<figure><img src="image_1.png" /></figure>

**HAL**
open science

# Recent Advances in Material Homogenization

Hajer Lamari, Amine Ammar, Patrice Cartraud, Francisco Chinesta, Frédéric Jacquemin, Grégory Legrain

► To cite this version:

Hajer Lamari, Amine Ammar, Patrice Cartraud, Francisco Chinesta, Frédéric Jacquemin, et al. Recent Advances in Material Homogenization. 13th ESAFORM Conference on Material Forming, 2010, Brechia, Italy. 10.1007/s12289-010-0913-y. hal-01008739

HAL Id: hal-01008739

https://hal.science/hal-01008739

Submitted on 25 Mar 2020

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# RECENT ADVANCES IN MATERIAL HOMOGENIZATION

H. Lamari<fnref n="1" />*, A. Ammar<fnref n="2" />, P. Cartraud<fnref n="1" />, F. Chinesta<fnref n="1" />, F. Jacquemin<fnref n="3" />, G. Legrain<fnref n="1" />

**ABSTRACT:** Heterogeneous materials involve different length scales in their mechanical properties. Obviously a mechanical description taking into account all the microscopic details is impossible from a computational point of view except for parts of very small dimensions. The main aim of material homogenization is defining macroscopic homogeneous properties able to represent at the macroscopic scale the real material and allowing for ignoring the microscopic scale in the numerical representation.

**KEYWORDS:** Homogenization, Model reduction, Separated representations, Proper Generalized Decompositions, Parametric models.

## 1 INTRODUCTION

The interest and issues of material homogenization are today well established. Two important recurrent issues concern the definition of homogenized properties when the microstructure evolves from one point to other at the macroscopic scale, that is, throughout the macroscopic part. This evolution of the microstructure needs the definition of a homogenized model at each node (or integration point) in the macroscopic part, even in the linear case. The second issue concerns precisely the eventual non-linearity of thermomechanical properties. In that case moreover, one should solve a homogenization problem at each location throughout the part as well as at each time step because the evolution of microscopic properties with the model solution itself.

In this work we present a technique able to alleviate the difficulty related to the evolution of the microstructure in the macroscopic scale. This technique is based on the consideration of a representative volume element consisting of different cells whose properties will be introduced as extra-coordinates of the thermo-mechanical model. Thus, the original thermo-mechanical model is transformed into a multi-dimensional parametric model suffering the redoubtable curse of dimensionality illness. The use of the proper generalized decomposition – PGD - allows circumventing this difficulty [1,2]. Section 2 revisits the application of PGD on parametric models. Some preliminary results involving heterogeneous linear thermal behaviours, perfectly defined at the microscopic scale, will be presented and discussed in section 3.

## 2 THERMAL PARAMETRIC MODEL

In this section we illustrate the application of the Proper Generalized Decomposition for solving a parametric thermal model involving a single extra-coordinate, the thermal conductivity characterizing an isotropic homogenous medium. This procedure will be extended later for including many extra-coordinates.

Consider the heat transfer equation:

$$ \frac{\partial u}{\partial t} - k \Delta u - f = 0 \qquad (1) $$

where $(x, t, k) \in \Omega \times I \times \mathfrak{T}$ and for the sake of simplicity the source term is assumed constant, i.e. $f = cte$. Because the conductivity is considered unknown, it is assumed as a new coordinate defined in the interval $\mathfrak{T}$. Thus, instead of solving the thermal model for different values of the conductivity parameter we prefer introducing it as a new coordinate. The price to be paid is the increase of the model dimensionality; however, as the complexity of the PGD scales linearly with the space dimension the consideration of the conductivity as a new coordinate allows for faster and cheaper solutions.

The solution of Eq. (1) is searched under the form:

$$ u(x, t, k) \approx \sum_{i=1}^{N} X_i(x) \cdot T_i(t) \cdot K_i(k) \qquad (2) $$

In what follows we are assuming that the approximation at iteration $n$ is already done:

<footnote n="1">¹GEM, Centrale Nantes, 1 rue de la Noe, BP 92101, F-44321 Nantes cedex 3, France</footnote>

<footnote n="2">²Laboratoire de Rhéologie, 13 rue de la Piscine, BP 53, F-38041 Grenoble cedex 9, France</footnote>

<footnote n="3">³GEM, Université de Nantes, 58, rue Michel Ange, BP 420, F-44606 Saint Nazaire cedex, France</footnote>

<footnote>* Hajer Lamari: GEM, Centrale nantes, 1 rue de la Noe, BP 92101, F-44321 Nantes cedex 3, France, Hajer.Lamari@ec-nantes.fr</footnote>
$$u^n(\mathbf{x}, t, k) = \sum_{i=1}^{n} X_i(\mathbf{x}) \cdot T_i(t) \cdot K_i(k) \quad (3)$$

and at present iteration we look for the next functional product $X_{n+1}(\mathbf{x}) \cdot T_{n+1}(t) \cdot K_{n+1}(k)$ that for alleviating the notation will be denoted by $R(\mathbf{x}) \cdot S(t) \cdot W(k)$. Prior to solve the resulting non linear model related to the calculation of these three functions a model linearization is compulsory. The simplest choice consists in using an alternating directions fixed point algorithm. It proceeds by assuming $S(t)$ and $W(k)$ given at the previous iteration of the non-linear solver and then computing $R(\mathbf{x})$. From the just updated $R(\mathbf{x})$ and $W(k)$ we can update $S(t)$, and finally from the just computed $R(\mathbf{x})$ and $S(t)$ we compute $W(k)$. The procedure continues until reaching convergence. The converged functions $R(\mathbf{x})$, $S(t)$ and $W(k)$ allow defining the searched functions:

$$X_{n+1}(\mathbf{x}) = R(\mathbf{x}), \quad T_{n+1}(t) = S(t) \quad \text{and} \\ K_{n+1}(k) = W(k).$$

We are illustrating each one of the just referred steps:

### I. Computing $R(\mathbf{x})$ from $S(t)$ and $W(k)$:

We consider the global weak form of Eq. (1):

$$\int_{\Omega \times I \times \mathfrak{T}} u^* \left( \frac{\partial u}{\partial t} - k \Delta u - f \right) d\mathbf{x} \, dt \, dk = 0 \quad (4)$$

where the trial and test functions write respectively:

$$u(\mathbf{x}, t, k) = \sum_{i=1}^{i=n} X_i(\mathbf{x}) \cdot T_i(t) \cdot K_i(k) + \\ + R(\mathbf{x}) \cdot S(t) \cdot W(k) \quad (5)$$

and

$$u^*(\mathbf{x}, t, k) = R^*(\mathbf{x}) \cdot S(t) \cdot W(k) \quad (6)$$

Introducing (5) and (6) into (4) it results:

$$\int_{\Omega \times I \times \mathfrak{T}} R^* \cdot S \cdot W \cdot \left( R \cdot \frac{\partial S}{\partial t} \cdot W - k \cdot \Delta R \cdot S \cdot W \right) d\mathbf{x} \, dt \, dk = \\ = - \int_{\Omega \times I \times \mathfrak{T}} R^* \cdot S \cdot W \cdot \mathfrak{R}^{(n)} d\mathbf{x} \, dt \, dk \quad (7)$$

$$\text{with } \mathfrak{R}^{(n)} = \sum_{i=1}^{i=n} X_i \cdot \frac{\partial T_i}{\partial t} \cdot K_i - \sum_{i=1}^{i=n} k \cdot \Delta X_i \cdot T_i \cdot K_i - f.$$

Now, being known all the functions involving the time and the parametric coordinate, we can integrate Eq. (7) in their respective domains $I \times \mathfrak{T}$. Integrating in $I \times \mathfrak{T}$ and taking into account the notation

$$\left[ \begin{array}{lll} w_1 = \int_{\mathfrak{T}} W^2 dk & s_1 = \int_I S^2 dt & r_1 = \int_{\Omega} R^2 d\mathbf{x} \\ w_2 = \int_{\mathfrak{T}} kW^2 dk & s_2 = \int_I S \cdot \frac{dS}{dt} dt & r_2 = \int_{\Omega} R \cdot \Delta R \, d\mathbf{x} \\ w_3 = \int_{\mathfrak{T}} W \, dk & s_3 = \int_I S \, dt & r_3 = \int_{\Omega} R \, d\mathbf{x} \\ w_4^j = \int_{\mathfrak{T}} W \cdot K_i \, dk & s_4^j = \int_I S \cdot \frac{dT_i}{dt} dt & r_4^j = \int_{\Omega} R \cdot \Delta X_i \, d\mathbf{x} \\ w_5^i = \int_{\mathfrak{T}} kW \cdot K_i \, dk & s_5^i = \int_I S \cdot T_i \, dt & r_5^i = \int_{\Omega} R \cdot X_i \, d\mathbf{x} \end{array} \right] \quad (8)$$

Eq. (7) reduces to:

$$\int_{\Omega} R^* \cdot (w_1 \cdot s_2 \cdot R - w_2 \cdot s_1 \cdot \Delta R) \, d\mathbf{x} = \\ = - \int_{\Omega} R^* \cdot \left( \sum_{i=1}^{i=n} w_4^i \cdot s_4^i \cdot X_i - \sum_{i=1}^{i=n} w_5^i \cdot s_5^i \cdot \Delta X_i - w_3 \cdot s_3 \cdot f \right) \, d\mathbf{x} \quad (9)$$

Eq. (9) defines an elliptic steady state boundary value problem that can be solved by using any discretization technique operating on the model weak form (finite elements, finite volumes ...). Another possibility consists in coming back to the strong form of Eq. (9):

$$w_1 \cdot s_2 \cdot R - w_2 \cdot s_1 \cdot \Delta R = \\ = - \left( \sum_{i=1}^{i=n} w_4^i \cdot s_4^i \cdot X_i - \sum_{i=1}^{i=n} w_5^i \cdot s_5^i \cdot \Delta X_i - w_3 \cdot s_3 \cdot f \right) \quad (10)$$

that could be solved by using any collocation technique (finite differences, SPH ...).

### II. Computing $S(t)$ from $R(\mathbf{x})$ and $W(k)$:

In the present case the test function writes:

$$u^*(\mathbf{x}, t, k) = S^*(t) \cdot R(\mathbf{x}) \cdot W(k) \quad (11)$$

Now, the weak form reads

$$\int_{\Omega \times I \times \mathfrak{T}} S^* \cdot R \cdot W \cdot \left( R \cdot \frac{\partial S}{\partial t} \cdot W - k \cdot \Delta R \cdot S \cdot W \right) d\mathbf{x} \, dt \, dk = \\ = - \int_{\Omega \times I \times \mathfrak{T}} S^* \cdot R \cdot W \cdot \mathfrak{R}^{(n)} d\mathbf{x} \, dt \, dk \quad (12)$$

that integrated in the domain $\Omega \times \mathfrak{T}$ and taking into account the notation (8) results:

$$\int_{I} S^* \cdot \left( w_1 \cdot r_1 \cdot \frac{dS}{dt} - w_2 \cdot r_2 \cdot S \right) dt = \\ = - \int_{I} S^* \cdot \left( \sum_{i=1}^{i=n} w_4^i \cdot r_5^i \cdot \frac{dT_i}{dt} - \sum_{i=1}^{i=n} w_5^i \cdot r_4^i \cdot T_i - w_3 \cdot r_3 \cdot f \right) dt \quad (13)$$

Eq. (13) represents the weak form of the ODE defining the time evolution of the field $S$ that can be solved by using any stabilized discretization technique (SU, Discontinuous Galerkin, ...). The strong form of Eq. (13) reads:
$$ 
\begin{aligned}
& w_1 \cdot r_1 \cdot \frac{dS}{dt} - w_2 \cdot r_2 \cdot S = \\
& = - \left( \sum_{i=1}^{i=n} w_4^i \cdot r_5^i \cdot \frac{dT_i}{dt} - \sum_{i=1}^{i=n} w_5^i \cdot r_4^i \cdot T_i - w_3 \cdot r_3 \cdot f \right)
\end{aligned} 
\quad (14) $$

than can be solved by using backward finite differences, or higher order Runge-Kutta schemes, among many other possibilities.

### III. Computing $W(k)$ from $R(x)$ and $S(t)$:

In the present case the test function writes:

$$ u^*(\mathbf{x}, t, k) = W^*(t) \cdot R(\mathbf{x}) \cdot S(k) \quad (15) $$

Now, the weak form reads:

$$ 
\begin{aligned}
& \int_{\Omega \times I \times \bar{\Omega}} W^* \cdot R \cdot S \cdot \left( R \cdot \frac{\partial S}{\partial t} \cdot W - k \cdot \Delta R \cdot S \cdot W \right) d\mathbf{x} \, dt \, dk = \\
& = - \int_{\Omega \times I \times \bar{\Omega}} W^* \cdot R \cdot S \cdot \mathfrak{R}^{(n)} d\mathbf{x} \, dt \, dk \quad (16)
\end{aligned}
 $$

that integrated in $\Omega \times I$ and taking into account the notation (8) results:

$$ 
\begin{aligned}
& \int_{\bar{\Omega}} W^* \cdot (r_1 \cdot s_2 \cdot W - r_2 \cdot s_1 \cdot k \cdot W) \, dk = \\
& = - \int_{\bar{\Omega}} W^* \cdot \left( \sum_{i=1}^{i=n} r_5^i \cdot s_4^i \cdot K_i - \sum_{i=1}^{i=n} r_4^i \cdot s_5^i \cdot k \cdot K_i - r_3 \cdot s_3 \cdot f \right) \, dk \quad (17)
\end{aligned}
 $$

Eq. (17) does not involve any differential operator. The strong form of Eq. (17) reads:

$$ 
\begin{aligned}
& (r_1 \cdot s_2 - r_2 \cdot s_1) \cdot W = \\
& = - \left( \sum_{i=1}^{i=n} (r_5^i \cdot s_4^i - r_4^i \cdot s_5^i) \cdot K_i - r_3 \cdot s_3 \cdot f \right) \quad (18)
\end{aligned}
 $$

that represents an algebraic equation. Thus, the introduction of parameters as additional model coordinates has not a noticeable effect in the computational cost, but can increase the number of sums involved in finite sums decomposition (2).

There are other minimization strategies more robust and exhibiting faster convergence for building-up the PGD (see [3]).

## 3 NUMERICAL EXAMPLES

In this section we are focusing in the thermal models defined in heterogeneous materials. Imagine a composite material involving a matrix and a fibrous reinforcement. A typical representative volume showing the microstructure heterogeneity is depicted in figure 1.

<figure><img src="image_3.png" /><figcaption>Figure 1. Microstructure of a composite material</figcaption></figure>

Now, in order to apply the PGD one should perform a separated representation of the thermal conductivity, allowing for an efficient thermal simulation, i.e.

$$ k(\mathbf{x}) \approx \sum_{i=1}^{P} K_x^i(x) \cdot K_y^i(y) \quad (19) $$

This separated representation can be performed by applying the SVD (singular value decomposition) to the matrix containing as entries the image pixels.

However, because the irregular distribution of the inclusions the number of sums in (19) becomes very high. Different possibilities exist to alleviate this representation some of them are being analyzed at present. A first possibility lies in considering only a reduced number of the modes of the SVD. The computational cost is significantly reduced and sometimes the precision is not too much degraded. Figure 2 depicts the image reconstruction for different number of sums in (20).

<figure><img src="image_4.png" /><figcaption>Figure 2. Reconstructed microstructure for $P = 23$ (left) and $P = 46$ in (19).</figcaption></figure>

Another possibility lies in the substitution of each circle by a square parallel to the coordinate axes, whose area equals the one of the associated circle. Eq. (20) represents a generic rectangle:

$$ \mathbf{x} \in [x_a, x_b] \times [y_a, y_b] \quad (20) $$

Now, we define the characteristic function:

$$ \chi_{a,b}(s) = \begin{cases} 0 & \text{if } s < a \\ 1 & \text{if } a \le s \le b \\ 0 & \text{if } s > b \end{cases} \quad (21) $$

Using this notation, the square defined in (20) can be represented in a separated form by:
$$ \chi_{x_a, x_b}(x) \times \chi_{y_a, y_b}(y). $$

Obviously, if the representative volume contains P inclusions, each one represented by a square, the separated representation (19) will contain P sums.

Obviously, all these approaches imply the solution of a thermal model for each representative volume (as soon as the microstrucuture evolve, the thermal conductivity is modified and then a new solution of the thermal model is required). If we consider a stochastic nature of the microstructure many realizations of the microstructure must be solved. One possibility for alleviating this task consists of considering the representative volume composed of a certain number of cells related to a grid of the representative volume as depicted in figure 3. Now, rather than solving a thermal model for each existing microstructure (represented by the different values of the thermal conductivity in each cell), we are introducing the thermal conductivity of each cell as an extra coordinate. In the example depicted in Fig. 3 the model will be defined in a space of dimension 27 (the x and y space coordinates and the 5×5 cells thermal conductivities).

<figure><img src="image_5.png" /><figcaption>Figure 3. Cellular microstructure of the representative volume element.</figcaption></figure>

Thus, the solution is searched under the form:

$$ u(x, y, k_{1,1}, L, k_{5,5}) \approx \sum_{i=1}^{N} X_i(x) \cdot Y_i(y) \cdot K_{1,1}^i(k_{1,1}) L \cdot K_{5,5}^i(k_{5,5}) \quad (22) $$

Thus, the thermal field for any possible microstructure only needs the solution of only one multidimensional problem, solution that can be performed efficiently by applying the PGD. As soon as the solution (22) is computed, the thermal field for any microstructure realization is obtained by substituting the known conductivities of each cell in Eq. (22). Thus, the thermal field related to the microstructure shown in Fig. 3 is depicted in Fig. 4 (a simple boundary condition

$$ u(\mathbf{x} \in \partial\Omega, k_{1,1}, L, k_{5,5}) = x \quad \text{was enforced} $$

<figure><img src="image_6.png" /><figcaption>Figure 4. Thermal field associated with the microstructure depicted in Fig. 3.</figcaption></figure>

## 4. CONCLUSIONS

This paper explores some possibilities related to the use of Proper Generalized Decomposition in the advanced simulation. This novel discretization technique allows the efficient solution of highly multidimensional models by circumventing the redoubtable curse of dimensionality that mesh based discretization techniques suffer.

As soon as an efficient solver for multidimensional models is available, many models of computational mechanics can be rewritten in higher dimensional spaces. Thus, one could for example compute the thermal field for any value of the thermal conductivity (as illustrated in this paper for homogeneous and heterogeneous materials) simplifying inverse identification of optimization.

## REFERENCES

[1] A. Ammar, B. Mokdad, F. Chinesta, R. Keunings, A new family of solvers for some classes of multidimensional partial differential equations encountered in kinetic theory modeling of complex fluids, J. Non-Newtonian Fluid Mech., 139: 153-176, 2006.

[2] A. Ammar, B. Mokdad, F. Chinesta, R. Keunings, A new family of solvers for some classes of multidimensional partial differential equations encountered in kinetic theory modeling of complex fluids. Part II: transient simulation using space-time separated representations, J. Non-Newtonian Fluid Mech., 144: 98-121, 2007.

[3] F. Chinesta, A. Ammar, E. Cueto Recent Advances in the Use of the Proper Generalized Decomposition for Solving Multidimensional Models. Archives of Computational Methods in Engineering - State of the Art Reviews, In press.