ISSN 1984-8218
OTIMIZANDO ELEMENTOS ESTRUTURAIS PLANOS
ATRAVES DE METODO ADAPTATIVO
Eliane Regina Flôres Oliveira
Universidade Federal de Uberlândia – UFU, Instituto de Física – INFIS
38400-000 Uberlândia, MG
e-mail: [email protected]
Elaine Gomes Assis
Universidade Federal de Uberlândia – UFU, Faculdade de Engenharia Mecânica – FEMEC
38400-000 Uberlândia, MG
e-mail: [email protected]
Resumo. Um processo adaptativo usando gradientes ou densidades da energia de deformação
pode levar a um bom modelo de elementos finitos, identificando-se de forma clara quais as
regiões requerem ou não refinamento de malha. O modelo deve conter uma estratégia de
adaptação para evitar a formação de elementos distorcidos, principalmente quando se usa o
método adaptativo de relocação nodal (Método r), que possui essa característica que lhe é
desfavorável. De maneira heurística é delineado o acoplamento do conceito de malha e forma
otimizadas, passando por uma remalhagem, que tem a função de atenuar no processo iterativo
o problema da distorção no uso do método r. Exemplos de aplicação são utilizados para
estudar a influência, do modelo proposto, na otimização de forma e de se usar o módulo
intermediário de remalhagem, o qual equilibra o desbalanceamento das áreas dos elementos na
discretização do modelo.
Introdução
A otimização de forma tem atraído grande atenção da comunidade científica e muitas
técnicas têm sido desenvolvidas e usadas com sucesso em análises e projetos de engenharia.
Geralmente, estas técnicas consistem em variar alguns contornos do modelo a ser projetado a
fim de melhorar seu comportamento mecânico, como por exemplo, reduzir as altas
concentrações de tensões que normalmente acontecem em locais de cantos e furos ou próximos
a eles. Este processo normalmente é feito impondo restrições e usando um método de
otimização selecionado, onde aspectos tais como definições geométricas, geração de malhas,
análise e processamento de resultados são normalmente envolvidos. Além disso, outros
elementos fazem um papel decisivo no processo de otimização, tal como a análise de
sensibilidade e o algoritmo de otimização utilizado.
O primeiro passo é definir os modelos físicos e geométrico onde no modelo geométrico as
variáveis de projeto são facilmente impostas e permitem uma explícita integração com outras
ferramentas de projeto, tal como os sistemas CAD ou CAM. O modelo físico é usado para obter
a resposta da parte estrutural, submetida a ações externas com uma análise numérica de
sensibilidade das variáveis predominantes e finalmente, um algoritmo de otimização apropriado
deve ser usado para resolver o problema de um modo eficaz e seguro.
Como o progresso da otimização muitas vezes depende do desenvolvimento de um bom
modelo de elementos finitos, este trabalho utiliza uma técnica para estimar erros numéricos na
discretização de problemas bidimensionais, utilizando-se o elemento finito triangular
hierárquico com função de forma quadrática. Verifica-se a existência da homogeneidade do
erro ao longo do domínio nos elementos; caso esta não ocorra, os nós livres pré estabelecidos,
serão realocados em função da diferença de gradiente do erro, reiniciando-se uma nova análise
do modelo físico otimizando a malha automaticamente. São avaliados dois métodos adaptativos
usando realocação nodal (método r): método heurístico e método geométrico. Com o objetivo de
atenuar o problema da distorção que ocorre no uso do processo iterativo do método r, este
trabalho desenvolve uma técnica de remalhagem automática de elementos finitos, usando o
elemento triangular hierárquico de grau 2, através do balanceamento da área dos elementos da
malha obtendo-se um método geométrico de relocação nodal. A otimização de forma é baseada
49
ISSN 1984-8218
no método de otimização de malha, Clapis [1999], sendo que a formulação da otimização de
forma leva em conta a homogeneização da densidade de energia de deformação por distorção
dos elementos finitos. O critério de otimalidade é sempre o da máxima densidade de energia de
deformação por distorção. Em linguagem Fortran foi implementado o acoplamento do
otimizador de forma com o remalhador e com o otimizador de malha de maneira iterativa e
automática. A avaliação da potencialidade dos métodos é feita por meio de uma estrutura
tradicional citada na literatura, como apresentado por Oliveira [2003].
Otimização da malha de elementos finitos.
Uma malha de elementos finitos é considerada ótima, geralmente, quando a norma do erro é a
mesma para todos os elementos, sendo assim como o método (FEM) é baseado na minimização
da energia total, torna-se natural estimar-se o erro de discretização em um modelo baseando-se
também na norma desta energia.
O modelo físico utilizado neste trabalho é a equação de elasticidade linear. A aproximação pelo
método de elementos finitos consiste em minimizar a energia potencial total, dada pela equação
(1)

1 eT 
T
 B D
u
e 1 2

n
f


T
B d u  f u
e
(1)
,

em relação as variáveis nodais e as coordenadas nodais. Na equação (1) n é o número total de
elementos da malha necessários para discretizar a estrutura, D é a matriz de constantes elásticas,
B é a matriz de transformação de deslocamentos em deformações, ue é o vetor de deslocamentos
associados ao elemento e, u é o vetor global de deslocamentos, f é o vetor carga aplicada e  é o
volume ocupado pelo elemento.
Método 1 – Método Heurístico
McNeic e Marcal (1971) consideraram o problema de distribuir os nós de uma malha de forma
melhorada para obter melhores resultados. A técnica utilizada consiste em minimizar a energia
potencial total, dada pela equação (1).em relação as variáveis nodais e as coordenadas nodais.
Se c é o vetor de coordenadas nodais, a condição para que a energia potencial seja mínima é que
as derivadas da equação (1) em relação a c e u se igualem a zero. Considerando a decomposição
do vetor u segundo Oliveira (2003), chegamos à seguinte equação
K
n 
K nh  u n  f n 

K u  f , K     BT D B d ou seja 
    
K
K
hn
hh

 u h  f h 
e 1 

O erro por elemento pode ser estimado, em termos da energia de deformação, correspondente às
variáveis hierárquicas como:
e e  u h  K hh  u h , u h   1 1 2  2 3  3 
2
T
(3)
onde: {, } representa os novos graus de liberdade hierárquicos introduzidos, conforme mostra
a Fig. 1 e Khh a matriz de rigidez dos parâmetros hierárquicos.
3
3
3
2
2
1
y
1
1
(2)
2
x
Figura 1 - Parâmetros hierárquicos do elemento triangular
50
ISSN 1984-8218
Considerando e
2
como uma pressão negativa agindo num elemento finito do modelo
e
geométrico da estrutura, Clapis (1999), o vetor das forças nodais equivalentes, De para um
elemento triangular linear é dado por,
De 
e
2
e
y32
 x32
 y31
x31
y 21  x21
T
(4)
Impondo-se a condição geométrica de que a forma do contorno não deve ser alterada e
levando-se em conta a superposição dos vetores De para todos os elementos, obtém-se um
vetor de carga nodal equivalente global, DG .
Tem-se, portanto, novas condições de contorno e consequentemente uma nova equação de
equilíbrio global. A direção d 0 do movimento dos nós com uma ou duas direções livres para
se moverem pode ser obtida de:
K 0 d 0  DG
(5)
sendo:
K 0
é a matriz de rigidez obtida do modelo geométrico, d 0 é a “direção de busca” onde se
procura homogeneizar a distribuição da densidade do erro de discretização e DG é o vetor
carga nodal equivalente global.
Através de uma busca unidimensional, utilizando-se o método da secante, calcula-se o fator
de escala  * , o qual minimizará a norma do vetor desbalanceamento, ou seja:
MIN

Dx  d 0 
(6)
onde {x} é o vetor de coordenadas dos nós móveis.
Devido a relocação nodal, a densidade do erro de discretização também se altera e o
processo de busca é realizado de modo iterativo com a análise de novos modelos físicos
discretizados. O programa executa automaticamente o balanceamento dos erros entre os
elementos, usando o método do gradiente conjugado para gerar novas direções de busca no caso
de não haver convergência numa determinada direção. A condição de ótimo para o problema é
que o erro de discretização seja uniformemente distribuído ao longo de todo o domínio.
Método 2 – Método Geométrico
Considere um nó e os elementos finitos circunvizinhos ao mesmo, Figura 2. Sendo o
centróide de cada elemento o ponto onde atuam forças atrativas, com origem no nó “comum”. A
intensidade de cada uma é o indicador de erro do elemento. A resultante destas forças nos
fornece uma nova posição do nó, caso exista desequilíbrio entre os indicadores dos erros de
discretização na região considerada, Cheng (1993).
y
x
Figura 2 – Esquema de reposicionamento nodal
51
ISSN 1984-8218
A soma das contribuições de todos os elementos circunvizinhos, para a coordenada x, vale:
xnm  xnm 1 
*
m 1
 e ( x*  xn )
*
e
(7)
onde m indica o número da iteração, xn e x* são, respectivamente, a posição do nó considerado e
dos centróides dos elementos circunvizinhos, e e* é o erro por elemento finito. Este pode ser
estimado, em termos da norma da energia de deformação, Eq. (1), onde o critério de parada é a
homogeneização do mesmo ao longo do domínio.
Analogamente, para a coordenada y, tem-se equação para y n como a apresentada na
equação (7).
Nota-se que apenas os nós dos vértices do elemento têm a permissão de moverem sendo que
os nós do meio do lado, se houver, não podem participar do processo de reposicionamento. Suas
coordenadas são localizadas no meio entre os nós da extremidade. Este processo tem natureza
iterativa onde os nós se movem em direção aos elementos de indicador de erro maior até não se
ter diferença nos indicadores de erro dos elementos ou a relocação dos nós tornar-se
suficientemente pequeno.
Algoritmo de otimização de forma
Propomos um método de otimização de forma com base no método de otimização de malha,
Clapis (1999). Enquant o Método 1 leva em consideração o desbalanceamento do erro de
discretização do modelo, a formulação da otimização de forma toma por base o
desbalanceamento da densidade de energia de deformação por distorção dos elementos finitos.
O critério de otimalidade utilizado é o da máxima densidade de energia de deformação por
distorção, de forma que os elementos com menor energia migrem para as áreas onde as
densidades de energia são maiores.
Seja UDe uma pressão negativa agindo num elemento finito do modelo geométrico da
estrutura. O vetor das forças nodais equivalentes, De , para um elemento triangular é dado pela
equação (4). Da superposição dos vetores De para todos os elementos obtém-se um vetor de
carga nodal equivalente global, DG . A direção d0 do movimento dos nós pode ser obtido
utilizando-se a equação (5).
Portanto, d0 , passa a ser uma direção de busca procurando homogeneizar a distribuição
da densidade de energia de distorção, conforme a equação (6), onde {x} é o vetor das
coordenadas nodais dos nós livres para realocarem-se nas duas ou em uma das direções.
Algoritmo de remalhagem
O código numérico de remalhagem permite reestruturar malhas de elementos finitos em
problemas bidimensionais, utilizando elementos triangulares hierárquicos com funções de forma
de grau 2. Automaticamente, o vetor desbalanceamento, na formulação geométrica de relocação
nodal, é função do diferencial entre as áreas dos elementos funcionando como cargas nodais
passando pelos centróides dos elementos circunvizinhos ao nó considerado, Cheug (1993).
Verifica-se se existe uma homogeneidade de áreas dos elementos da malha na região
considerada, caso esta não ocorra, os nós livres, pré-estabelecidos, serão deslocados em função
da resultante das cargas e a nova coordenada x, vale:
m
m 1

xn  xn

 A* x*  x n
 A*

m 1
(8)
52
ISSN 1984-8218
onde m indica o número da iteração, xn e x* são, respectivamente, a posição do nó considerado e
dos centróides dos elementos circunvizinhos, e A* é a área de cada elemento
finito.Analogamente para a coordenada y, tem-se equação análoga a equação (8).
Este processo tem natureza iterativa onde os nós se movem até diminuir a distorção entre as


m 1 tornar-se suficientemente
áreas dos elementos ou a relocação dos nós até o valor x m
n  xn
pequena.
Exemplo de aplicação - Coluna de auto estrada
O projeto de otimização conjunta de forma e malha de uma das colunas da pista de uma
determinada auto estrada é considerado neste exemplo. Oda (1977), Kikuchi (1986), Rossow
(1976) e Clápis (1999) utilizaram deste mesmo tipo de estrutura para avaliarem suas propostas
na busca da forma ideal ou melhorada, do elemento estrutural considerado.
A Figura 3 mostra as características de projeto, considerando-se para efeito de somente
avaliar o método proposto, que a estrutura é de aço em vez de concreto armado o que
usualmente é utilizado para confeccionar um componente estrutural deste porte. Os dados
adotados são: módulo de elasticidade 205000 MPa, coeficiente de Poisson 0.30, espessura
constante da estrutura de 1000 mm e limite de escoamento à tração igual a 7 x 107 N/m2.
Convém salientar que a escolha do material baseia-se no critério da energia de distorção
(von Mises) levando em conta o escoamento, o que não é apropriado para o caso quando se usa
um material frágil como o concreto, onde o mais apropriado seria utilizar um critério de ruptura.
Devido à simetria em relação ao eixo vertical, utiliza-se somente a metade do modelo físico
na otimização conjunta de forma e malha,. As Figuras 4 (b).e 4 (c) mostram, respectivamente, o
modelo geométrico utilizado na otimização de forma com 156 elementos triangulares lineares e
o modelo geométrico utilizado na remalhagem e na otimização da malha com 156 elementos
triangulares hierárquicos em expansão quadrática na função de forma.
16 m
1000 N/m
1m
4m
8m
6m
Figura 3 – Seção transversal da coluna de uma auto-estrada
(b)
(c)
(a)
Figura 4 – Discretizações da coluna: (a) modelo físico, (b) modelo geométrico da forma e (c)
modelo geométrico utilizado na remalhagem e na otimização da malha
53
ISSN 1984-8218
Após 4 iterações, considerando-se o estado plano de tensões, chegou-se à forma otimizada
utilizando-se somente o programa otimizador de forma e o critério de parada pré-estabelecido,
Figura 5. Nota-se, neste caso, que o volume da seção transversal passou de um volume inicial de
48,5m3 para um volume final de 38,9m3 (redução de 24,68%), bem como se observa uma
melhor concordância do contorno na região onde existe uma alta concentração de tensões. Notase também nesta região uma grande distorção dos elementos finitos ocasionado pela relocação
nodal no balanceamento da energia de distorção.
Fazendo-se a mudança de forma com a intercalação dos módulos de remalhagem e de
otimização de malha no processamento numérico e estimando-se para a análise em questão o
número de iterações para cada módulo, tem-se para o método heurístico: 1 iteração para a
forma; 6 iterações para a remalhagem; 10 iterações para a otimização de malha. Para os erros de
discretização após a otimização da malha. O processo é refeito 8 vezes até que o critério de
parada seja atingido.
Através do método geométrico tem-se: 1 iteração para a forma; 5 iterações para a
remalhagem, 7 iterações para a otimização de malha. Para os erros de discretização após a
otimização da malha o processo é refeito 9 vezes utilizando-se o critério de parada anterior.
Observa-se que no processo de otimização de forma intercalando-se a otimização de malha,
a redução do volume da seção transversal é muito maior e a distorção dos elementos finitos é
bem menor do que no caso onde foi utilizado somente o otimizador de forma, obtendo-se uma
discretização final de qualidades física e geométrica melhorada que se deve a influência da
introdução do módulo de remalhagem. Usando-se a otimização conjunta forma e malha, pelo
método heurístico, o volume passou de um valor inicial de 48,5 m3 para um volume final de
32,47 m3 (redução de 49,37%). Pelo método geométrico, o volume final foi para 33,06m3
(redução de 46,7%), enquanto que se utilizando o módulo somente do otimizador de forma o
volume final é de 38,9 m3 (redução de 24,68%).
Figura 5 – Tensões equivalentes após otimização de forma – 4 iterações
Figura 6 – (a) Otimização da forma / tensões equivalentes – método heurístico (1 iteração, 8ª vez), (b) Remalhagem –
método heurístico (7 iterações, 8ª vez), (c) Otimização da malha – método heurístico (10 iterações, 8ª vez).
(a
))
(b
)
(c
)
Fazendo-se a mudança de forma com a intercalação dos módulos de remalhagem e de
otimização de malha no processamento numérico e estimando-se para a análise em questão o
número de iterações para cada módulo, tem-se para o método heurístico: 1 iteração para a
54
ISSN 1984-8218
forma, Figura 6: (a); 7 iterações para a remalhagem; 6(b); 10 iterações para a otimização de
malha, Figura 6(c).
A sensibilidade do uso do programa quanto ao número de iterações a ser utilizada em
cada um dos três módulos (forma, remalha e malha), depende da intuição e experiência do
projetista. No módulo de otimização de forma o número de iterações pode ser aumentado até
que não ocorra grande distorção dos elementos finitos. No módulo de remalhagem o critério de
parada utilizado foi em função do balanceamento das áreas dos elementos finitos e no módulo
de otimização de malha o número de iterações foi fixado em função da obtenção de erros de
discretização dos elementos finitos mais homogêneos sendo que na maioria dos casos fez-se
necessário alguns testes iniciais para se obter a melhor solução do problema.
Conclusões
Desenvolveu-se um programa de elementos finitos para otimizar a forma de um elemento
estrutural baseado na homogeneização da energia de distorção entre os elementos, adaptado
iterativamente com um remalhador e com um otimizador de malha, a cada iteração ou a partir
de certo número de iterações previamente fixado.
Conclui-se que o otimizador de forma deve ser acoplado ao otimizador de malha desde que
o modulo de remalhagem seja introduzido para corrigir a distorção dos elementos após a
otimização da forma e a convergência do sistema seja alcançada.
Na formulação numérica do modulo otimizador de forma optamos por usar o elemento
triangular linear, pois usando-se o elemento triangular hierárquico o programa se torna lento e
de difícil convergência.
A utilização do programa em estruturas onde havia alta concentração de tensões nos cantos,
a técnica de intercalar o otimizador de forma com o remalhador e com o otimizador de malha
levou a uma forma ótima na qual houve um decréscimo no valor das tensões equivalentes e a
uma melhor distribuição de tensões ao longo do contorno móvel quando comparado com o uso
somente do programa otimizador de forma.
Pode-se observar ainda que, na medida em que o procedimento aqui apresentado tende a
suavizar a distribuição de tensões no contorno, o mesmo pode resultar em economia de material,
pois o volume final é menor que o inicial, o que abre a perspectiva de utilizá-lo no projeto ótimo
de componentes estruturais.
O programa utiliza o método de elementos finitos, que é uma eficiente técnica para obter
soluções aproximadas para problemas onde a solução exata é difícil de ser obtida, a
discretização inicial da malha, o tipo e o número de elementos influenciam nos resultados finais
da forma ótima.
A partir do exemplo mostrado, pode-se concluir que o algoritmo desenvolvido apresenta
resultados satisfatórios, quando comparado com a bibliografia encontrada.
Referências
[1] Clápis, A. P., 1999, “Um método heurístico de otimização de forma de componentes
estruturais no estado plano de elasticidade linear”, Tese de Doutorado, Universidade Estadual de
Campinas, Campinas, SP, Brasil.
[2] Cheng, Jung-Ho, 1993, “Adaptative grid optimization for structural analysis-geometry-based
approach”, Comp. Meth. in Applied Mech. and Eng., 107,1-22
[3] Kikuchi N., Chug, K. Y., Torigaki, J. E., Taylor, J. E., 1986, “Adaptative finite elements
methods for shape optimization of linearly elastic structures”, Comp. Meth. in Applied Mech.
and Eng., 67-89.
[4] McNeic B.M. and Marcal P.V., 1971, “Optimization of Finite Element Grids Based on
Maximum Potential Energy”, Technical Report no. 7, Brown University, Providence.
[5] Oda, J., Yamazaki, K., 1977, “On a technique to obtain an optimum stength shape by the
finite element method”, Bulletin of JSME, vol. 20, pp. 1524-1532.
[6] Oliveira, E. R. F., 2003, “Otimização de forma de elementos estruturais planos”, Tese de
Doutorado, Universidade Federal de Uberlândia, MG, Brasil
[7] Rossow, M. P., Taylor, J. E., 1976, “An optimal structural design algorithm using opmality
criteria”, Society Engeneering Science, 13th Annual Meeting, Hampton, VA.
55