Análise da resposta constitutiva de materiais heterogêneos pelo Método dos Elementos de Contorno, considerando-se diferentes microestruturas para o EVR

Analysis of the constitutive response of heterogeneousmaterials by the Boundary Element Method,considering different RVE microstructures

RESUMO

São feitas análises numéricas para modelar a resposta constitutiva de materiais heterogêneos utilizando-se uma formulação do Método dos Elementos de Contorno desenvolvida no contexto de uma teoria multi-escala baseada em conceito de EVR. Nesse caso, a microestrutura do material é definida como EVR (Elemento deVolume Representativo), sendo modelada por uma placa em sub-regiões, onde cada sub-região representa a matriz ou uma inclusão, podendo ainda serem definidos vazios no interior da matriz, a fim de modelar a micro-estrutura de materiais porosos. O fenômeno dissipativo no EVR é levado em conta através da considera-ção de um campo de esforços iniciais que representam as forças dissipativas. Após impor ao EVR um vetor de deformações constante ao seu contorno, obtém-se sua resposta constitutiva homogeneizada. Mas para isso, deve-se antes resolver o problema de equilíbrio do EVR, que é definido em termos de flutuação dos deslocamentos. Nos exemplos numéricos consideram-se diferentes micro-estruturas para o EVR, a fim de mostrar como a resposta constitutiva se modifica em função da fração volumétrica adotada para as inclusões ou vazios. Além disso, para um EVR onde são definidas cinco inclusões de forma aleatória, analisa-se sua resposta constitutiva devido a diferentes vetores de deformações impostos ao seu contorno. Os valores homogeneizados das tensões e do tensor constitutivo são comparados com a formulação desenvolvida com o Método dos Elementos Finitos para validar a formulação do Método dos Elementos de Contorno.

Palavras-chave
elementos de contorno; problema bi-dimensional; técnicas de homogeneização; EVR; modelagem multi-escala

ABSTRACT

Numerical analyses to model the constitutive response of heterogeneous materials by a Boundary Element Formulation, developed in the context of a RVE-based multi-scale theory, are performed. In this case, the material microstructure is denoted as RVE (Representative Volume Element) and it is modelled by a zoned plate, where each sub-region represents the matrix or an inclusion. Besides, voids can also be defined inside the matrix to model the microstructure of porous materials. The dissipative phenomenon in the RVE is taken into account by considering an initial forces field which represents the dissipative forces. After imposing a constant strain vector to the RVE boundary, its homogenized constitutive response must be obtained. For that, the RVE equilibrium problem must be solved first, which is defined in terms of displacements fluctuations. In the numerical examples different microstructures are adopted to show how the RVE constitutive response changes when different volume fractions of inclusions or voids are considered. Besides, the constitutive response of a RVE where five inclusions are disposed aleatory is evaluated when different strains vectors are imposed to its boundary. The homogenized values for the stress vector and the constitutive tensor are com pared to the formulation developed by the Finite Element Method in order to validate the Boundary ElementFormulation.

Keywords
Boundary elements; two-dimensional problem; homogenization techniques; RVE; multi-scale modeling

1. INTRODUÇÃO

Esse trabalho trata da modelagem da resposta constitutiva de materiais heterogêneos através do MEC (Método dos Elementos de Contorno), adotando-se uma teoria multi-escala baseada no conceito de EVR (Elemento de Volume Representativo). A formulação aqui desenvolvida foi inicialmente proposta em [¹], sendo que no presente trabalho novas simulações numéricas são discutidas a fim auxiliar na validação dessa formulação.

Em análises não-lineares convencionais, o comportamento do material é governado por um modelo constitutivo, onde dada uma deformação, o mesmo fornece as tensões e o tensor constitutivo. Nesse caso, para representar de forma mais precisa o comportamento complexo de certos materiais, é necessário um modelo constitutivo muito sofisticado, onde se definem muitos parâmetros, sendo alguns de difícil identificação. No presente trabalho, a resposta constitutiva do material heterogêneo é dada pela resposta homogeneizada do EVR, que somente pode ser determinada após resolver seu problema de equilíbrio, que é definido em termos de flutuação dos deslocamentos. Note que surge um campo de flutuações de deslocamentos no EVR apenas quando sua microestrutura é heterogênea, ou seja, para EVRs uniformes sua resposta homogeneizada é a mesma daquela quando a deformação é considerada diretamente no modelo constitutivo, sem necessidade de definição de EVR. Nos trabalhos [²] a [⁹] tem-se outros diferentes tipos de modelagem do comportamento constitutivo do material através de uma abordagem multi-escala. Assim, no domínio do EVR , que representa a micro-estrutura do material estudado, podem-se definir vazios ou inclusões (elásticas ou elasto-plásticas) dentro de sua matriz, sendo a formação e propagação de micro-fissuras ou o processo de plastificação no EVR, estudados após impor ao contorno do mesmo um vetor de deformações, o que é feito impondo-se deslocamentos lineares aos nós do contorno do EVR. Pode-se considerar diferentes propriedades elásticas e modelos constitutivos nas fases do EVR, sendo cada fase a representação de um material diferente de sua microestrutura. Assim, após impor os deslocamentos lineares como condição de contorno do problema, deve-se resolver o problema de equilíbrio do EVR, obtendo-se o campo de tensões e os tensores constitutivos nas células do EVR que satisfazem sua equação de equilíbrio. Em seguida, através de princípios de homogeneidade e conceito de média volumétrica, obtém-se a resposta constitutiva homogeneizada do EVR.

Apesar de nesse trabalho não ser considerada uma análise multi-escala de uma estrutura, é importante definir a presente formulação considerando-se uma abordagem multi-escala, a fim de que seja possível definir algumas variáveis e também para facilitar no futuro a realização de análises multi-escala de uma estrutura, que deverá ser obtida através do acoplamento da presente formulação com um modelo numérico que representa o problema do macro-contínuo a ser estudado. Observe que a formulação aqui apresentada para modelar a micro-estrutura do material, segue os mesmos princípios daquela proposta em [¹⁰]-[¹³], onde o problema de equilíbrio do EVR é modelado pelo MEF (Método dos Elementos Finitos). Note que o simples fato de se adotar o MEC ao invés do MEF na modelagem da microestrutura requer uma solução numérica muito diferente para se obter a resposta elástica necessária para se iniciar o processo iterativo no EVR e também para calcular as forças no contorno depois que o problema de equilíbrio do EVR é solucionado, a partir das quais se obtém o tensor de tensões homogeneizado. Além disso, no modelo com o MEC não há necessidade de aproximar as forças elásticas no domínio do EVR, apenas as forças plásticas devem ser aproximadas. Esses fatores levam a obter uma resposta homogeneizada apenas similar com os dois métodos numéricos, mas não idêntica.

De acordo com a teoria multi-escala proposta em [¹⁰]-[¹³], o campo de deslocamento no EVR é dividido em duas parcelas: uma referente à deformação imposta ao EVR e outra denominada de flutuação dos deslocamentos, que surge quando a microestrutura do material é heterogênea, isto é, composta por mais de uma fase (matriz e inclusões) ou quando se definem vazios no interior da matriz. Assim, a solução do EVR é obtida quando se encontra o campo de flutuações de deslocamentos que satisfaz sua equação de equilíbrio, que por sua vez é satisfeita quando se encontra um campo de tensões no EVR que esteja auto-equilibrado. Porém, para que o problema de equilíbrio do EVR seja resolvido é necessário impor condições de contorno ao mesmo em termos de flutuação dos deslocamentos, sendo que cada condição de contorno leva a uma resposta numérica diferente, ou seja, a um modelo multi-escala diferente.

O surgimento do MEC e a sua utilização como uma alternativa para a obtenção de soluções numéricas, representa um avanço significativo que ocorreu nessa área do conhecimento nos últimos anos. Em muitos problemas, comprovadamente, esse método é uma alternativa mais precisa e que permite a obtenção de respostas mais confiáveis quando comparadas com as dos métodos usuais (veja [¹⁴]-[¹⁶]). Nos trabalhos [¹⁷]-[²⁰] desenvolvem-se modelos multi-escala para análise de estruturas, onde se acopla a formulação do MEC para o macro-contínuo proposta em [²¹] com a formulação do MEF para o modelar o EVR proposta em [¹⁰], [¹²] e [¹³]. Em [²²] uma formulação do MEC para análise de flexão composta de placas enrijecidas é proposta onde o problema de flexão simples é acoplado com o problema de chapa, podendo as vigas e lajes ser definidas com diferentes propriedades elásticas. Em [²³] a formulação elástica para o problema bidimensional proposta em [²²] é modificada a fim de obter uma formulação do MEC para modelar a resposta constitutiva elástica de materiais heterogêneos. Em seguida, em [¹] a formulação proposta em [²³] é novamente modificada a fim de considerar fenômenos dissipativos no EVR e assim possibilitar modelar a resposta constitutiva de materiais heterogêneos levando-se em conta deformações residuais.

No presente trabalho é apresentada a mesma formulação proposta em [¹], porém diferentes exemplos numéricos são discutidos, a fim de auxiliar na validação dessa formulação. Assim, nos exemplos numéricos inicialmente estuda-se como a variação do volume de fração de uma inclusão centrada no EVR influencia seu comportamento constitutivo homogeneizado. Em seguida, esse mesmo estudo é feito considerando-se um vazio no lugar da inclusão. Então, adotando-se um EVR com cinco inclusões aleatoriamente posicionadas, estuda-se como a resposta constitutiva homogeneizada do EVR modifica em função da condição de contorno adotada em termos de flutuação dos deslocamentos, ou seja, estuda-se como a resposta constitutiva do material estudado é modificada em função do modelo multi-escala adotado. A fim de validar o modelo aqui desenvolvido, a resposta constitutiva do EVR, ou seja, os valores homogeneizados do vetor de tensão e do tensor constitutivo são comparados com o modelo do MEF proposto em [¹⁰], [¹²] e [¹³].

MATERIAIS E MÉTODOS

A seguir, inicialmente serão apresentadas as equações integrais e algébricas do MEC para resolver o problema da chapa em sub-regiões composta de diferentes materiais. Então, serão discutidos conceitos relacionados à modelagem multi-escala como: definição do campo de deslocamentos no EVR, equação de equilíbrio do EVR, solução do seu problema de equilíbrio e obtenção da sua resposta constitutiva homogeneizada. Esses tópicos referentes à modelagem multi-escala seguem os princípios da teoria multi-escala proposta em [¹⁰], [¹²] e [¹³], que serão adotados nesse trabalho. Maiores detalhes sobre a implementação da formulação desenvolvida podem ser encontrados em [¹].

2.1 Equações Integrais e Algébricas do MEC

O EVR será modelado por uma chapa dividida em sub-regiões, sendo cada sub-região a representação de uma inclusão ou da matriz, como se pode ver na Figura 1, onde se tem um EVR com 1 inclusão (domínio Ω₂) e 1 vazio (domínio Ω₃) definidos dentro da matriz (domínio Ω₁), totalizando 2 sub-regiões que formam o domínio sólido. Além disso, na Figura 1 Γjk representa a interface entre as sub-regiões adjacentes Ω_j and Ω_k, sendo o sistema de coordenadas cartesianas (X₁ e X₂) definidos na superfície da chapa. Note que no modelo proposto, o contorno externo sempre coincidirá com o contorno da matriz Γ₁, ou seja, as inclusões ou vazios serão sempre dispostas no interior da matriz.

Figura 1
EVR com uma inclusão e um vazio definidos na matriz

A partir do Teorema de Reciprocidade de Betti pode-se obter a equação integral que representa o problema bi-dimensional na chapa. Como descrito em maiores detalhes em [²³], para uma chapa composta de Ns sub-regiões no domínio sólido, considerando-se o caso de estado plano de tensões, chega-se a:

$${\int }_{\Omega }{\epsilon }_{kij}^{*}\stackrel{.}{N_{ij}^e}d\Omega =\underset{m=1}{\overset{N_S}{\Sigma }}\left[\frac{{{\displaystyle \stackrel{=}{E}}}_m}{{\displaystyle \stackrel{=}{E}}}\frac{v_m}{v}{\int }_{\Omega }{\epsilon }_{kij}^{*}\stackrel{.}{N_{ij}^e}d\Omega +{\stackrel{=}{E}}_m\left(1-\frac{v_m}{v}\right){\int }_{\Omega }{\epsilon }_{kij}^{*}\stackrel{.}{N_{ij}^e}d{\Omega }_m\right]\mathrm{k},\mathrm{i},\mathrm{j}=1,2$$ (1)

onde k é a direção da carga fundamental, os termos em * são referentes ao problema fundamental, v é o coeficiente de Poisson, E é o módulo de elasticidade, t a espessura,$\stackrel{.}{N_{ij}^e}$ representa as forças de membrana elásticas;${\epsilon }_{kij}^{*},\overline{E},ve{N}_{kij}^{*}$são valores referentes à sub-região onde está o ponto de colocação; para uma sub-região Ω_m, tem-se:${\overline{E}}_m=\frac{{\overline{E}}_m}{{\left(1-{v_m}^2\right)}^{.}}$.

Porém, quando se considera fenômenos dissipativos, as forças de membrana elásticas (obtidas a partir das deformações totais) são definidas como a soma das forças N_ij calculadas a partir das tensões que satisfazem o modelo constitutivo e das forças plásticas N_ij^p (ver [²¹] e [²⁴]), ou seja:

$${\dot{N}}_{ij}^e=\frac{\overline{E}}{\left(1-v^2\right)}\left[v{\dot{\epsilon }}_{kk}{\delta }_{ij}+\left(1-v\right){\dot{\epsilon }}_{ij}\right]={\dot{N}}_{ij}+{\dot{N}}_{ij}^p$$ (2)

Substituindo-se a Eq. (2) na Eq. (1), esta se torna:

$${\int }_{\Omega }{\dot{\epsilon }}_{kij}^{*}\left({\dot{N}}_{ij}+{\dot{N}}_{ij}^p\right)d\Omega =\underset{m=1}{\overset{N_S}{\Sigma }}\left[\frac{{{\displaystyle \stackrel{=}{E}}}_m}{{\displaystyle \stackrel{=}{E}}}\frac{v_m}{v}{\int }_{\Omega }{\dot{\epsilon }}_{kij}^{*}\left({\dot{N}}_{ij}+{\dot{N}}_{ij}^p\right)d\Omega +{\stackrel{=}{E}}_m\left(1-\frac{v_m}{v}\right){\int }_{\Omega }{\dot{\epsilon }}_{ij}{\epsilon }_{kij}^{*}d{\Omega }_m\right]$$ (3)

Deve-se agora fazer a integração por partes da equação (3) a fim de se obter a equação integral dos deslocamentos no plano da chapa escrita em função de integrais ao longo do contorno externo da placa e interfaces. Para um caso geral em que se tenham N_vazios vazios e N_inc inclusões, a equação integral do deslocamento de um ponto do domínio se torna:

$$\begin{gathered} C_{k1}{\dot{u}}_1\left(s\right)+C_{k2}{\dot{u}}_2\left(s\right)=\frac{{{\displaystyle \stackrel{=}{E}}}_1{v}_1}{{\stackrel{=}{E}}_v}{\int }_{{\Gamma }_1}\left({\dot{u}}_1\left(P\right)p_{k1}^{*}\left(s,P\right)+{\dot{u}}_2\left(P\right)p_{k2}^{*}\left(s,P\right)\right)d\Gamma + \\ -\underset{m=1}{\overset{N_{inc}}{\Sigma }}\left(\frac{{\stackrel{=}{E}}_m{v}_m-{\stackrel{=}{E}}_1{v}_1}{{\stackrel{=}{E}}_v}\right){\int }_{{\Gamma }_{m1}}\left({\dot{u}}_1\left(P\right)p_{k1}^{*}\left(s,P\right)+{\dot{u}}_2\left(P\right)p_{k2}^{*}\left(s,P\right)\right)d{\Gamma }_{m1} \\ -\underset{m=1}{\overset{N_{vazios}}{\Sigma }}\frac{{\displaystyle {\stackrel{=}{E}}_1{v}_1}}{{\stackrel{=}{E}}_v}{\int }_{{\Gamma }_{1m}}\left({\dot{u}}_1\left(P\right)p_{k1}^{*}\left(s,P\right)+{\dot{u}}_2\left(P\right)p_{k2}^{*}\left(s,P\right)\right)d{\Gamma }_{1m}+ \\ {\int }_{{\Gamma }_1}\left(u_{k1}^{*}\left(P\right){\dot{p}}_1\left(s,P\right)+u_{k2}^{*}\left(P\right){\dot{p}}_2\left(s,P\right)\right)d\Gamma + \\ -{\stackrel{=}{E}}_1\left(1-\frac{v_1}{v}\right){\int }_{{\Gamma }_1}\left({\dot{u}}_1\left(P\right){\epsilon }_{k1}^{*}\left(s,P\right)+{\dot{u}}_2\left(P\right){\epsilon }_{k2}^{*}\left(s,P\right)\right)d\Gamma + \\ \underset{m=1}{\overset{N_S}{\Sigma }}{\stackrel{=}{E}}_m\left(1-\frac{v_m}{v}\right){\int }_{{\Omega }_m}{\dot{u}}_1\left(P\right){\epsilon }_{kij,j}^{*}\left(s,p\right)d{\Omega }_m \\ -\underset{m=1}{\overset{N_{inc}}{\Sigma }}\left[{\stackrel{=}{E}}_m\left(1-\frac{v_m}{v}\right)-{\stackrel{=}{E}}_1\left(1-\frac{v_1}{v}\right)\right]{\int }_{{\Gamma }_{m1}}\left[{\dot{u}}_1\left(P\right){\epsilon }_{k1}^{*}\left(s,P\right)+{\dot{u}}_2\left(P\right){\epsilon }_{k2}^{*}\left(s,P\right)\right]d{\Gamma }_{m1} \\ +{\int }_{\Omega }{\epsilon }_{kij}^{*}\left(s,p\right)+{\dot{N}}_{ij}^p\left(p\right)d\Omega k,i,j=1,2 \end{gathered}$$ (4)

onde ${\stackrel{=}{E}}_1$e v₁se referem às propriedades da matriz, k é a direção da carga fundamental, u₁ e u₂ são os deslocamentos no plano da chapa, p1 e p2 são as forças no plano da chapa, os termos livres são definidos em [²³], s é o ponto de colocação e p o ponto campo.

Na equação (4) têm-se duas integrais de domínios a ser resolvidas: uma integral envolvendo os deslocamentos e outra envolvendo os esforços iniciais ou plásticos. Para resolvê-las, o domínio será discretizado através de células triangulares, onde os deslocamentos terão variação linear e os esforços iniciais serão considerados constantes (maiores detalhes sobre a implementação numérica do modelo ver em [¹]). Assim, após a discretização da chapa em elementos de contorno, elementos de interface e células no domínio, nessa formulação ter-se-á quatro variáveis nos nós do contorno (u₁, u₂, p₁, p₂), das quais apenas duas são incógnitas, pois duas são dadas como condição de contorno. Além disso, têm-se duas incógnitas (u₁, u₂) nos nós de interface e nos nós de domínio e as forças plásticas nos nós centrais das células. Portanto, na montagem do sistema de equações serão escritas duas equações de deslocamentos nos nós sobre o contorno externo e nos nós de interface e de domínio. Como no problema de multi-escala são prescritos deslocamentos em todos os nós do contorno, todas as incógnitas no contorno externo são relacionadas às forças de superfície. Portanto, escrevem-se duas equações de deslocamentos em todos os nós do RVE, chegando-se ao seguinte sistema de equações algébricas para um determinado incremento de carga:

$$\left[\begin{array}{cc}{\left[H\right]}_{cc}& {\left[H\right]}_{ci}\\ {\left[H\right]}_{ic}& {\left[H\right]}_{ii}\end{array}\right]\left\{\begin{array}{c}{\left\{\Delta U\right\}}_c\\ {\left\{\Delta U\right\}}_i\end{array}\right\}=\left[\begin{array}{c}{\left[G\right]}_{cc}\\ {\left[G\right]}_{ic}\end{array}\right]{\left\{\Delta P\right\}}_c+\left[E\right]\left\{\Delta N^p\right\}$$ (5)

onde os índices c e i se referem, respectivamente, a nós sobre o contorno e internos (onde nesses estão incluídos também aqueles sobre as interfaces); {U} traz os deslocamentos nodais; {N^p} é o vetor das forças iniciais ou plásticas, de dimensão 3 x nº de células (pois N^p tem variação constante na célula); {P}_c é o vetor das forças nodais no contorno, que nesse caso, serão forças reativas devido aos deslocamentos nodais do contorno do RVE.

Impondo-se as condições de contorno, referentes aos deslocamentos lineares prescritos nos nós do contorno, obtêm-se os incrementos das incógnitas:

$$\left\{\Delta X\right\}=\left\{\Delta L\right\}+\left[R\right]\left\{\Delta N^p\right\}$$ (6)

onde:

$$\left\{\Delta X\right\}=\left\{\begin{array}{c}{\left\{\Delta P\right\}}_c\\ {\left\{\Delta U\right\}}_i\end{array}\right\}$$ (7)

$$\left\{\Delta L\right\}={\left[A\right]}^{-1}\left\{\Delta B\right\}$$ (8)

$$\left\{\Delta B\right\}=\left[\begin{array}{c}-{\left[H\right]}_{cc}\\ -{\left[H\right]}_{ic}\end{array}\right]{\left\{\Delta U\right\}}_c$$ (9)

$$\left[R\right]={\left[A\right]}^{-1}\left[E\right]$$ (10)

$$\left[A\right]=\left[\begin{array}{cc}-{\left[G\right]}_{cc}& {\left[H\right]}_{ci}\\ -{\left[G\right]}_{ic}& {\left[H\right]}_{ii}\end{array}\right]$$ (11)

Para resolver o problema de equilíbrio do EVR, a ser discutido na próxima seção, devem-se ter as forças de membrana elásticas nos nós centrais das células. Essas forças são obtidas a partir das deformações (ver equação 2), as quais, por sua vez, são obtidas a partir das derivadas dos deslocamentos. Assim, após escrever a equação da derivada dos deslocamentos em todos nós centrais das células e aplicar a lei de Hooke, chega-se à seguinte equação algébrica para os esforços de membrana elásticos, após já ter imposto as condições de contorno:

$$\left\{\Delta N^e\right\}=-\left[A\text{'}\right]\left\{\Delta X\right\}+\left\{\Delta B\text{'}\right\}+\left[E\text{'}\right]\left\{\Delta N^P\right\}$$ (12)

onde:

$$\left[A\text{'}\right]=\left[{\left[G\text{'}\right]}_c-{\left[H\text{'}\right]}_i\right]$$ (13)

$$\left\{\Delta B\text{'}\right\}=-{\left[H\text{'}\right]}_c{\left\{\Delta U\right\}}_c$$ (14)

ou substituindo-se a Eq. (7) na Eq. (12), chega-se a:

$$\left\{\Delta N^e\right\}=\left\{\Delta K_n\right\}+\left[S\right]\left\{\Delta N^P\right\}$$ (15)

onde:

$$\left\{\Delta K_N\right\}=\left\{\Delta B\text{'}\right\}+\left[A\text{'}\right]\left\{\Delta L\right\}$$ (16)

$$\left[S\right]=\left[E\text{'}\right]-\left[A\text{'}\right]\left[R\right]$$ (17)

Observe que nas equações (6) e (15) as parcelas ΔL e ΔK_N se referem à solução elástica e as parcelas que envolvem os esforços plásticos apenas serão levadas em conta após a solução do problema de equilíbrio do EVR, a ser discutido na próxima seção.

2.2 Modelagem multi-escala: campo de deslocamento no EVR, problema de equilíbrio do EVR, resposta constitutiva homogeneizada.

Note que apesar de nesse trabalho as análises serem feitas apenas no nível da micro-estrutura, a formulação é desenvolvida no contexto da modelagem multi-escala e por isso, a seguir definem-se variáveis referentes ao macro-contínuo. De acordo com a formulação multi-escala proposta em [¹⁰], [¹²] e [¹³], o campo de deslocamento no EVR é composto pela parcela u^ε obtida a partir da macro deformação ε(x) imposta de forma constante ao contorno do EVR e da parcela ũ denominada flutuação de deslocamentos. Essa segunda parcela surge quando a microestrutura é heterogênea, isto é, quando se definem vazios e/ou inclusões no interior da matriz. Assim, considerando-se um ponto genérico y do EVR, seu campo de deslocamento é definido como:

$$u_{\mu }\left(y\right)=u_{\mu }^{\epsilon }\left(y\right)+{\tilde{u}}_{\mu }\left(y\right)$$ (18)

onde o sub-índice µ indica valores referente à microestrutura.

Analogamente, o campo de deformações no EVR é dado por:

$${\epsilon }_{\mu }\left(y\right)=\epsilon \left(x\right)+{\tilde{\epsilon }}_{\mu }\left(y\right)$$ (19)

onde:

$${\tilde{\epsilon }}_{\mu }=B{\tilde{u}}_{\mu }$$ (20)

sendo x o ponto do macro-contínuo representado pelo EVR, ε(x) é a deformação imposta pelo macro-contínuo de forma constante ao contorno do EVR e B a matriz que relaciona deslocamentos com deformações na célula.

De acordo com a formulação proposta em [¹⁰], [¹²], e [¹³], o EVR encontra-se em equilíbrio quando se obtém um campo de tensões que esteja auto equilibrado, ou seja, sua equação de equilíbrio é definida como:

$$R={\int }_{{\Omega }_{\mu }}{\sigma }_{\mu }\left(y\right)dV=0$$ (21)

ou, após a discretização do EVR em células, a equação de equilíbrio pode ser escrita em termos de flutuação dos deslocamentos como (ver mais detalhes em [¹] e [²³]):

$$R=\underset{e=1}{\overset{N_{cel}}{\Sigma }}{B}_e^T\Delta N^e{A}_e=\underset{e=1}{\overset{N_{cel}}{\Sigma }}{\left[B\right]}_e^T\left[C^{ep\left(e\right)}\right]{\left(\left\{\Delta \epsilon \right\}+{\left[B\right]}_e\left\{\widetilde{\Delta u}\right\}\right)}_e{A}_{e}t=0$$ (22)

onde A_e é a área da célula e, t a espessura da chapa, Δε é o vetor de deformação imposto pelo macro (no contorno essa deformação é constante, mas nos pontos internos ela não será constante se a microestrutura for heterogênea), C^ep é o tensor constitutivo obtido a partir do modelo constitutivo adotado na célula e.

Observe que para obter Δε, primeiro deve-se calcular ΔK_N (definido na equação 16), devido aos deslocamentos com variação linear no contorno impostos pela deformação constante do macro nos nós do contorno. Então, usando-se a inversa da lei de Hooke, obtém-se Δε. Assim, para satisfazer a equação de equilíbrio (22) deve-se encontrar o incremento de flutuações que anula a mesma. Porém, para resolver o problema de equilíbrio do EVR devem ser impostas condições de contorno em termos de flutuação dos deslocamentos. Nesse trabalho, serão consideradas as mesmas condições de contorno propostas em [¹⁰], [¹²] e [¹³] que são: deslocamento linear no contorno, flutuações periódicas no contorno e forças uniformes no contorno. Cada uma dessas condições leva a resultados numéricos diferentes, ou seja, gera um modelo multi-escala diferente (para maiores detalhes ver [¹], [¹⁷] e [¹⁸]).

Como já comentado, o incremento de deformações Δε não é constante em todos os nós do EVR, como ocorria com a formulação em elementos finitos proposta em [¹⁰], [¹²] e [¹³], ou seja, no modelo do MEC aqui desenvolvido, a deformação imposta pelo macro gera automaticamente um campo de flutuação no interior do EVR. Assim, a deformação imposta pelo macro é constante apenas para os pontos do contorno, pois como os deslocamentos dos nós internos são incógnitas do problema, na solução elástica linear do EVR, os deslocamentos obtidos nos nós internos não serão iguais àqueles referentes a uma deformação Δε constante, se forem definidas diferentes fases ou vazios. Note que no caso de o EVR ser uniforme (composto apenas por uma matriz, sem inclusões ou vazios) a deformação Δε imposta pelo macro gerará um campo de deformações constante no EVR, já que todos os nós de célula têm as mesmas propriedades elásticas. Além disso, se o RVE é uniforme, a equação de equilíbrio (22) é satisfeita para um campo de flutuações nulo, mesmos se todos os nós tiverem plastificado, ou seja, não é necessário resolver o processo iterativo, sendo a resposta da estrutura idêntica a uma análise não-linear convencional.

A matriz tangente do problema bidimensional é obtida linearizando-se a equação (22) de resíduos, através do Método de Newton-Raphson, ou seja, para uma determinada iteração i a correção no estado de flutuação dos deslocamentos (ẟΔ(ũ)ⁱ_n) é calculada a partir da seguinte equação:

$${\left\{\delta \Delta \tilde{u}\right\}}_n^i={\left[-\frac{\partial {\left\{R\right\}}_n^{i-1}}{\partial {\left\{\tilde{u}\right\}}_n^{i-1}}\right]}^{-1}{\left\{R\right\}}_n^{i-1}$$ (23)

ou seja, a matriz tangente consistente é dada por:

$$K^i=-\frac{\partial {\left\{R\right\}}_n^{i-1}}{\partial {\left\{\tilde{u}\right\}}_n^{i-1}}=\underset{e=1}{\overset{N_{cel}}{\Sigma }}{B}_e^T{C}^{ep\left(e\right)}{B}_e{A}_{e}t$$ (24)

O sistema de equações definido na equação (23) pode ser reduzido de acordo com as condições de contorno adotadas em termos de flutuação dos deslocamentos (para maiores detalhes ver [¹], [¹⁷] e [¹⁸]). Assim, na forma reduzida, a eq. (23) pode ser escrita na forma da eq. (25) de onde se obtêm as flutuações incógnitas ẟũ_Rⁱ⁺¹:

$$\delta {{\tilde{u}}_R}^{i+1}=-K_R^{i^{-1}}{F_R}^i$$ (25)

onde os dois vetores e a matriz K_R variam de acordo com o modelo multi-escala (deslocamento linear no contorno, flutuações periódicas no contorno ou forças uniformes no contorno).

Agora será mostrado como se obtém a resposta constitutiva homogeneizada do EVR, ou seja, como se faz a passagem da microestrutura para a macroestrutura. Para isso, admite-se que a tensão σ=σ(x,t) de um ponto x do macro-contínuo é definida como sendo a homogeneização do campo de tensões microscópicas σ_µ=σ_µ(y,t) que agem no EVR associado a esse ponto x, ou seja:

$$\sigma \left(x,t\right)=\frac{1}{V_{\mu }}{\int }_{{\Omega }_{\mu }}{\sigma }_{\mu }\left(y,t\right)dV$$ (26)

onde σ=σ(x,t) é a tensão macroscópica ou homogeneizada.

Na definição de σ (equação 26) o EVR é considerado como um meio contínuo e, portanto, o conceito de tensão permanece válido na escala microscópica. Assim, o EVR deve ser suficientemente grande para que a representação dele como meio contínuo faça sentido. No entanto, suas dimensões são completamente independentes daquelas definidas para estrutura estudada na macro-escala. Observe que há outros tipos de modelos para análise em multi-escala (que não será considerado aqui) em que a tensão macroscópica é relacionada diretamente às forças médias moleculares ou atômicas. Nesse caso, o EVR não é considerado como meio contínuo.

A tensão homogeneizada σ é calculada em termos das forças de superfície do contorno do EVR, através da seguinte equação (maiores detalhes ver [¹⁰], [¹²], [¹⁷], [¹⁸]):

$$\sigma =\frac{1}{2V_{\mu }}\left(\overline{\sigma }+{\overline{\sigma }}^T\right)$$ (27)

onde considerando-se que o contorno foi discretizado em N_b e definindo F como as forças dos nós sobre o contorno e y o vetor das coordenadas x₁ e x₂ de um ponto y do EVR, tem-se:

$$\overline{\sigma }=\underset{i=1}{\overset{Nb}{\Sigma }}{\left\{F\right\}}_i{\left\{y\right\}}_i^T$$ (28)

As forças F são obtidas a partir das forças de contorno P_B (por unidade de comprimento). Assim, adotando-se L_e-1 e L_e como sendo, respectivamente, os comprimentos dos elementos aos quais o nó genérico i do contorno pertence, tem-se:

$${\left\{F\right\}}_i={\left\{P_c\right\}}_i\frac{1}{2}\left(L_{e-1}+L_e\right)$$ (29)

onde para um incremento n, tem-se:

$${\left\{P_c\right\}}_n={\left\{P_c\right\}}_{n-1}+{\left\{\Delta P_c\right\}}_n$$ (30)

sendo as forças ΔP_c calculadas pela equação (6), após resolver o problema de equilíbrio do EVR.

O tensor constitutivo homogeneizado também é obtido considerando-se um processo de homogeneização (ver mais detalhes em [¹⁰], [¹²], [¹⁷] e [¹⁸]), sendo calculado através da seguinte equação:

$$C^{ep}\left(x\right)=\frac{1}{V_{\mu }}{\int }_{{\Omega }_{{\mu }^S}}\frac{\partial {\sigma }_{\mu }\left(y\right)}{\partial \epsilon \left(x\right)}dV=\frac{1}{V_{\mu }}{\int }_{{\Omega }_{{\mu }^S}}{C}_{\mu }^{ep}dV+\frac{1}{V_{\mu }}{\int }_{{\Omega }_{{\mu }^S}}{C}_{\mu }^{ep}\frac{\partial {\nabla }^S{\tilde{u}}_{\mu }}{\partial \epsilon \left(x\right)}dV=C^{ep\left(Taylor\right)}+{\tilde{C}}^{ep}$$ (31)

onde:

$$C^{ep\left(Taylor\right)}=\frac{1}{V_{\mu }}{\int }_{{\Omega }_{{\mu }^S}}{C}_{\mu }^{ep}dV=\underset{e=1}{\overset{N_{cel}}{\Sigma }}\frac{V_e}{V_{\mu }}{C}_{\mu }^{ep\left(e\right)}$$ (32)

$${\tilde{C}}^{ep}=-\frac{1}{V_{\mu }}{G}_R{K_R}^{-1}{G_R}^T$$ (33)

$$G=\underset{e=1}{\overset{Ncel}{\Sigma }}{C}_{\mu }^{ep\left(e\right)}{B}_e{V}_e$$ (34)

Nas equações (32), (33) e (34) $C_{\mu }^{ep}$ é o operador tangente consistente, obtido com a lei constitutiva incremental microscópica, ou seja, nesse caso C^ep(Taylor) representa a média volumétrica do tensor constitutivo incremental. As matrizes G_R e K_R são formas reduzidas das matrizes G (eq. 34) e K (eq. 24).

Maiores detalhes da implementação e o algoritmo do modelo aqui desenvolvido podem ser encontrados em [¹].

3. RESULTADOS E DISCUSSÃO

Nesta seção diferentes microestruturas (EVRs) foram adotadas, cujas respostas constitutivas homogeneizadas são comparadas com aquelas obtidas a partir do modelo do MEF para a micro-estrutura proposto em [¹⁰], [¹²], e [¹³] a fim de validar a formulação aqui desenvolvida. Note que a formulação multi-escala proposta em [¹⁰], [¹²] e [¹³] e implementada em [¹⁷], [¹⁸], [¹⁹] e [²⁰], já foi publicada em diversos trabalhos que têm coautoria do Prof. Eduardo A. de Souza Neto e seus colaboradores, por isso foi escolhida para validar os resultados no presente trabalho. Nessa comparação, adotou-se sempre a mesma malha para ambos os modelos do MEC e do MEF, ou seja, para o modelo em elementos finitos, os elementos triangulares foram adotados coincidentes com as células triangulares do modelo desenvolvido com o MEC. Vale ressaltar que apesar de não serem apresentados aqui, testes de convergência dos resultados com o refinamento da malha foram testados tanto para o modelo do MEC quanto do MEF. Como os dois modelos numéricos são bem diferentes, adotou-se uma malha que apresentasse bons resultados para ambos os métodos, a fim de que se pudesse fazer a comparação dos resultados. É importante lembrar novamente que todas as análises são feitas apenas no nível da microestrutura, isto é, uma análise multi-escala totalmente acoplada não é considerada nesse trabalho. Dessa forma, para simular o processo incremental de carga referente ao problema do macro-contínuo numa análise multi-escala, uma deformação qualquer é imposta ao EVR em incrementos. Então, após a solução do EVR para cada incremento de deformação, os valores homogeneizados do vetor de tensão e do tensor constitutivo são obtidos. É importante observar que uma análise sobre a convergência dos resultados com o refinamento da malha foi feita em todos os exemplos. Assim, quando esse estudo não é apresentado, utilizou-se sempre uma malha cujos resultados já haviam alcançado essa convergência.

3.1 Influência da fração de volumes das inclusões

Nesse primeiro conjunto de EVRs, uma inclusão elástica foi definida no centro do EVR, sendo o comportamento do material da matriz governado pelo modelo constitutivo de von Mises. Diferentes frações de volume foram consideradas para essa inclusão, a fim de verificar como isso influenciava na resposta constitutiva do EVR. Adotou-se em todos os cálculos o modelo de flutuações periódicas no contorno do EVR. As seguintes frações de volume foram adotadas: fv =10%, fv =30% and fv =37% (veja Figura 2). As seguintes propriedades foram adotadas para as inclusões, que são adotadas como elásticas: Módulo de elasticidade E=200GPa e coeficiente de Poisson ν=0,2. Para a matriz adotou-se: E=70GPa, ν=0,3, tensão de escoamento σ_y=243MPa e modulo de encruamento K=2,24GPa. O seguinte vetor de deformação foi imposto ao EVR em 25 incrementos: {ε}={ε₁ ε₂ ε₁₂}={-0,0015 0,0048 0}

Figura 2
EVRs com diferentes frações de volume (fv) de inclusões: a) fv =10%, b) fv =30% e c) fv =37%

As malhas consideradas para o EVR com fv =10% , fv =30% and fv =37% são apresentadas na Figura (2) e contêm, respectivamente: a) 220 células e 131 nós (40 elementos de contorno e 12 elementos de interface), b) 204 células e 123 nós (40 elementos de contorno e 20 elementos de interface) and c) 236 células e 139 nós (40 elementos de contorno e 24 elementos de interface). Na Figura (3) tem-se a tensão homogeneizada na direção x₂ ao longo do processo incremental e considerando-se os três diferentes EVRs. Como esperado, a resposta mais rígida se refere ao EVR com fv =37% e a resposta mais flexível é referente ao EVR com fv=10%. Na Tabela 1, apresentam-se o tensor constitutivo e o vetor de tensão para o EVR com fv =37% referentes ao último incremento de deformação. Além disso, tanto na Figura (3) como na Tabela (1) nota-se que os resultados do MEC e do MEF são bem similares.

Figura 3
Tensão Homogeneizada na direção x2, considerando-se diferentes frações de volume para as inclusões

Tabela 1
Resposta Homogeneizada para o EVR com fv =37%

VETOR DE TENSÃO (MPA) - MEC	VETOR DE TENSÃO (MPA) - MEF
$\left[\begin{array}{ccc}14,64& 273,89& -0,009\end{array}\right]$	$\left[\begin{array}{ccc}14,25& 273,89& -0,012\end{array}\right]$
TENSOR CONSTITUTIVO (GPA) - MEC	TENSOR CONSTITUTIVO (GPA) - MEF
$\left[\begin{array}{ccc}82,434& 34,797& 0,0054\\ 34,797& 20,748& -0,0033\\ 0,0054& -0,0033& 26,389\end{array}\right]$	$\left[\begin{array}{ccc}82,533& 34,662& 0,0147\\ 34,662& 20,646& -0,0026\\ 0,0147& -0,0026& 26,785\end{array}\right]$

3.2 Influência da fração de volumes dos vazios

Nesse conjunto de exemplos numéricos pretende-se analisar o comportamento constitutivo de materiais porosos, onde um vazio é adotado no centro do EVR. Diferentes frações de volume (fv) foram adotadas para o vazio, a fim de verificar como o tamanho do vazio altera a reposta constitutiva do material estudado. Na Figura 4 têm-se os EVRS adotados com as seguintes frações de volume para o vazio: fv =10%, fv =30% e fv =37%, os quais possuem, respectivamente, as malhas: a) 540 células e 312 nós (64 elementos de contorno e 20 elementos de interface); b) 384 células e 240 nós (64 elementos de contorno e 32 elementos de interface); c) 356 células e 228 nós (64 elementos de contorno e 36 elementos de interface).

Figura 4
EVRs com diferentes frações de volume de vazios: a) fv =10%, b) fv=30% e c) fv =37%

O comportamento do material da matriz é governado pelo critério de Von Mises, onde adotam-se: E=70GPa, ν=0,3, σ_y=243MPa e K=2.24GPa. Além disso, o seguinte vetor de deformação foi imposto ao EVR em 25 incrementos: {ε}={ε₁ ε₂ ε₁₂}={-0,0015 0,0048 0,0006} A fim de estudar a convergência dos resultados com o refinamento da malha, quatro malhas foram adotadas para cada EVR (ver Figura 5 para o EVR com fv =37%). Porém, nesse estudo é mostrado apenas o EVR com fv =37%, pois o mesmo tipo de comportamento foi observado para os outros EVRs. Como pode ser observado na Figura 6, mesmo a malha mais pobre já apresenta bons resultados.

Figura 5
Diferentes discretizações para o EVR com fv =37%

Figura 6
Convergência dos resultados com o refinamento da malha para o EVR com fv =37%

Na Figura 7 apresenta-se a tensão homogeneizada na direção x₂, ao longo do processo incremental de deformação, para os três EVRs definidos na Figura 4. Observa-se que o EVR com fv =10% apresenta o comportamento mais rígido, como esperado, já que representa a microestrutura de um material com menor porosidade. É importante notar o quanto a tensão σ₂ diminuiu ao se comparar a Figura 7 com a Figura 3, que é referente ao EVR com inclusões elásticas. Isso foi devido à substituição de uma inclusão elástica por um vazio, diminuindo a rigidez do EVR. Ainda comparando as Figuras 7 e 3, pode-se reparar o salto de tensão que acontece entre os EVRs (a) e (b). No caso do EVR com inclusão elástica, aquele que possui apenas 10% de inclusão (EVR a) apresenta resposta mais flexível que o EVRs (b) (com 30% de inclusão) e o EVR (c) (com 37% de inclusão). Porém a diferença entre os resultados é bem menor se comparado com os resultados dos EVRs (a), (b) e (c) da figura 7, que possuem vazios nos locais das inclusões. Isso mostra que a variação da fração de volume de vazios influenciou mais na rigidez do EVR do que a variação da fração de volume de inclusão. Além disso, na Tabela 2 tem-se o tensor constitutivo e o vetor de tensão para o EVR com fv =37%. Nota-se que na Figura 7 e na Tabela 2 os valores referentes ao modelo do MEC são semelhantes com aqueles obtidos com o MEF.

Figura 7
Tensão Homogeneizada na direção x₂, considerando-se diferentes frações de volume para os vazios

Tabela 2
Resposta Homogeneizada para o EVR com fv =37%.

VETOR DE TENSÃO (MPA) - MEC	VETOR DE TENSÃO (MPA) - MEF
$\left[\begin{array}{ccc}-23,01& 89,54& 2,019\end{array}\right]$	$\left[\begin{array}{ccc}-28,08& 89,31& 1,139\end{array}\right]$
TENSOR CONSTITUTIVO (GPA) - MEC	TENSOR CONSTITUTIVO (GPA) - MEF
$\left[\begin{array}{ccc}21,100& 1,5059& -0,0338\\ 1,5059& 3,0702& -0,4360\\ -0,0338& -0,4360& 2,3492\end{array}\right]$	$\left[\begin{array}{ccc}21,570& 1,4072& -0,0429\\ 1,4072& 2,9924& -0,4109\\ -0,0429& -0,4109& 2,3052\end{array}\right]$

Note que a maior componente do vetor de tensão, assim como o tensor constitutivo do modelo MEF são muito similares ao modelo proposto. Já a primeira componente de tensão apresenta maior diferença, mas não se espera que essa diferença resulte em diferenças significativas quando uma análise multi-escala acoplada de uma estrutura for realizada. Como os dois modelos numéricos são muito diferentes, não se espera que eles reproduzam exatamente o mesmo resultado, mas que tenham apenas resultados similares. As principais diferenças entre os modelos são: i) no MEC as forças elásticas não são aproximadas no domínio, apenas as forças inelásticas, ii) no modelo do MEC a mesma ordem de aproximação é adotada para os deslocamentos e forças no contorno, enquanto no MEF a função que aproxima as forças tem ordem menor que aquela dos deslocamentos; iii) para obter a solução elástica do MEC os deslocamentos devem ser aproximados no domínio apenas no caso de se ter fases com diferentes coeficientes de Poisson iv) no modelo proposto a deformação é imposta ao contorno do EVR, o que gera um campo de flutuações no domínio mesmo antes de resolver o problema de equilíbrio do EVR, enquanto que no modelo do MEF a deformação constante é imposta a todos os nós do EVR o que leva a ter um campo de flutuações inicialmente nulo.

3.3 Influência da condição de contorno em termos de flutuação dos deslocamentos

Considere agora o EVR com 5 inclusões elásticas representado na Figura 8, onde tem-se definidas quatro diferentes discretizações a fim de estudar a convergência dos resultados com o refinamento da malha. Na ordem da esquerda para a direita, essas discretização são as seguintes: : 204 células e 123 nós (40 elementos de contorno e 48 elementos de interface); 488 células e 277 nós (64 elementos de contorno e 72 elementos de interface); 880 células e 481 nós (80 elementos de contorno e 96 elementos de interface); 1222 células e 658 nós (92 elementos de contorno e 116 elementos de interface).

Figura 8
Discretizações adotadas para o EVR com cinco inclusões

Nesse item será estudado como a resposta constitutiva do EVR se modifica de acordo com o modelo multi-escala adotado, isto é, de acordo com as condições de contorno adotadas em termos de flutuações dos deslocamentos, que podem ser: deslocamentos lineares, flutuações periódicas e forças uniformes. Como já mostrado em diversos outros trabalhos, como por exemplo [¹¹], [²³] e [²⁴], as forças uniformes resultam na resposta mais flexível, enquanto deslocamentos lineares produzem a resposta mais rígida, já que se refere ao modelo mais restringido.

Para as inclusões elásticas são adotadas as seguintes propriedades: módulo de elasticidade E=200GPa e coeficiente de Poisson ν=0,2. Para o material da matriz considerou-se um comportamento perfeitamente plástico e governado pelo critério de von Mises, além das seguintes propriedades: E=70GPa, ν=0,3 e tensão de escoamento σy=243MPa. Além disso, foi imposto ao EVR, em 25 incrementos, o seguinte vetor de deformações: {ε}={ε₁ ε₂ ε₁₂}={-0,0015 0,0048 0,0006}

Figura 9
Convergência dos resultados com refinamento da malha para o EVR com 5 inclusões, considerando-se flutuações periódicas no contorno

Observa-se na Figura 9 que a malha mais pobre já apresenta bons resultados e que as duas malhas mais refinadas fornecem resultados praticamente idênticos. Em análises numéricas, quanto mais refinadas as malhas, mais precisos os resultados. Assim, tais testes de convergência são importantes para mostrar que ao se adotarem malhas mais finas os valores convergem para um determinado valor que permanece inalterado mesmo se malhas mais finas forem consideradas. Isso demonstra que o modelo é estável. A tensão homogeneizada ao longo do processo incremental de deformações é apresentada na Figura 10, onde observa-se que para as três condições de contorno os resultados do MEC se comparam muito bem com os do MEF. É também interessante notar que, como o material da matriz tem comportamento perfeitamente plástico, a resposta homogeneizada do EVR também se assemelha a um comportamento perfeitamente plástico, tendo sido obtidas, respectivamente, para deslocamento linear, flutuação periódica e forças uniformes as seguintes tensões últimas: 259MPa, 252,5MPa e 239,5 MPa. Como esperado, o modelo com forças uniformes apresenta o menor valor de tensão limite e aquele com deslocamento linear o maior valor para essa tensão. No caso do presente exemplo, observou-se diferença entre os três diferentes tipos de modelos multi-escala apenas em relação à carga limite, ou seja, houve diferença apenas na resistência final do material, mas não houve diferença na sua rigidez, já que as curvas da Figura 10 se encontram praticamente superpostas. Porém, há casos, como por exemplo se no EVR adotado nesse exemplo fossem considerados vazios no lugar das inclusões (ver [¹]), que se tem também uma diferença significativa na rigidez do EVR de acordo com o modelo multi-escala adotado.

Figura 10
Tensão homogeneizada na direção x₂ para o EVR com 5 inclusões considerando-se diferentes condições de contorno em termos de flutuações de deslocamentos

Na Tabela 3 é apresentada a resposta homogeneizada para a deformação ε₂=0,00288 e considerando-se cada uma das três condições de contorno, onde nota-se que os resultados do MEF são muito similares àqueles referentes ao modelo aqui apresentado do MEC.

Tabela 3
Resposta homogeneizada para o EVR com 5 inclusões (referente à deformação ε₂=0,00288), considerando-se diferentes condições de contorno em termos de flutuações dos deslocamentos

DESLOCAMENTOS LINEARES	FLUTUAÇÕES PERIÓDICAS	FORÇAS UNIFORMES
TENSOR CONSTITUTIVO (GPA)
	MEC
$\left[\begin{array}{ccc}100,83& 32,303& 0,1035\\ 32,303& 42,301& -2,9089\\ 0,1035& -2,9089& 35,049\end{array}\right]$	$\left[\begin{array}{ccc}99,318& 33,585& -0,0025\\ 33,585& 40,543& -2,8085\\ -0,0025& -2,8085& 34,559\end{array}\right]$	$\left[\begin{array}{ccc}94,532& 40,157& 0,4716\\ 40,157& 33,813& -3,9532\\ 0,4716& -3,9532& 34,073\end{array}\right]$
	MEF
$\left[\begin{array}{ccc}100,86& 32,085& 0,0542\\ 32,085& 41,257& -2,9410\\ 0,0542& -2,9410& 35,096\end{array}\right]$	$\left[\begin{array}{ccc}99,191& 33,556& 0,0105\\ 33,556& 40,465& -2,7136\\ 0,0105& -2,7136& 34,660\end{array}\right]$	$\left[\begin{array}{ccc}92,963& 41,691& 0,6173\\ 41,691& 30,517& -4,1411\\ 0,6173& -4,1411& 33,951\end{array}\right]$
VETOR DE TENSÃO (MPA)
	MEC
$\left[\begin{array}{ccc}-11,66& 240,44& 11,91\end{array}\right]$	$\left[\begin{array}{ccc}-10,24& 239,39& 11,67\end{array}\right]$	$\left[\begin{array}{ccc}-4,25& 239,52& 11,58\end{array}\right]$
	MEF
$\left[\begin{array}{ccc}-12,30& 240,77& 11,98\end{array}\right]$	$\left[\begin{array}{ccc}-10,76& 239,78& 11,75\end{array}\right]$	$\left[\begin{array}{ccc}-4,42& 239,72& 11,60\end{array}\right]$

4. CONCLUSÕES

Nos exemplos numéricos, diferentes EVRs sujeitos a deformações representativas foram analisados, onde o desempenho da formulação proposta foi confirmado comparando o vetor de tensão e o tensor constitutivo homogeneizado com a formulação implementada em [¹⁷], [¹⁸], [¹⁹] e [²⁰], pelo Método dos Elementos Finitos (MEF). Em um trabalho futuro, acoplando-se essa formulação a um modelo do MEC para análise do macrocontínuo, a análise em multi-escala de estruturas composta por materiais heterogêneos pode ser realizada usando somente o Método dos Elementos de Contorno. Como a análise em multi-escala é muito cara computacionalmente, é muito importante investigar diferentes formulações, a fim de obter um modelo numérico menos caro computacionalmente, mas ainda preciso.

O custo computacional da modelagem proposta foi comparado ao modelo de elementos finitos desenvolvido em [¹⁷] e [¹⁸], para evidenciar uma vantagem importante do modelo de elementos de contorno. Foi observado que o tempo de processamento para resolver o processo iterativo do EVR pelo MEC foi menor do que o MEF. Como comparação pode-se evidenciar o EVR com 880 células da Figura 8. Definindo t_MEC como o tempo necessário para resolver o processo iterativo pelo MEC e t_MEF de forma análoga para o MEF, foi obtida a razão t_MEC / t_MEF = 0,55. Essa diferença é muito significativa, principalmente quando em trabalhos futuros for considerada uma análise multi-escala completa, pois nesse caso um problema de equilíbrio do EVR deve ser solucionado para cada ponto do macro-contínuo em uma determinada iteração do processo incremental de carga da estrutura. É importante dizer que essa análise do tempo de processamento foi feita para diferentes EVRs, havendo pouca diferença quando diferentes micro-estrutras eram adotadas. Como comentado mais detalhadamente nos exemplos numéricos, os dois modelos numéricos são muito diferentes o que influencia significativamente no tempo de processamento.

Nos exemplos numéricos, com relação à influência do volume de inclusões e vazios no comportamento do EVR, percebe-se que a diferença de tensões homogeneizadas é mais significativa quando se aumenta o volume de vazios, ou seja, os vazios influenciam mais na rigidez do EVR do que as inclusões. Como esperado, o aumento do volume de inclusão leva a um aumento da rigidez e da resistência do EVR, enquanto o aumento do volume de vazios o torna menos rígido e resistente.

Com relação à adoção de diferentes modelos multi-escala, verificou-se como esperado, que aquele que considera forças uniformes no contorno do EVR resulta num EVR com menor resistência. Porém, para o exemplo estudado a diferença não foi significativa para as outras duas condições de contorno, apesar do modelo que considera deslocamento linear no contorno conseguir absorver maiores valores de deformação.

Portanto, de maneira geral, pode-se concluir que a formulação do MEC foi capaz de simular o comportamento mecânico de microestruturas heterogêneas considerando-se não-linearidade física, com resultados próximos ao MEF, um método já consolidado, o que nos auxilia a validar a formulação do MEC aqui apresentada.

BIBLIOGRAFIA

• ¹
  FERNANDES G. R., CROZARIOL L. H. R., SOARES A. F., et al, “A 2D Boundary Element Formulation to Model the Constitutive Behaviour of Heterogeneous Microstructures Considering Dissipative Phenomena”, Engineering Analysis of Boundary Elements, v. 99, pp. 1-22, Feb. 2019.
• ²
  GAL, E., KRYVORUK, R. “Fiber reinforced concrete properties – a multiscale approach”, Computers and Concrete, v.8, n. 5, pp. 525-539, Oct. 2011.
• ³
  TERADA, K., HORI, M., KYOYA, T., et al, “Simulation of the multiscale convergence in computational homogenization approaches”, International Journal of Solids and Structures, v.37, pp. 2285–2311, Apr. 2000.
• ⁴
  NGUYEN, V.P., LLOBERAS VALLS, O., STROEVEN, M., et al, “Homogenization-based multiscale crack modelling: from micro-diffusive damage to macro-cracks”, Computer Methods in Applied Mechanics and Engineering, v.200, pp.1220–1236, Feb. 2011.
• ⁵
  KOUZNETSOVA, V., GEERS, M.G.D., BREKELMANS, W.A.M., “Multi-scale second-order computational homogenization of multi-phase materials: a nested finite element solution strategy”, Computer Methods in Applied Mechanics and Engineering, v.193, pp.5525–5550, Dec. 2004.
• ⁶
  NEMAT-NASSER, S., HORI, M., Micromechanics: Overall Properties of Heterogeneous Materials, Amsterdam, Elsevier Science, 1999.
• ⁷
  SANTOSH, M., GHOSH, M. A. “Multi-scale identification of concrete material parameters.”, Theoretical and Applied Fracture Mechanics, v.75, pp. 8-15, Feb. 2015.
• ⁸
  TORO, S., SÁNCHEZ, P. J., BLANCO, P. J., et al, “Multiscale formulation for material failure accounting for cohesive cracks at the macro and micro scales.”, International Journal of Plasticity, V.76, pp. 75-110, Jan. 2016.
• ⁹
  EDMANS, B. D., ALFANO, G., BAHAI, H., “Nonlinear multi-scale homogenization with different structural models at different scales.”, International Journal for Numerical Methods in Engineering v. 94, pp. 355-373, Apr. 2013.
• ¹⁰
  PERIC D., SOUZA NETO E. A., FEIJÓO R., et al, “On Micro-to-Macro Transitions for Multiscale Analysis of Heterogeneous Materials: Unified Variational Basis and Finite Element Implementation.”, International Journal for Numerical Methods in Engineering v.87, pp. 149-170, Jul. 2011.
• ¹¹
  GIUSTI S.M., BLANCO P.J., SOUZA NETO, E.A., et al, “An assessment of the Gurson yield criterion by a computational multi-scale approach.”, Engineering Computations, v. 26, n. 3, pp. 281-301, Apr. 2009.
• ¹²
  SOUZA NETO, E.A. & FEIJÓO, R.A. “Variational foundations of multi-scale constitutive models of solid: Small and large strain kinematical formulation”, In: Internal Research & Development Report, National Laboratory for Scientific Computing (LNCC/MCT), Brazil, No. 16, 2006.
• ¹³
  SOMER D.D., SOUZA NETO E.A., DETTMER W.G., et al, “A sub-stepping scheme for multi-scale analysis of solids.”, Computer Methods in Applied Mechanics and Engineering, v. 198, pp. 1006–1016, Feb. 2009.
• ¹⁴
  ALIABADI MH. The Boundary Element Method: applications in solids and structures, v.2, Wiley, 577, pp., 2002.
• ¹⁵
  BESKOS D.E., Boundary element analysis of plates and shells Berlin, Springer Verlag, 1991.
• ¹⁶
  BREBBIA C.A., TELLES J.C.F., WROBEL L.C. Boundary element techniques. Theory and applications in engineering, Berlin and New York, Springer-Verlag, 1984.
• ¹⁷
  FERNANDES G. R., PITUBA J. J. C., SOUZA NETO E. A. “Multi-Scale Modelling For Bending Analysis of Heteregeneous Plates by Coupling BEM AND FEM.”, Engineering Analysis of Boundary Elements v.51 pp.1-13, Feb. 2015.
• ¹⁸
  FERNANDES G. R., PITUBA J. J. C., SOUZA NETO E. A. “FEM/BEM formulation for multi-scale analysis of stretched plates.”, Engineering Analysis of Boundary Elements v. 54 pp.47-59, May. 2015.
• ¹⁹
  FURTADO, A. S., PITUBA, J. J. C., FERNANDES, G. R. “Analysis of the heterogeneous microstructures influence on the macromechanical response of the two-dimensional problem of plates.” Materia-Rio de Janeiro, v.22, e11830 - , Jun. 2017.
• ²⁰
  SOUZA JUNIOR, N. S., FERNANDES, G. R., PITUBA, J. J. C., “Multi-scale formulation for analysis of the plate bending problem considering dissipative processes in the microstructure and coupling BEM and FEM.”, Materia-Rio de Janeiro, v.22, e11821 - , June. 2017.
• ²¹
  FERNANDES, G.R., SOUZA NETO, E.A. “Self-consistent linearization of non-linear BEM formulations with quadratic convergence”, Computational Mechanics, v.52, pp.1125-1139, Nov. 2013.
• ²²
  FERNANDES G. R., ROSA NETO J. “Analysis of stiffened plates composed by different materials by the boundary element method.”, Structural Engineering and Mechanics, an International Journal, v.56, n.4 pp.605-623, Nov. 2015.
• ²³
  OHLAND, G. A., VIEIRA, J. F., FERNANDES, G. R. “A Boundary Element Formulation to perform elastic analysis of heterogeneous microstructures.”, Engineering Analysis of Boundary Elements, v.87, pp.47-65, Feb. 2018.
• ²⁴
  DORO, V. C., FERNANDES, G. R. “Análise não linear de chapas através de uma formulação do método dos elementos de contorno com convergência quadrática.” Materia-Rio de Janeiro, v.21, pp. 27-48, Mar. 2016.

Datas de Publicação

• Publicação nesta coleção
  06 Abr 2020
• Data do Fascículo
  2020

Histórico

• Recebido
  14 Jan 2019
• Aceito
  27 Ago 2019