XXIV IBERIAN LATIN-AMERICANCONGRESS ON COMPUTATIONALMETHODS IN ENGINEERING

2003

CIL

AM

CE

Ouro Preto/MG - Brazil

Geração de Superfícies de Interação pelo Método da Regressão Linear Múltipla com oModelo de Dano em Vigas de Timoshenko 3D - Aspectos Teóricos.

Pedro C. S. VieiraWilliam T. M. SilvaUniversidade de Brasília, Departamento de Engenharia Civil e Ambiental,PECC/FT, UnB/DF, Brasilpcvieira@unb.br, taylor@unb.brA. D. HanganuCentro Internacional de Métodos Numéricos en Ingeniería, Edificio C1,Campus Norte UPC, Gran Capitã, s/n, 08014, Barcelona - Españaalex@cimne.upc.es

Resumo. Na literatura técnica, são encontradas dificuldades para obtenção de superfícies deinteração em resultantes de tensões que abordem mais de três esforços seccionais, como porexemplo: momentos fletores, torçores, axiais e cortantes. Estes tipos de superfícies apresen-tam simplificações importantes nas análises de plasticidade, evitando o processo de integraçãonumérica das tensões nas seções transversais do elemento. Geralmente, as funções de escoa-mento f trabalham no espaço de tensões e dentro deste escopo, vê-se que a interação entreas tensões normal e tangencial pelo critério de Mises, aplicadas para os principais pontos detensão numa seção metálica, é usualmente considerada como um limite para projetos elásti-cos de elementos resistentes. Expressões em tensões, que dependem dos esforços dados pelaResistência dos Materiais, permitem aplicações de condições limites de forma direta. Quandoesta forma de critério é dada, a interação de surpefícies limites para trios de esforços aplicadosresulta em planos, quádricas, surperfícies mais complexas, ou uma mistura destas. Técnicasque usam formulações análiticas são mais ou menos complexas e dependem de características,como por exemplo: combinação de tensões ou de esforços seccionais, e o tipo de seção anal-isada. A formulação usada, neste trabalho, para a obtenção destas superfícies leva em conta aanálise de elementos sólidos de viga de Timoshenko 3D, no qual, os esforços seccionais usadoscomo base para gerar a superfície são retirados da seção mais plastificada de um elementoengastado livre. Assim, estes esforços são obtidos em função da combinação de diversos car-regamentos para a obtenção da superfície de interação entre eles. O modelo de regressão linearmúltipla permite o tratamento dos esforços resultantes de diversas análises para a geração dasuperfície que leva em conta a combinação de esforços tratada. Esta parte, aborda os aspectosteóricos da formulação.

Keywords: Funções de Escoamento, Curvas de Interação, Vigas de Timoshenko 3D, RegressãoLinear Múltipla, Esforços Seccionais.

1. INTRODUÇÃO

Pode-se ver, na literatura, que o enfoque dado para a análise não-linear de estruturas comvigas 3D, usando superfícies de interação, leva em conta, somente, os três esforços seccionais,devido as dificuldades de trabalhar com hiper-superfícies e modelar, experimentalmente, os seisesforços seccionais existentes na análise de pórticos espaciais. O uso de elementos sólidosaumenta o custo computacional da análise de forma importante quando são empregados paraanálises reais e geralmente os softwares comerciais para projetos de pórticos espaciais são feitospor meios elementos finitos de barra com formulações que simplificam a análise. Sabe-se que acapacidade de resistência de um elemento estrutural depende do tipo de esforços com o qual estásendo solicitado, ou seja, um sistema estrutural que trabalha só com esforços axiais tem umacapacidade de resistência maior do que quando trabalha com flexo-compressão, flexo-torção,etc. Estendendo o conceito, vê-se que os edifícios que trabalham com cargas de vento, concen-tradas, distribuídas, etc possuem um sistema que trabalha com os seis esforços seccionais; seanalisarmos a influência que um tipo de esforço tem em alterar a capacidade resistente da estru-tura, pode-se ver a necessidade de compreender melhor a interação que há entre os esforços, anível de influenciar a estabilidade global ou local de uma estrutura.

Quando é feita uma análise elasto-plástica com modelos de viga 3D, necessita-se de umafunção da superfície de escoamento que controlará o término da fase elástica e o estado plásticoda estrutura. O limite entre a zona elástica e a plástica se estabelece mediante a superfície defluência ou superfície de descontinuidade, e a partir deste limite esta superfície adquire mo-bilidade no espaço de tensões, seguindo a evolução do processo plástico, transformando-se nadenominada superfície de carga plástica. Para estabelecer, durante o processo de carga, o iníciodo comportamento inelástico e a posterior evolução das fronteiras do domínio elástico dentrodo espaço, adota-se o critério de fluência ou descontinuidade (Oller, 2001).

Os modelos que foram utilizados neste trabalho levam em conta a necessidade de umasolução de esforços seccionais como resultado de uma análise 3-D. São apresentados os fun-damentos da teoria de viga de Timoshenko 3D, das relações constitutivas, da regressão linearmúltipla e das curvas de interação, os quais foram utilizados no presente trabalho.

2. FORMULAÇÃO DE VIGA DE TIMOSHENKO 3D

O elemento finito de barra 3D foi desenvolvido a partir do elemento de viga de Timoshenko.Este elemento permite modelar o comportamento de um barra prismática sobre qualquer car-regamento (Hanganu, 1997). Se trata de um elemento finito lagrangiano de continuidade C0,de três nós e seis graus de liberdade por nó. O fato de que o modelo constitutivo necessita deinformação a nível de tensão-deformação, faz necessária uma discretização da seção transversalem uma malha retangular (ver fig. 1).

A formulação seguinte descreve as relações existentes entre as variáveis pontuais e sec-cionais. A barra é considerada no sistema de coordenadas locais, com seu eixo x0 formandocom os restantes do eixos um triedo direito. Considera-se que os eixos y0 e z0 são os eixosprincipais de inércia de cada seção. O convênio dos sinais para os deslocamentos e rotações é oda mecânica clássica (Hanganu, 1997). Em coordenadas locais, os campos de deslocamentos ede deformações são:

vw

u

x

y

z

Figura 1: Barra 3D com seção dividida mediante uma malha retangular. Eixos globais e locais.

up =

up

vp

wp

=

u0 + z0θy0 − y0θz0v0 − z0θx0w0 + y0θx0

=

1 0 0 0 z0 −y00 1 0 −z0 0 00 0 1 y0 0 0

u0

v0

w0

θx0θy0θz0

= Su0 (1)

ε =

εx0γx0y0γx0z0

=

∂up

∂x0∂up

∂y0+

∂vp

∂x0∂up

∂z0+

∂wp

∂x0

=

du0

dx0+ z

dθy0

dx0− y

dθz0

dx0dv0

dx0− θz0 − z

dθx0

dx0dw0

dx0+ θy0 + y

dθx0

dx0

= Sε (2)

onde ε =½

du0

dx0dv0

dx0− θz0

dw0

dx0+ θy0

dθx0

dx0dθy0

dx0dθz0

dx0

¾T

.

As variáveis das equações anteriores tem o seguinte significado:

• up− vetor de deslocamentos locais de um ponto qualquer da seção;

• ε− vetor de deformações;• u0− vetor de deslocamentos em coordenadas locais do elemento finito de barra 3D, cor-respondente ao eixo da seção;

• ε− vetor de deformações seccionais;• S− matriz geométrica de relação ponto-seção.

Utilizando o princípio dos deslocamentos virtuais para escrever as equações de equilíbrio,o trabalho interno Lint , na equação(3), correspondente a um campo de deformação virtual δεtoma a seguinte forma:

Lint =

ZV

δεTσdV =

ZV

δεTSTσdV =

Z l

0

δεT·Z

A

STσdA

¸dx (3)

=

Z l

0

δεT σdx0=

Z l

0

δεT Cεdx0 (4)

Na relação (3), é visto o vetor de esforços seccionais σ como o conjugado energético dovetor de deformações seccionais. A matriz S = S (y, z) varia nas duas direções da seção,sendo isto, de interesse para este trabalho porque permite fazer uma análise tridimensional doprocesso de plastificação:

σ =

ZA

STσdA =

ZA

1 0 0 0 z0 −y0

0 1 0 −z0 0 00 0 1 y

00 0

σx0

τx0y0

τx0z0

dA

=

ZA

©σx0 τx0y0 τx0z0 −z0τx0y0 + y

0τx0z0 z

0σx0 −y0σx0

ªTdA

=©Nx0 Qy0 Qz0 Tx0 My0 Mz0

ªT (5)

Onde σ = Cε, sendo C o tensor constitutivo de rigidez local, cujo valor é:

C =

E 0 00 G 00 0 G

(6)

Empregando as relações da equação (2) junto com a equação (5) , obtem-se:

σ =

ZA

STσdA =

ZA

STCεdA =

ZA

STCSεdA = Cε (7)

Baseando-se na equação (7) , apresenta-se o valor da matriz constitutiva seccional C querelaciona as deformações e os esforços seccionais.

C =

ZA

STCSdA (8)

A equação (8) é resolvida com a integração sobre uma quadrícula de suporte 3D. O processode maneira detalhado é apresentado, a seguir:

C=

ZA

E 0 0 0 z

0E −y0E

0 G 0 −z0G 0 00 0 G y

0G 0 0

0 −z0G y0G

¡z02 + y

02¢G 0 0

z0E 0 0 0 z

02E −y0z0E−y0E 0 0 0 −y0z0E y

02E

dA

=nXi=1

Z y0i+1

y0i

Z z0i+1

z0i

Ei 0 0 0 z

0Ei −y0Ei

0 Gi 0 −z0Gi 0 00 0 Gi y

0Gi 0 0

0 −z0Gi y0Gi

¡z02 + y

02¢Gi 0 0

z0Ei 0 0 0 z

02Ei −y0z0Ei

−y0Ei 0 0 0 −y0z0Ei y02Ei

dy0dz

0

(9)

As quadrículas do modelo 3D, (ver fig. 1) , são retângulos com os lados paralelos aos eixosde inércia, as integrais duplas podem ser integradas de maneira independente, em função decada variável. O elemento finito está definido por 3 nós com 6 graus de liberdade cada um(Hanganu, 1997). São empregadas as funções de forma e o vetor de deslocamentos nodaisseguintes:

N0i = NiI6, com N

0 =©N

01 N

02 N

03

ªa0i =

©u0i v

0i w

0i θx0 i θy0 i θz0 i

ªT , com a0 = © a01 a02 a

03

ª(10)

Onde I6 é uma matriz identidade de fila (rank) 6,N0 são as funções de forma e a0i é o vetor

de deslocamentos nodais.

A matriz B0 é apresentada a seguir:

B0i =

dNi

dx00 0 0 0 0

0dNi

dx00 0 0 −Ni

0 0dNi

dx00 Ni 0

0 0 0dNi

dx00 0

0 0 0 0dNi

dx00

0 0 0 0 0dNi

dx0

(11)

Sendo que: B0=©B

01 B

02 B

03

ª.

A transformação do sistema local ao global de coordenadas é obtida através da matriz detransformação T0, definida a seguir:

a0i = Tai ;

T =

·T0 0303 T0

¸;

T0 =

cos ¡x0 , x¢ cos¡x0, y¢cos¡x0, z¢

cos¡y0, x¢cos¡y0, y¢cos¡y0, z¢

cos¡z0, x¢cos¡z0, y¢cos¡z0, z¢ (12)

Onde: cos¡x0, x¢é o coseno do ângulo entre a direção local x0 e a direção global x, e assim

sucessivamente para os demais (Hanganu, 1997).

As funções de forma empregadas na teoria 3D definem o campo contínuo elementar, in-terpolando os valores nodais. São utilizadas funções quadráticas lagrangianas correspondentesa um elemento de barra de três nós (ver fig. 2). Cada nó tem uma função de forma associada,de maneira que esta vale 1 no nó e 0 nos demais. Expressam-se como função de uma variávelnormalizada ζ, que varia de -1 a 1 (Hanganu, 1997). A seguir, são apresentados seus valores:

Figura 2: Representação das funções de forma

N1 (ζ) =1

2ζ (ζ − 1) ; N2 (ζ) = 1− ζ2 ; N3 (ζ) =

1

2ζ (ζ + 1) (13)

As matrizesN0i eB

0i são reapresentadas em função da matriz transformaçãoT ( ver eq. 12):

Ni = N0iT ; Bi = B

0iT (14)

O vetor de deslocamentos seccionais da barra nos eixos locais u0 e as deformações sec-cionais ε são apresentados em função da matriz de formaN e de B, respectivamente:

u0= Na ; ε = Ba (15)

A expressão desenvolvida para o vetor de forças internas elásticas Fe é vista a seguir:

Fe =

Zl

BT σdx =

Zl

BT Cεdx0 (16)

Onde l é o comprimento do elemento finito.Com esta útlima transformação a força elástica toma a seguinte forma:

Fe =

Zl

BT CBdx0a = Ka (17)

Pode-se ver, na equação (17), a matriz de rigidez global do elemento, K. O processode integração é feito com uma quadratura gaussiana reduzida de 2 pontos para os termos decortante, para evitar o efeito de bloqueo de cortante (Hanganu, 1997).

O processo de cálculo não-linear consiste na avaliação das deformações seccionais ε corre-spondentes aos deslocamentos a, tal como se pode ver na equação (15). Desta maneira avalia-seas deformações pontuais εmediante a equação (2) e as tensões correspondentes que são corregi-das dentro do modelo constitutivo para depois integrá-las sobre a seção mediante a relação (5).No final, são obtidos os esforços seccionais correspondentes de maneira que pode-se calcularas forças residuais com algoritmos normais (Hanganu, 1997).

A seção da viga está discretizada mediante uma quadrícula ortogonal (ver fig. 1). Os eixosda quadrícula devem ser paralelos as direções principais de inércia da seção. Cada retângulo daquadrícula pode estar caracterizado por um material e dimensões geométricas distintas, sendoque neste trabalho o material definido é homogêneo. Os quatros cantos de cada retângulo sãoos pontos de cálculo das deformações e tensões. Para integrar as tensões seccionais a partirdas tensões do modelo constitutivo, considera-se que todas as tensões envolvidas tenham umavariação linear dentro de uma célula da quadrícula (Hanganu, 1997).

3. MODELO DE DANO ISOTRÓPICO

O programa de análise 3D usa o modelo constitutivo de dano isótropico. O modelo éusado para problemas térmicamente estáveis, na configuração material lagrangiana com peque-nas deformações e deslocamentos (Hanganu, 1997). Define-se, a seguir, a variável de dano dnassociada a uma superfície elementar com um volume de material degradado.

dn =Sn − Sn

Sn= 1− Sn

Sn(18)

Onde: Sn é a área total da seção, Sn é a área resistente efetiva, e¡Sn − Sn

¢é a área ocupada

pelas aberturas.Na relação anterior, dn representa a densidade dos defeitos do material e terá o valor zero

no estado inicial, sem dano. A medida que a fissuração avança, dn tenderá a um valor crítico,próximo da unidade que corresponde a completa falta de área resistente Sn. A relação deequilíbrio entre a tensão de Cauchy σ e a tensão efetiva σ é vista a seguir:

σSn = σSn (19)

Usando as equações (18) e (19), obtem-se:

σ = (1− dn) σ = (1− dn)Eε (20)

Durante um processo de degradação em desenvolvimento, é a área efetiva que suportaa carga externa, sendo assim, σ é um parâmetro fisicamente mais representativo que σ. Osmodelos de dano descrevem o comportamento não-linear mediante uma ou várias variáveisinternas de dano, que medem a perda de rigidez secante do material e que se normalizam comrespeito a unidade, a cual corresponde o dano máximo. O efeito do dano se traduz na diminuiçãodo módulo de rigidez secante. O modelo considerado na tese doutoral de Hanganu foi o de danobaseado na mecânica do sólido deformável com somente uma variável interna (Hanganu, 1997).

Dentro destes conceitos foram feitas adequações para a aplicação da teoria na análise demateriais homogêneos, como os metais. Definiu-se a energia de fratura com um valor alto parapoder obter o comportamento dos metais. Foram feitos testes para os esforços axiais, momen-tos fletores e torçores (ver fig. 3) para testar a viabilidade do modelo de dano em materiaishomogêneos, especificamente o aço. Chegou-se a conclusão que a teoria podia ser aplicadaperfeitamente nas análises. Nos exemplos numéricos, serão explicados, com mais detalhes, ashipóteses adotadas nas análises através do elemento de viga 3D.

Figura 3: Gráfico da relação carga × deslocamento normalizados.

4. REGRESSÃO LINEARMÚLTIPLA

Muitas aplicações da análise de regressão envolvem situações em que há mais de umavariável de regressão. Um modelo de regressão que contém mais de um regressor recebe onome de modelo de regressão múltipla (Montegomery et. all, 1998). Um modelo de regressãomúltipla pode ser escrito como a relação seguinte:

Y = β0 + β1x1 + β2x2 + · · ·βkxk+ ∈ (21)

Este é um modelo de regressão com várias variáveis, sendo que Y é a variável depen-dente ou resposta, e pode estar relacionada com k variáveis independentes ou regressores. Osparâmetros βj , j = 0, 1, · · · , k, se conhecem como coeficientes de regressão. Este modelo de-screve um hiperplano no espaço de dimensão k formado pelas variáveis de regressão {xj}. Oparâmetro βj representa a variação esperada na reposta Y por unidade de variação em xj quandotodos os demais regressores xi (i 6= j) se mantém constantes. Frequentemente estes modelos seempregam como funções de aproximação e se desconhece a verdadeira relação funcional entreY e x1, x2, . . . , xk. Sobre certos tipos de variáveis independentes o modelo de regressão linearconstitue uma aproximação adequada (Montegomery et. all, 1998). Os modelos que tem umaestrutura mais complexa que a dada pela equação (21) com frequência, também, podem seranalisadas com as técnicas da regressão linear múltipla. Por exemplo, considerando um modelopolinomial cúbico com uma variável de regressão.

Y = β0 + β1x+ β2x2 + β3x

3+ ∈ (22)

Tomando-se x1 = x, x2 = x2, x3 = x3, então a equação (22) pode ser escrita da formausual do modelo de regresão múltipla.

Os modelos que incluem efeitos de interação, que é o caso deste trabalho, também podemser analisados com os métodos da regressão linear múltipla. Uma interação entre duas variáveispode ser representada como um produto entre variáveis, tal como

Y = β0 + β1x1 + β2x2 + β12x1x2+ ∈ (23)

Faz-se as seguintes modificações: x3 = x1x2 e β3 = β12, então a equação (23) pode serescrita como

Y = β0 + β1x1 + β2x2 + β3x3+ ∈ (24)

que é um modelo de regressão linear múltipla.Note-se que, ainda que este seja um modelo de regressão linear, a forma da superfície

gerada pelo modelo não é linear. Em geral, qualquer modelo de regressão que é linear nosparâmetros (β) é um modelo de regressão linear, sem importar a forma de superfície que estegera (Montegomery et. all, 1998).

4.1 Enfoque matricial para a regressão

E mais conveniente expressar o modelo com operações matemáticas em forma matricial.Suponha-se que existem k variáveis de regressão e n observações (xi1, xi2, . . . , xik, yi), i =1, 2, . . . , n, e que o modelo que relaciona os regressores com a resposta seja:

Y = β0 + β1xi1 + β2xi2 + · · ·βkxik+ ∈i, i = 1, 2, . . . , n (25)

Este modelo é um sistema de n equações que pode expressar-se em notação matricial como

y = Xβ + ∈ (26)

onde:

y =

y1y2...yn

X =1 x11 x12 · · · x1k1 x21 x22 · · · x2k...

......

...1 xn1 xn2 · · · xnk

(27)

β =

β0β1...βk

e ∈ =∈1∈2...∈n

(28)

Em geral, y é um vetor de observações de (n× 1), X é um tensor (matriz) de (n× p)dos níveis das varíaveis independentes, β é um vetor de (p× 1) formado pelos coeficientes deregressão e ∈ é um vetor (n× 1) dos erros aleatórios.

Deve-se encontrar o vetor dos estimadores dos mínimos quadrados, β, que minimiza

L =nXi=1

∈2i= ∈T∈ = (y−Xβ)T (y−Xβ) (29)

O estimador de mínimos quadrados β é a solução para β nas equações

∂L

∂β= 0 (30)

Desevolvendo-se a equação (30) chega-se a:

XTXβ = XTy (31)

As equações (31) são as equações normais dos mínimos quadrados em forma matricial, esão identicas a forma escalar, como é apresentado a seguir:

nβ0 + β1

nXi=1

xi1 + β2

nXi=1

xi2 + · · ·+ βk

nXi=1

xik =nXi=1

yi

β0

nXi=1

xi1 + β1

nXi=1

x2i1 + β2

nXi=1

xi1xi2 + · · ·+ βk

nXi=1

xi1xik =nXi=1

xi1yi

......

......

...

β0

nXi=1

xik + β1

nXi=1

xikxi1 + β2

nXi=1

xikxi2 + · · ·+ βk

nXi=1

x2ik =nXi=1

xikyi (32)

Para resolverem as equacões normais, multiplicam-se ambos membros da equação (31)pela inversa deXTX. Por conseguinte, o estimador de mínimos quadrados de β é:

β =¡XTX

¢−1XTy (33)

Note-se que existem p = k + 1 equações normais e p = k + 1 incógnitas, ou seja, osvalores de β0, β1, · · · , βk. Por outro lado, a matrizXTX não é singular, de modo que podem-seempregar os métodos de inversão de matrizes que existem na literatura.

A forma matricial das equações normais de (32) é apresentada a seguir:

n

Pni=1 xi1

Pni=1 xi2 · · · Pn

i=1 xikPni=1 xi1

Pni=1 x

2i1

Pni=1 xi1xi2 · · · Pn

i=1 xi1xik...

......

...Pni=1 xik

Pni=1 xikxi1

Pni=1 xikxi2 · · · Pn

i=1 x2ik

β0β1...βk

=

Pn

i=1 yiPni=1 xi1yi...Pn

i=1 xikyi

(34)

Pode-se observar que a matriz XTX é uma matriz simétrica de (p× p), e que XTy é umvetor coluna de (p× 1). Os elementos da matriz diagonal deXTX são as somas dos quadradosdos elementos nas colunas deX, enquanto que os elementos que estão fora da diagonal principalsão as somas dos produtos cruzados dos elementos das colunas de X (Montegomery et. all,1998). Os elementos de XTy são as somas dos produtos cruzados das colunas de X e asobservações de y.

O modelo de regressão ajustado tem a seguinte forma:

yi = β0 +nX

j=1

βjxij , i = 1, 2, · · · , n (35)

A forma matricial do modelo é:

y = Xβ (36)

A diferença entre a observação yi e o valor ajustado yi é um resíduo, ei = yi − yi. O vetorde resíduos de (n× 1) se denota como:

e = y− y (37)

4.2 Propriedades dos estimadores de mínimos quadrados

O resíduo quadrático médio σ2 está dado pelo erro ou resíduo quadrático médio

σ2 =MSE =SSEn− p

(38)

Onde SSE é a soma dos quadrados dos erros, sendo representado a seguir:

SSE =nXi=1

(yi − yi)2 =

nXi=1

e2i = eTe (39)

Substituindo-se e = y− y = y−Xβ, chega-se a:

SSE =³y−Xβ

´T ³y−Xβ

´= yTy− βT

XTy− yTXβ + βTXTXβ

= yTy−2βTXTy + β

TXTXβ (40)

Sendo queXTXβ =XTy.

A equação anterior se converte em:

SSE = yTy− βT

XTy (41)

Por conseguinte, outra maneira de escrever a equação (38) é apresentada a seguir:

σ2 =SSEn− p

=yTy− βT

XTy

n− p(42)

4.3 Prova de hipótese na regressão

Nos problemas de regressão linear múltipla, existem certas provas de hipótese sobre osparâmetros do modelo que são utéis para medir a adequação do modelo. A prova de hipóteserequer que os termos de erro ∈i do modelo de regressão tenham distribuições normais e inde-pendentes com média zero, e variância σ2. A seguir, é apresentada uma tabela com o modelopara a análise da variância para a prova de significância da regressão.

Tabela 1: Tabela da prova de significância

Fonte Soma de Graus de Média dede variação quadrados liberdade quadrados F0 P > F0Regresão SSR k MSRErro ou resíduo SSE n− p MSE

MSRMSE

Total Syy n− 1

A prova de significância para a regressão é uma prova para determinar se existe uma relaçãolinear entre a variável linear e a variável de resposta e um sub-conjunto de variáveis de regressãox1, x2, . . . , xk. As hipóteses apropriadas são:

H0 : β1 = β2 = · · · = βk = 0

H1 : βj 6= 0, pelo menos para um j (43)

Se H0 : β1 = β2 = · · · = βk = 0 for rejeitado implica que pelo menos umas dasvariáveis de regressão x1, x2, . . . , xk tem contribuição significativa no modelo. A soma totaldos quadrados Syy divide-se em: uma soma de quadrados devida a regressão e uma soma dequadrados devido ao erro.

Syy = SSR + SSE (44)

Se H0 : β1 = β2 = · · · = βk = 0 é verdadeira, então SSRσ2é uma variável aleatória ji-

quadrada com k graus de liberdade. O número de graus de liberdade para esta variável aleatóriaji-quadrada é igual ao número de variáveis de regressão do modelo. Pode-se demonstrar queSSEσ2é uma variável aleatória con n− p graus de liberdade, e que SSE e SSR são independentes

(Montegomery et. all, 1998). O estatístico de prova para H0 : β1 = β2 = · · · = βk = 0 é:

F0 =SSRk

SSE(n−p)

(45)

Deve-se rejeitar H0, si o valor calculado do estatístico de prova da equação (45), F0, émaior do que fα,k,n−p.

As fórmulas computacionais empregadas são:

SSE = yTy− βT

XTy (46)

SSR = βTXTy−

µnPi=1

yi

¶2n

(47)

4.4 Prova sobre os coeficientes individuais de regressão

Existe, também, o interesse de fazer provas de hipóteses sobre os coeficientes de regressão.Tais provas são úteis para determinar o valor potencial de cada uma das variáveis de regressão.O modelo pode se tornar mais eficaz com a inclusão de outras variáveis ou a eliminação de umou mais regressores presentes no modelo. A adição de uma variável ao modelo de regressãosempre faz com que a soma dos quadrados da regressão aumente e que a soma dos quadradosdo erro diminua. Portanto, deve-se decidir se o aumento na soma dos quadrados da regressão ésuficientemente grande para justificar o uso de uma variável a mais no modelo (Montegomeryet. all, 1998). Por outra parte a adição de uma variável sem importância pode aumentar o erroquadrático médio, indicando que a varíavel diminue a qualidade com que o modelo ajusta osdados.

A hipótese para a prova de significância de qualquer coeficiente de regressão individual,por exemplo βj , é:

H0 : βj = 0

H1 : βj 6= 0 (48)

A não é rejeição da opção H0 : βj = 0, indica que o regressor xj pode ser eliminado domodelo. O estatístico de prova para esta hipótese é:

T0 =βj

2pσ2Cjj

(49)

Onde Cjj é o elemento da diagonal de¡XTX

¢−1que corresponde a βj , e o denominadorda equação (49) é o erro normalizado do coeficiente de regressão βj .

A hipótese nulaH0 : βj = 0 é rejeitada se |t0| > tα2,n−p. Isto é chamado de prova parcial o

marginal porque o coeficiente de regressão βj depende de todas as demais variáveis de regressãoxi (i 6= j) que estão no modelo (Montegomery et. all, 1998).

4.5 Coeficiente de determinação múltipla

Para medir a adequação do modelo, podem ser empregadas várias técnicas. Dentre estas,se apresenta o coeficiente de determinação múltipla R2 definido como:

R2 =SSRSyy

= 1− SSESyy

(50)

R2 é uma medida da magnitude da redução na variabilidade de y obtida mediante o em-prego das variáveis de regressão x1, x2, · · · , xk, com 0 ≤ R2≤ 1. Um valor grande de R2 não

indica necessáriamente que seja um bom modelo. A adição de uma variável ao modelo sempreaumenta o R2, sem importar se a variável é ou não estatísticamente significativa. Portanto, sãonecessárias análises conjuntas de outras informações para determinar a competência do modelo.A raiz quadrada positiva de R2 recebe o nome de coeficiente de correlação múltipla entre y e oconjunto de variáveis de regressão x1, x2, · · · , xk, ou seja,R é uma medida da associação linearque existe entre y e x1, x2, · · · , xk (Montegomery et. all, 1998).

Outro critério similar ao R2 é o coeficiente R2 ajustado que leva em conta o número devariáveis do modelo. Este coeficiente é definido como:

R2 = 1− n− 1n− p

¡1−R2

¢(51)

Reapresentado a equação (51) tem-se que:

R2 = 1− n− 1n− p

¡1−R2

¢= 1− n− 1

n− p

µSSESyy

¶= 1− n− 1

Syy(MSE) (52)

Pode-se perceber que R2 pode diminuir a medida que p aumenta se a redução em (n−1)(1−R2) não

é compensada pela perda de um grau de liberdade em n− p. O normal é que o experimentadorselecione o modelo de regressão que tenha o valor máximo de R2. Entretanto, ao fazer isto éequivalente ao modelo que minimizaMSE (ver 52 ) .

5. CURVAS DE INTERAÇÃO

A obtenção de curvas de interação em resultantes de tensões facilita a análise, de maneiraque evita o processo de integração numérica. Para a obtenção das curvas de interação, em re-sultantes de tensões, foram utilizados os resultados dos esforços seccionais da análise numérica3D, apresentada anteriormente. Dentro do processo, foram feitas várias combinações de car-regamentos de forma a ter um grupo de pontos para gerar a superfície proposta, ou seja, pontosque tenham alcançados a superfície de escoamento. Para um dado carregamento, obtem-seum ponto, como por exemplo o ponto 1, coordenadas (n1,m1) , na figura (4) .As coordenadasrepresentadas são as combinações dos esforços adimensionalisados, como é visto a seguir naforma matricial:

Figura 4: Pontos gerados para criar a função de escoamento (caso uniaxial).

x1j =

"1

õN

Nu

¶α1!12

õMx

Mxu

¶α2!13

õMy

Myu

¶α3!14

õMz

Mzu

¶α4!15

· · ·#

=£1 (nα1)12 (mα2

x )13¡mα3

y

¢14

(mα4z )15 · · · ¤ (53)

, com j = 1, 2, · · · , kUm modelo que leva em conta combinações dos esforços seccionais independentes para

pórticos espaciais é apresentado a seguir:

xj1 =

µN

Nu

¶α1

, xj2 =

µMx

Mxu

¶α2

, xj3 =

µMy

Myu

¶α3

, xj4 =

µMz

Mzu

¶α4

,

xj5 =

µN

Nu

¶α5µMx

Mxu

¶α6

, xj6 =

µN

Nu

¶α7µMy

Myu

¶α8

,

xj7 =

µN

Nu

¶α9µMz

Mzu

¶α10

, xj8 =

µMx

Mxu

¶α11µMy

Myu

¶α12

,

xj9 =

µMx

Mxu

¶α13µMz

Mzu

¶α14

e xj10 =µMy

Myu

¶α15µMz

Mzu

¶α16

(54)

As observações de (54) , apresentadas para a regressão linear múltipla (27) , na forma ma-tricial são:

X =

1 x11 x12 · · · x1k1 x21 x22 · · · x2k...

......

...1 xn1 xn2 · · · xnk

(55)

Onde os termos xij são os esforços seccionais adimensionalisados, (54), vistos anterior-mente. Os αi são as constantes que determinam o grau da função;

©N, Mx, My, Mz

ªe©Nu, Mxu, Myu, Mzu

ªsão os esforços de cálculo e límite elasto-plástico, respectiva-

mente.A seguir, tem-se um modelo de equação geral para as curvas de interação levando em conta

os esforços seccionais independentes:

1 = β0 + β1

µN

Nu

¶α1

+ β2

µMx

Mxu

¶α2

+ β3

µMy

Myu

¶α3

+ β4

µMz

Mzu

¶α4

+ β5

µN

Nu

¶α5µMx

Mxu

¶α6

+ β6

µN

Nu

¶α7µMy

Myu

¶α8

+ β7

µN

Nu

¶α9µMz

Mzu

¶α10

+ β8

µMx

Mxu

¶α11µMy

Myu

¶α12

+ β9

µMx

Mxu

¶α13µMz

Mzu

¶α14

+ β10

µMy

Myu

¶α15µMz

Mzu

¶α16

(56)

Reapresentando a equação (56) e adotando os seguintes processos, chega-se a equação(58):

x1 =

µN

Nu

¶α1

, x2 =

µMx

Mxu

¶α2

, x3 =

µMy

Myu

¶α3

, x4 =

µMz

Mzu

¶α4

,

x5 =

µN

Nu

¶α5µMx

Mxu

¶α6

, x6 =

µN

Nu

¶α7µMy

Myu

¶α8

,

x7 =

µN

Nu

¶α9µMz

Mzu

¶α10

, x8 =

µMx

Mxu

¶α11µMy

Myu

¶α12

,

x9 =

µMx

Mxu

¶α13µMz

Mzu

¶α14

e x10 =

µMy

Myu

¶α15µMz

Mzu

¶α16

(57)

comµN

Nu

¶= n,

µMx

Mxu

¶= mx,

µMy

Myu

¶= my e

µMz

Mzu

¶= mz

Então, chega-se a:

1 = β0 + β1n+ β2mx + β3my + β4mz + β5nmx

+β6nmy + β7nmz + β8mxmy + β9mxmz + β10mymz

sendo reapresentada na forma corrente de regressão, a seguir:

1 = β0 + β1x1 + β2x2 + β3x3 + β4x4 + β5x5

+β6x6 + β7x7 + β8x8 + β9x9 + β10x10 (58)

A equação anterior mostra a viabilidade do método de regressão (ver eq. 35), de maneiraque este modelo pode ser extendido para análises mais gerais com os seis esforços seccionais evários tipos de curvas de interação que atendam aos critérios de uma superfície de escoamento.

5.1 Limites plásticos dos esforços seccionais

Os limites plásticos foram adotados, para uma seção retangular, em função das fórmulasexistentes na literatura. São apresentados, a seguir, as fórmulas empregadas:

Axial

Nu = Aσu (Tração) ;

Nu = Aσuω(Compressão) . (59)

Onde ω é o coeficiente de flambagem, A é a área da seção e σu é a resistência de cálculo.

σu =σeγa

(60)

Onde σe é límite elástico do aço e γa é o coeficiente de minoração com os seguintes valores;

• γa = 1 para aços com límite elástico mínimo garantido, e

• γa = 1, 1 para aços cujo límite elástico seja determinado por métodos estatísticos.

Os coeficientes de flambagem ω dependem da função de esbeltez λ =l

ique podem ser

vistas na literatura, como por exemplo: (NBE EA-95, 2001), com l, sendo o comprimento doelemento e i o raio de giração.

Momento fletor

Miu = σu Wiψi; i = y, z (61)

ψi =SiWi

(62)

OndeWi é o momento resistente mínimo da seção no plano de flexão, ψi é o fator de forma(coeficiente que depende da forma da seção) e Si é a soma dos momentos estáticos com respeitoao eixo neutro plástico. Todos estes dados são retirados de tabelas existentes na literatura, comopor exemplo: (NBE EA-95, 2001).

Momento torçor

Mxu =1

6kb2 (3h− b) (63)

Onde k = σu, sendo que b é a largura da seção e h a altura.

5.2 Funções de escoamento em resultantes de tensão

Na literatura, são apresentadas algumas funções de escoamento em resultantes de tensãocom procedimentos aproximados para alguns casos e em outros análiticos. Dentro deste escopo,será apresentado um resumo de algumas funções existentes, sendo que as definições vistas nasequações de (53) a (57) serão empregadas na maneira de apresentar as funções de escoamento.Existem outras funções, para outros tipos de seções com combinações de esforços seccionais,porém, neste trabalho, os exemplos numéricos serão comparados com as funções apresentadas,de maneira a comprovar a eficácia do método proposto.

Interação momento fletor e força axial Neste caso, a deformação é assumida de maneiraque a tensão de escoamento é alcançada para toda seção, com uma relação tensão-deformaçãoperfeitamente plástica (Crisfield, 1990). Com isto, chega-se as seguintes equações:

n =N

Nu=

N

σut= 1− 2η (64)

Tomando-se os momentos sobre o centro da viga (ver fig. 5), chega-se a:

mi =Mi

Miu=4Mi

σut2= 4

¡η − η2

¢; i = y, z (65)

Eliminando-se a altura adimensional, η, nas equações (64) e (65), obtem-se:

Figura 5: Tensão e deformação para o caso uniaxial. (a) deformação; (b) tensão; (c) resultantes detensão; (d) deformação reversa; (e) deformação no plano dominante.

f = n2 +mi − 1 = 0, i = y, z (66)

Pode-se ver que a função de escoamento da equação (66) é perfeitamente aplicada na re-gressão linear múltipla, ver equação (56). Alterando a função anterior para a deformação re-versa (Crisfield, 1990), obtem-se:

f = n2 −mi − 1 = 0, i = y, z (67)

Combinando-se as duas funções (66) e (67), chega-se a:

f = n2 + smi − 1 = 0 (68)

Onde: s =Mi

|Miu| ; i = y, z

A seguir, são vistas funções de aproximação para a equação (68), dadas por:

f = n2 +s√3min+m2

i − 1 = 0 (69)

f = n2 + 3smin+9

4m2

i − 1 = 0 (70)

com i = y, z.Estas funções, apresentadas por Crisfield, são melhores para reservados critérios de escoa-

mentos de seção cheia com rápidas análises aproximadas (Crisfield, 1990).Lubliner (Lubliner, 1990), também, apresenta para o caso de uma viga retangular com

largura b e altura h, e chegam-se as seguintes equações:

Mi = σyb

·h2

4− y0

¸, i = y, z;

N = 2σyby0. (71)

Onde: σy= tensão última e y0= coordenada onde se anula o centro elástico.A seguir, apresenta-se a equação (71) na forma adimensional:

mi = 1− η2;

n = η (72)

com η = 2y0h.

E, também, na forma explicíta:

mi = 1− n2 (73)

Onde: mi =Mi

Miucom i = y, z; e n =

N

Nu.

Interação momento fletor, força axial e torçor Uma barra com seção retangular com umacombinação de força axial, momento fletor e torçor apresenta a seguinte função de escoamento(Lubliner, 1990):

my =p1−m2

x −n2xp1−m2

x

(74)

Reapresentando a equação anterior, chega-se a:

m2y + 2myn+ n2 +m2

x = 1 (75)

Interação de momentos fletores A função apresentada, a seguir, é para o caso de uma barracom seção retangular e foi definida por Lubliner da seguinte forma:

mz +3

4m2

y = 1,my

mz≤ 1; (76)

my +3

4m2

z = 1,my

mz≥ 1. (77)

Esta parte enfoca a formulação teórica dos modelos usados no trabalho, sendo que osexemplos numéricos com os resultados e conclusões são apresentados na segunda parte, notrabalho: Geração de Superfícies de Interação pelo Método da Regressão Linear Múltipla como Modelo de Dano em Vigas de Timoshenko 3D - Exemplos Numéricos.

Referências

Argyris, J., 1982, "An Excursion into Large Rotations", Comp. Methods Appl. Mech. Engrg.,vol. 32, pp. 85-155.

Atluri, S. N.,1983, "Alternate Stress and Conjugate Strain Measures and Mixed Variational For-mulations Involving Rigid Rotations, for Computational Analyses of Finitely DeformedSolids, with Application to Plates and Shells. Theory ", Computer.&Structures, vol.18(1),pp. 93-106.

Chen, W. F.; Astuta, T., 1977, "Theory of Beam-Columns", Space Behaviour and Design,McGraw-Hill, New York, vol 2.

Crisfield, M. A. , 1990, "A Consistent Co-rotational Formulation for Non-linear, Three-dimensional, Beam-elements", Comp.Methods Appl. Mech. Engrg., vol. 81, pp. 131-150.

Crivelli, L.A., 1991, "A Total-lagragian Beam Element for Analysis of Nonlinear Space Struc-tures", PhD Tesis, CSSC, College of Engineering University of Colorado, Boulder, Col-orado.

Faria, H.P., 1998, "Análise Não-linear de Instabilidade Elástica de Pórticos Planos", Dissertaçãode Mestrado, ENC/FT/UnB.

Felippa, C. A.; Militello, C. , 1990, "Variational Formulation of High-performance Finite-elements Parametrized Variational-principles", Comp. Structures, vol. 36 (1), pp. 1-11.

Felippa, C.A., 2001,"Customizing High Performance Elements by Fourier Methods", Trends inComputational Structural Mechanics, CIMNE, Barcelona, Spain.

Gerardin, M; Cardona, A., 1987, "Kinematics and Dynamics of Rigid and Flexible MechanicsUsing Finite Elements and Quartenion Algebra", Comp. Mech. In Press.

Gere, J. M.; Weaver, W. JR., 1987, "Análise de Estruturas Reticulares", Guanabara S.A.

Hanganu, A. D., 1997, "Metodologia de Evaluación del Deterioro en Estructuras de HormigónArmado", Monografia CIMNE, no 39, Barcelona, Spain.

IMSL, 1997, "IMSL Fortran and C Application Development Tools", Visual Numerics, Inc.

Irles, R.M.; Irles, F.M., 2000, "Elastic Interaction Graphs for Steel H-sections Subjected toBending, Shear and Axial Forces", International Journal of Solids and Structures., vol.37, pp. 1327-1337.

Irles, R.M.; Irles, F.M., 2001, "Biaxial Bending-axial Force Elastic Interaction Diagrams in Hol-low Steel Sections", International Journal of Solids and Structures., vol. 38, pp. 423-433.

Kozar, I.; Ibrahimbegovic, A., 1995, "Finite Element Formulation of The Rotation Solid Ele-ment", Finite Elements in Analysis and Design., vol. 20, pp. 101-126.

Kondoh, K.; Tanaka, K.; Atluri, S. N. , 1986, "An Explicit Expression for Tangent-stiffness ofa Finitely Deformed 3-D Beam and its Use Analisys of Space Frames", Comp. Struct.,vol. 24, pp. 253-272.

Kondoh, K.; Atluri, S. N., 1987, "Large-deformation, Elasto-plastic Analysis of Frames UnderNonconservative Loading, Using Explicit Derived Tangent Stiffness Based on AssumedStress", Computational Mechanics, vol. 2, pp. 1-25.

Krenk, S.; Vissing-Jorgensen; Thesbjerg, C. L., 1999, "Efficient Collapse Analysis Tecniquesfor Framed Structures", Comp.Structures, vol. 72, pp. 481-496.

Lekhnitskii, S.G., 1963, "Theory of Elasticity of Anisotropic Elastic Body", Hoden Day, SanFrancisco.

Li, M., 1997, "The Finite Deformation Theory for Beam, Plate and Shell: Part I. TheTwo-dimensional Beam Theory", Comput. Methods Appl. Mech. Engrg., vol. 146, pp. 53-63.

Li, M., 1998, "The Finite Deformation Theory for Beam, Plane and Shell Part. III. The Three-dimensional Beam Theory and the FE Formulation", Comput. Methods Appl. Mech.Engrg., vol. 162, pp. 287-300.

Lubliner, J., 1990, "Plasticity Theory", Macmillan Publishing Company., New York, USA.

Menezes, L.F.; Teodosiu, C., 2000, "Three-dimensional Numerical Simulation of The Deep-drawing Process Using Solid Finite Elements",Journal of Materials Processing Technol-ogy., vol. 97, pp. 100-106.

Montegomery, D. C.; Runger, G. C., 1998, "Probabilidad y Estadística Aplicadas a la Inge-niería", McGraw-Hill Interamericana Editores, S. A., D. F., México.

NBE EA-95, 2001, "Norma Básica de la Edificación NBE EA-95: Estructuras de Acero enEdificación", Dirección General de la Vivienda, la Arquitectura y el Urbanismo, Madrid.

Oller, S., 2001, "Fractura Mecánica. Un enfoque global", CIMNE, Barcelona, Spain.

Oñate, E., 1992, "Cálculo de Estructuras por el Método de Elementos Finitos-Análisis estáticolineal", CIMNE, Barcelona, Spain.

Orbison, J.G.; McGuire, W. ; Abel, J.F., 1982, "Yield Surface Aplications in Non-linear SteelFrame Analysis", Comp. Methods Appl. Mech. Engrg., vol. 33, pp. 557-573.

Park, M.S.; Lee, B.C., 1996, "Geometrically Non-linear and Elastoplastic Three-dimensionalShear Flexible Beam Element of Von-mises-type Hardening Material", Int. j. numer.Methods eng., vol. 39, pp. 383-408.

Rathod, H. T. ; Sridevi, K. , 2001, "General Complete Lagrange Interpolations with Appli-cations to Three-dimensional Finite Element Analysis", Comp. Methods Appl. Mech.Engrg., vol. 190, pp. 3325-3368.

Saje, M., 1990, "A Variational Principle for Finite Planar Deformation of Straight Slender Elas-tic Beams", Internat. J. Solids Structures, vol. 26, pp. 887-900.

Saje, M.; G. Jeleni , 1995, "A Kinematically Exact Space Finite Strain Beam Model-finiteElement Formulation by Generalized Virtual Work Principle", Comp. Methods Appl.Mech. Engrg., vol. 120, pp. 131-161.

Shi, G.; Atluri, S. N., 1988, "Elasto-plastic Large Deformation Analisys of Space-frames: APlastic-hinge and Stress-based Explicit Derivation of Tangent Stiffness", Int. j. numer.Methods eng., vol. 26, pp. 589-615.

Simo, J. C.; L. Vu-Quoc, 1991, "A Geometrically-exact Rod Model Incorporating Ar andTorsion-warping Deformation", Internat. J. Solids Structures, vol. 27 (3), pp. 371-393.

Teh, L. H., Murray J. C., 1998, "Co-rotational and Lagragian Formulations for Elastic Three-dimensional Beam Finite Elements", Journal of Constructional Steel Research, vol. 48,pp. 123-144.

Timoshenko, S. P.; Goodier, J. N., 1968, "Teoria de la Elasticidad", Edic. Urmo.

Wells, G. N. ; Sluys, L. J., 2001, "Analysis of Slip in Three-dimensional Solids", Comp. Meth-ods Appl. Mech. Engrg., vol. 190, pp. 3591-3606.

Wood, R.D.; Zienkiewicz, O.C., 1977, "Geometrically Nonlinear Finite Element Analysis ofBeams, Frames, Arches and Axisymetric Shells", Computer.&Structures, vol. 7, pp. 725-735.