# A Manifold Learning Approach for Integrated Computational Materials Engineering

Elena Lopez, David Gonzalez, Jose Aguado, Emmanuelle Abisset-Chavanne, Elias Cueto, Christophe Binetruy, Francisco Chinesta

► To cite this version:

Elena Lopez, David Gonzalez, Jose Aguado, Emmanuelle Abisset-Chavanne, Elias Cueto, et al. A Manifold Learning Approach for Integrated Computational Materials Engineering. Archives of Computational Methods in Engineering, 2018, 25 (1), pp.59-68. 10.1007/s11831-016-9172-5. hal-02265692

HAL Id: hal-02265692

https://hal.science/hal-02265692

Submitted on 11 Aug 2019

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_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# A manifold learning approach for Integrated Computational Materials Engineering

E. Lopez¹, D. Gonzalez², J.V. Aguado¹, E. Abisset-Chavanne¹, F. Lebel³, R. Upadhyay³, E. Cueto², C. Binetruy¹ F. Chinesta¹,*

Received: October 2012 / Accepted: date

**Abstract** Image-based simulation is becoming an appealing technique to homogenize properties of real microstructures of heterogeneous materials. However fast computation techniques are needed to take decisions in a limited time-scale. Techniques based on standard computational homogenization are seriously compromised by the real-time constraint. The combination of model reduction techniques and high performance computing contribute to alleviate such a constraint but the amount of computation remains excessive in many cases. In this paper we consider an alternative route that makes use of techniques traditionally considered for machine learning purposes in order to extract the manifold in which data and fields can be interpolated accurately and in real-time and with minimum amount of online computation. Locally Linear Embedding – LLE - is considered in this work for the real-time thermal homogenization of heterogeneous microstructures.

This work has been partially supported by the Spanish Ministry of Science and Competitiveness, through grants number CICYT-DPI2014-51844-C2-1-R. Professor Chinesta is also supported by the Institut Universitaire de France.

E. Lopez & J.V. Aguado & E. Abisset-Chavanne & C. Binetruy
GEM, UMR CNRS - Centrale Nantes

1 rue de la Noe, BP 92101, F-44321 Nantes cedex 3, France

E-mail: {Elena.Lopez,Jose.Aguado-Lopez,Emmanuelle.Abisset-Chavanne,Christophe.Binetruy}@ec-nantes.fr

E. Cueto & D. Gonzalez

I3A, Universidad de Zaragoza

Maria de Luna s/n, 50018 Zaragoza, Spain

E-mail: {ecueto,gonzal}@unizar.es

F. Lebel & R. Upadhyay

Composite Technologies, GE Global Research, One research circle, Niskayuna, NY 12309, USA

E-mail: {Francois.upadhyay}ge.com

F. Chinesta

ESI Chair "Advanced Computational Manufacturing Processes"

GEM, UMR CNRS - Centrale Nantes

Institut Universitaire de France

1 rue de la Noe, BP 92101, F-44321 Nantes cedex 3, France

E-mail: Francisco.Chinesta@ec-nantes.fr

*Corresponding author: Francisco Chinesta, E-mail: Francisco.Chinesta@ec-nantes.fr
**Keywords:** Real time thermal simulation, Composite materials, Model Order Reduction, Computational homogenization, Locally Linear Embedding, Machine learning, manifold learning

# 1 Introduction

Image-based simulation is an appealing technique to homogenize properties of real microstructures of heterogeneous materials such as composite materials. Macroscale homogenized properties constitute useful information to adapt manufacturing parameters to the actual microstructure of the material to be processed. This method is known as adaptive manufacturing and is acknowledged as a practical way to reduce defects in parts. However fast computation techniques are needed to take decisions in a limited time-scale.

Moreover, the possibility of analyzing one part in real time reduces significantly the uncertainty and consequently allows to make better predictions, limiting the necessity of quantifying and propagating uncertainty all along the process.

When proceeding with heterogeneous materials consisting of disordered inclusions (fibers) into a matrix (polymer) the final thermal or mechanical homogenized properties strongly depend on the fraction of inclusions, their shape and their spatial distribution. When all of them vary significantly within the part, simulations at the macroscopic scale must extract homogenized material properties locally by performing efficient microscopic calculations at different representative volumes.

This procedure implies the necessity of solving complex numerical models in real-time. While such calculations can be envisaged by using nowadays powerful computational capabilities, the objective here is to deploy such capabilities in the production plant, where accurate solutions are needed in real-time and by using computing devices as cheap as possible. There have been a number of model order reduction attempts in this direction, see among many others [9, 15, 13, 8, 1] and references therein. However, none of them seem to have accomplished the challenge of real-time analysis and decision-making because of the multiscale nature of the problem to be solved as discussed below.

In this paper we consider a linear, homogeneous and isotropic material occupying the domain $\Omega \in \mathbb{R}^2$ (the extension to 3D is straightforward) containing a series of circular inclusions composed of a different linear, homogeneous and isotropic material, in which a linear thermal problem is defined. It consists of determining the temperature field $T(\mathbf{X}, t)$, with $\mathbf{X} \in \Omega$, at each time $t \in (0, T]$, fulfilling heat conduction equation

$$ \rho C \frac{\partial T(\mathbf{X}, t)}{\partial t} = \nabla \cdot (\mathbf{K}(\mathbf{X})\nabla T(\mathbf{X}, t)) = 0, \quad t \in (0, T], \quad \mathbf{X} \in \Omega, \quad (1) $$

with appropriate initial and boundary conditions. In Eq. (1) $\mathbf{K}(\mathbf{X})$ refers to the conductivity tensor at position $\mathbf{X}$, $\rho$ is the material density and $C$ the specific heat.

The description of heterogeneous materials, when this heterogeneity involves the fine scale, needs a resolution, high enough for capturing both the distribution of the material constituents (heterogeneity) and its effects on the temperature field. However, this high-resolution description requires too many degrees of freedom even for nowadays computational capabilities.

One possibility for circumventing this issue consists of separating the scales (whenever possible) and associating to each point $\mathbf{X} \in \Omega$ a representative volume at the finest scale $\omega(\mathbf{X})$, where the homogenized conductivity $\mathbf{K}(\mathbf{X})$ is calculated.
If the material is assumed to be composed of inclusions immersed into a matrix, both having a negligible evolution of their conductivities with temperature, then the resulting thermal model can be assumed linear. However, if the matrix is undergoing thermochemical or thermophyiscal transformations, the dependence of the conductiv- ity on the temperature cannot be ignored anymore, and the resulting thermal model becomes nonlinear.

In the linear case, the calculation of the homogenized conductivity tensor K(X) can
be performed as soon as the spacial distribution and thermal properties of constituents
(matrix and inclusions) are known, as detailed in A, requiring the solution of two (three
in the 3D case) boundary value problems – BVP – in ω(X) at each characteristic
microscopic representative volume Xᵢ, i = 1, ..., M. This homogenization procedure
was addressed using model reduction techniques [9, 45, 26, 18, 23, 27]

Efficient solutions of the thermal problem (1) can be accomplished by considering
any of the available model order reduction technique based on the Proper Orthogonal
Decomposition – POD – [6, 32, 36, 37, 38, 5, 7, 16, 43], Proper Generalized Decomposition
– PGD – [2, 3, 10, 11, 12, 22, 24, 25] or Reduced Bases – RB – [28, 29, 30, 35], that allow
real time simulations with the possibility of integrating them in Dynamic Data Driven
Application Systems – DDDAS – . When the microstructure is perfectly defined ho-
mogenization can be performed a priori and off-line. However, when the microstructure
is only available on-the-fly the situation is radically different because homogenization
must be performed in real-time for the microstructure just acquired, and as just men-
tioned and detailed in A, the calculation of the associated conductivity requires the
solution of two (or three) BVP. The main issue in those solutions is that as the in-
clusion distribution differs from one microstructure to other, the coefficients matrix
associated with the discrete thermal problem, can change substantially, implying its
assemblage and its solution for each new acquired microstructure. In [27] this issue was
deeply addressed and different approaches for alleviate the computational cost related
to these calculations proposed, in both the linear and nonlinear case. However all those
strategies were based in combining model reduction techniques with optimized data
and algebraic manipulation.

In this paper we consider a radically different route. POD, that is equivalent to PCA – Principal Components Analysis - can be viewed as an information extractor from a data set that attempts to find a linear subspace of lower dimensionality than the original space. If the data has more complicated structures which cannot be well represented in a linear subspace, standard PCA will not be very helpful. Fortunately, kernel PCA allows us to generalize standard PCA to nonlinear dimensionality reduction [44, 40, 41]. Locally Linear Embedding – LLE – [34] results from a particular choice of the kernel within the kPCA framework [46].

Our main objective in this work is to analyze the possibility of using these kind
of strategies widely employed for machine learning purposes, for performing suitable
interpolations on the data manifold from the last constructed from available offline
information, and like this infer in real-time and with a minimum amount of calculation
homogenized properties in heterogeneous microstructures.

In the next section we revisit LLE for data reconstruction based on the offline ex-
traction of the data manifold in which data, fields and PGD-based parametric solutions
can be safely interpolated. In section 3 different numerical solutions will be described
and discussed, proving the potentiality of the proposed approach.
## 2 Robust data reconstructors

It is well known that microstructures do not allow simple reduced descriptions. Imagine a domain containing a single inclusion. The system can be described by means of the phase field or the inclusion characteristic function (that takes a unit value at points located inside the inclusion vanishing outside). Imagine the same domain but now with the single inclusion located in a different region (with empty intersection with the region occupied by the one considered in the first system). This last situation can be described again by using the characteristic function related to the inclusion. Third, and finally, we consider the same domain containing again a single inclusion but in a different location that in both previous cases (the regions occupied by the inclusion in the three cases have null intersection). It is clear that the characteristic function related to the third inclusion cannot be written, in general, as a linear combination of the functions associated to the first two inclusions.

When instead of studying the microstructure itself we consider a field (e.g. the temperature field) the situation is similar. The different fields associated with different microstructures show some similarity but it is difficult to define an interpolated field in a new microstructure from the known fields and microstructures. In fact the main concern is how quantifying similarities or resemblances, and how to take profit of them.

In [27] we tried to extract significant information by considering different microstructures, solving the thermal models involved in computational homogenization in each one of them, and then applying a POD on these calculated solutions. When the number of microstructures is large enough there is only one significative mode that, as expected, corresponds to the Hill's mode (affine temperature). The remaining modes have all similar significance, making impossible the definition of a reduced basis accurate enough for expressing the solution of a new thermal problem defined in a different microstructure. For this reason in [27] more than using a reduced basis for solving the thermal problems involved in the new microstructures we simply improved the algebraic strategies for speeding-up the solution of those models in the new microstructures. In the nonlinear case the difficulty was twofold, because one must add to the previous one the one related to the nonlinear behavior, both making impossible the solution under the real time constraint.

In this paper we consider and analyze an alternative route based on the use of the Localy Linear Embedding – LLE – technique [34], a member of the large family of the so-called machine learning techniques. This technique has been widely used for performing proper data reconstruction, ensuring that the interpolated solution lives in the manifold defined by the solutions from which it is reconstructed.

We are describing the procedure directly in the problem we are interested in. The idea is that from some amount of calculation performed offline on a series of microstructures, we would like, for a new microstructure, to infer online its homogenized conductivity without solving any thermal problem in it, under the real-time constraint. The procedure consists of two steps: the offline analysis of many samples followed by an online data reconstruction, both described in the sequel.

### 2.1 Offline construction of the data manifold

First we assume the existence of $M$ microstructures $\mathcal{M}_m, m = 1, \dots, M$, defined in the RVE $\omega$. In what follows and without loss of generality we consider $2D$ microstructures
and temperature fields. Moreover, and again without loss of generality, we assume the existence of two phases, the circular inclusions and the continuous phase, occupying the domains $\omega_f^m$ and $\omega_p^m$ respectively, with $\omega_f^m \cup \omega_p^m = \omega$, $m = 1, \dots, M$. A regular mesh is associated to each RVE consisting in $N$ nodes ($N = N_n^2$, with $N_n$ the number of nodes along the $x$ and $y$ directions). The coordinates of each node are $\mathbf{x}_i$, $i = 1, \dots, N$ ($\mathbf{x}_i^T = (x_i, y_i)$).

For each microstructure $\mathcal{M}_m$ we define the phase field $\chi(\mathbf{x}; \mathcal{M}_m)$:

$$ \chi(\mathbf{x}; \mathcal{M}_m) = \begin{cases} 1 & \text{if } \mathbf{x} \in \omega_f^m \\ 0 & \text{if } \mathbf{x} \in \omega_p^m \end{cases} \qquad (2) $$

The microstructures can be represented in a discrete way from vectors $\chi^m$, whose $i$-th component writes $\chi_i^m = \chi(\mathbf{x}_i; \mathcal{M}_m)$. Vectors $\chi^m$ are defined in $\mathbb{R}^N$, i.e. the dimension coincides with the number of nodes considered in the discrete microscopic description.

Thus each vector $\chi^m$ defines a point in a space of dimension $N$, and then, the set of microstructures, represents a set of $M$ points in $\mathbb{R}^N$. The question that arises is: Do all these points belong to a certain low-dimensional manifold embedded in the high-dimensional space $\mathbb{R}^N$? Imagine that despite the impressive space dimension $N$, the $M$ points belong to a curve, a surface or a hyper-surface of dimension $d \ll N$. When $N = 3$ a simple observation suffices for checking if these points are located on a curve (one-dimensional manifold) or on a surface (two-dimensional manifold). However, when dealing with spaces of thousands of dimensions simple visual observation is unsuitable.

Instead, appropriate techniques are needed to extract the underlying manifold (when it exists) when proceeding in extremely multidimensional spaces. There is a variety of techniques for accomplishing this task. The interested reader can refer to [42, 34, 33, 44, 4]. In this work we focus on the LLE – Locally Linear Embedding – technology [34]. It proceeds as follows.

- Each point $\chi^m$, $m = 1, \dots, M$ is linearly reconstructed from its $K$-nearest neighbors. In principle $K$ should be greater than the expected dimension $d$ of the underlying manifold and the points should be close enough to ensure the validity of the linear approximation. In general, a large-enough number of neighbors $K$ and a large-enough sampling $M$ ensures a satisfactory reconstruction as proved later. For each point $\chi^m$ we can write the locally linear data reconstruction as:

$$ \chi^m = \sum_{i \in S_m} W_{mi} \chi^i, \qquad (3) $$

where $W_{mi}$ are the unknown weights and $S_m$ the set of the $K$-nearest neighbors of $\chi^m$.

As the same weights appears in different locally linear reconstructions, the best compromise is searched by looking for the weights, all them grouped in vector $\mathbf{W}$, that minimize the functional

$$ \mathcal{F}(\mathbf{W}) = \sum_{m=1}^{M} \left\| \chi^m - \sum_{i=1}^{M} W_{mi} \chi^i \right\|^2 \qquad (4) $$

where here $W_{mi}$ is zero if $\chi^i$ does not belong to the set of $K$-nearest neighbors of $\chi^m$.

The minimization of $\mathcal{F}(\mathbf{W})$ allows to determine all the weights involved in all the locally linear data reconstruction.
- We suppose now that each linear patch around $\chi^m$, $\forall m$, is mapped into a lower dimensional embedding space of dimension $d$, $d \ll N$. Because of the linear mapping of each patch, weights remain unchanged. The problem becomes the determination of the coordinates of each point $\chi^m$ when it is mapped into the low dimensional space, $\xi^m \in \mathbb{R}^d$.

For this purpose a new functional $\mathcal{G}$ is introduced, that depends on the searched coordinates $\xi^1, \dots, \xi^M$:

$$ \mathcal{G}(\xi^1, \dots, \xi^M) = \sum_{m=1}^{M} \left\| \xi^m - \sum_{i=1}^{M} W_{mi} \xi^i \right\|^2, \quad (5) $$

where now the weights are known and the reduced coordinates $\xi^m$ are unknown. The minimization of functional $\mathcal{G}$ results in a $M \times M$ eigenvalue problem whose $d$-bottom non-zero eigenvalues define the set of orthogonal coordinates in which the manifold is mapped.

Thus, considering a new point $\xi$ in the reduced space $\mathbb{R}^d$ after identifying its neighbors set $\mathcal{S}(\xi)$ and calculating the locally linear approximation weights, we can come back to $\mathbb{R}^N$ and reconstruct the phase field $\chi$ from its neighbors $\chi^i$, $i \in \mathcal{S}(\xi)$. However the reconstructed microstructure $\chi$ is no longer binary because it involves interpolation of binary vectors. Thus, reconstructed microstructures show an amount of spurious smoothing.

The reconstruction could be improved by substituting the phase field function by a more regular function able to identify the phase distribution with the same accuracy but from a continuous description. The simplest choice consists of using a level-set function $\phi(\mathbf{x})$, whose zero level coincides with the inclusion-matrix interfaces and its value in any other point $\mathbf{x}$ is simply the distance to the nearest interface. We proved in our numerical experiments that instead of using the phase field, when a level set-based description is employed, microstructures can be interpolated with better accuracy. However, if the objective is not to reconstruct geometrical data but smoother fields, a phase-field based description could be accurate enough.

## 2.2 Temperature field reconstruction

As previously discussed and detailed in A the homogenization process involves two thermal problems defined in $\omega(\mathbf{X})$ to obtain the temperature fields $\Theta^1(\mathbf{x})$ and $\Theta^2(\mathbf{x})$, associated with two particular choices of the boundary conditions, from which the homogenized thermal conductivity $\mathbf{K}(\mathbf{X})$ is defined. Suppose that both thermal problems have been solved offline for each microstructure $M_m$, $m = 1, \dots, M$. Now, from $\Theta_m^1(\mathbf{x})$ and $\Theta_m^2(\mathbf{x})$, $m = 1, \dots, M$ we can calculate the localization tensor for each microstructure and then the resulting homogenized thermal conductivity $\mathbf{K}^m$.

The new concern is: for a new microstructure $M$ is it possible to infer, without further thermal calculations, the homogenized conductivity tensor? Obviously as soon as the associated temperature fields $\Theta^1(\mathbf{x}); M)$ and $\Theta^2(\mathbf{x}; M)$ are available, as accurately as possible, the homogenized thermal conductivity $\mathbf{K}^M$ can be calculated.

First, we must prove that the new microstructure belongs to the manifold defined by the $M$ microstructures previously analyzed. For that purpose it suffices to determine $\chi$ and its reduced image $\xi$ and check if $\xi$ belongs to the $d$-dimensional manifold. If it
is the case, as previously explained, we can reconstruct $\xi$ from its nearest neighbors in the set $S(\xi)$ from a local and linear approximation. The reconstruction involves the weights $W_{\xi i}$, $i \in S(\xi)$. It is supposed here that the temperature field related to the new microstructure $\Theta^1(\mathbf{x}; \mathcal{M})$ and $\Theta^2(\mathbf{x}; \mathcal{M})$ can be linearly interpolated from the temperature fields $\Theta_i^1(\mathbf{x})$ and $\Theta_i^2(\mathbf{x})$, with $i \in S(\xi)$ by using the weights $W_{\xi i}$.

## 2.3 Checking the accuracy

The following test is proposed to check the accuracy of the whole procedure. We compare the discrete temperature fields known at each microstructure $\Theta_m^1$ and $\Theta_m^2$, $m = 1, \dots, M$, with the ones linearly reconstructed from its neighbors by using the mapping weights $W_{mi}$. For that purpose we define the errors $\mathcal{E}^1(M)$ and $\mathcal{E}^2(M)$

$$ \mathcal{E}^1(M) = \left( \frac{1}{M} \sum_{m=1}^{M} \left\| \Theta_m^1 - \sum_{i \in S_m} W_{mi} \Theta_i^1 \right\|^2 \right)^{\frac{1}{2}} \quad (6) $$

$$ \mathcal{E}^2(M) = \left( \frac{1}{M} \sum_{m=1}^{M} \left\| \Theta_m^2 - \sum_{i \in S_m} W_{mi} \Theta_i^2 \right\|^2 \right)^{\frac{1}{2}} \quad (7) $$

with the total error $\mathcal{E} = \mathcal{E}^1 + \mathcal{E}^2$.

If the error is low enough then it can be stressed that as soon as $\xi$ falls inside the convex hull of $\xi^m$, $m = 1, \dots, M$, the locally linear data reconstruction can be considered accurate enough.

## 2.4 Addressing nonlinearities

As described in [27] when addressing nonlinear models, linear homogenization must be carried out by freezing the thermal conductivity in each phase in $\omega(\mathbf{X})$ to the macroscopic temperature existing at the position $\mathbf{X}$ to which the RVE is attached $\omega(\mathbf{X})$. To avoid the solution of the numerous thermal problems that the nonlinearity involves until reaching convergence, authors proposed in [26, 27] the calculation of the parametric solutions $\Theta^1(\mathbf{x}, T(\mathbf{X}))$ and $\Theta^2(\mathbf{x}, T(\mathbf{X}))$ where $T(\mathbf{X})$ refers to the macroscopic temperature existing at point $\mathbf{X}$ where the RVE $\omega(\mathbf{X})$ is located. This parametric solution was performed within the PGD framework using the rationale described in [11]. The resulting strategy consists of solving both temperature fields for any macroscopic temperature within a certain predefined interval $T$.

As soon as both parametric solutions are known, the homogenized conductivity can be calculated in real time for a given macroscopic temperature $T(\mathbf{X})$. Thus, all the difficulties related to the nonlinearity can be efficiently circumvented, however a major issue remains. Both parametric solutions cannot be computed in real-time for any possible microstructure.

Following the rationale described in the previous sections we could compute $2 \times M$ parametric solutions related to the $M$ microstructures: $\Theta_m^1(\mathbf{x}, T)$ and $\Theta_m^2(\mathbf{x}, T)$, $m = 1, \dots, M$. Given a new microstructure $\mathcal{M}$, expressed as $\chi$, and after checking it
belongs to the samples manifold by calculating its image $\xi$. Then, weights $W_{\xi i}$ with respect to its K-nearest neighbors in the set $S(\xi)$ are calculated and from them the parametric solution is reconstructed (locally and linearly) according to

$$ \Theta^1(\mathbf{x}, T; \mathcal{M}) = \sum_{i \in S(\xi)} W_{\xi i} \Theta_i^1(\mathbf{x}, T), \quad (8) $$

and

$$ \Theta^2(\mathbf{x}, T; \mathcal{M}) = \sum_{i \in S(\xi)} W_{\xi i} \Theta_i^2(\mathbf{x}, T), \quad (9) $$

that define an appropriate interpolation of the parametric solution on the manifold. This procedure was successfully employed in [17] for interpolating parametric solutions in nonlinear elastic problems.

## 3 Numerical results

In the numerical tests we consider first $M = 30$ microstructures consisting of 100 circular inclusions randomly distributed in $\omega$. The inclusions volume fraction being $r = \frac{|\omega_f|}{|\omega_p|} \approx 0.5$. The temperature fields $\Theta_m^1(\mathbf{x})$ and $\Theta_m^2(\mathbf{x})$ were calculated by discretizing the thermal problems using a regular finite element mesh consisting of $N = 10^4$ nodes.

The $M$ samples were described by the phase field vectors $\chi^m$, all of them defined in $\mathbb{R}^N$. Then the weights involved in the linear data reconstruction were calculated as well as the mapping to the reduced space. The dimension of the manifold where the different microstructures $\chi^m$ were mapped results in $d = 2$. Figure 1 depicts the resulting points $\xi^m \in \mathbb{R}^3$, $m = 1, \dots, M$, from which one can realize that the dimension of the manifold is $d = 2$ since it corresponds to a plane. Fig. 2 represents these points on the plane. In order to check that by considering more microstructures the manifold is better filled $M = 100$ microstructures are now considered. Fig. 3 depicts the resulting reduced points in $\mathbb{R}^3$ and on the plane that constitutes the manifold.

Now, one of these points is randomly chosen, for instance the one surrounded by a red square in Figure 4. Fig. 5 shows the selected microstructure. As this point is related to one of the microstructures that were analyzed, the reference temperature field $\Theta^1(\mathbf{x})$ depicted in Figure 6 is actually known. The temperature at this microstructure is linearly interpolated from the K-nearest neighbors by considering the weights associated to the microstructure mapping

$$ \tilde{\Theta}^1(\mathbf{x}) = \sum_{i \in S(\xi)} W_{\xi i} \Theta_i^1(\mathbf{x}). \quad (10) $$

Fig. 7 shows the K-nearest microstructures (in our case we considered $K = 9$). A very good match between both solutions is found: $||\tilde{\Theta}^1(\mathbf{x}) - \Theta^1(\mathbf{x})|| = 0.0033$.

The same procedure was repeated many times in order to evaluate if the resulting error varies significantly when considering the different points on the manifold. Fig. 8 shows the errors (defined as previously) when reconstructing the temperature fields at indicated positions with respect to the reference thermal solutions. As it can be noticed the error does not evolve significantly even if for $M = 30$ microstructures the distribution of samples is rather heterogeneous.
<figure><img src="image_2.png" /><figcaption>Fig. 1 Two views of points $\xi^m$ in $\mathbb{R}^3$</figcaption></figure>

Finally we consider a new microstructure different to all the previously analyzed. As expected, its reduced coordinates $\xi$ belong to the previously reconstructed manifold, being its coordinates $\xi^T = (-1.14, -0.82, -1)$. This safety check, to verify that the new microstructure belongs to the manifold making possible further calculations, in particular the temperature fields interpolation, was carried out in about 0.3 seconds in a standard laptop using Matlab. The interpolation itself was performed in about 6 milliseconds. Both thermal problems were then solved in order to evaluate the accuracy of the interpolated fields, and the errors were again quite small, of around $10^{-3}$.

## 4 Conclusions

The proposed method employs manifold learning techniques for a near-optimal characterization of the microstructure composed of circular inclusions randomly distributed
<figure><img src="image_3.png" /><figcaption>Fig. 2 Microstructures $\xi^m$ on the two-dimensional manifold</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Fig. 3 Two views of points $\xi^m$ in $\mathbb{R}^3$ (left) and on the plane (right) when considering $M = 100$ microstructures</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Fig. 4 Selected microstructure were the temperature field is being reconstructed</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>Fig. 5 Microstructures $\xi^m$ on the two-dimensional manifold</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Fig. 6 Temperature reference solution for the selected microstructure</figcaption></figure>

in an homogeneous matrix. We have demonstrated that the employment of locally linear embedding techniques allows for a very efficient embedding of the microstructures in a way in which it is possible to properly interpolate between carefully selected "reference microstructures". It is thus avoided a costly and detailed analysis of every possible microstructure found by imaging techniques.
<figure><img src="image_8.png" /><figcaption>**Fig. 7** K-nearest microstructures ($K = 9$) from which the temperature field related to microstructure depicted in Fig. 5 is reconstructed</figcaption></figure>

This method proves to be reasonably efficient and to render accurate enough results with minimal computational cost. In fact, the application of LLE techniques to costly 2D or 3D images has proven to run in a few seconds, much less than the usual cost of a full, non-linear, finite element simulation.

# 5 References

## References

1. J.V. Aguado, A. Huerta, F. Chinesta, E. Cueto. Real-time monitoring of thermal processes by reduced order modelling. International Journal for Numerical Methods in Engineering, 102/5, 991-1017, 2015.

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. Journal of Non-Newtonian Fluid Mechanics, 139, 153-176, 2006.

3. 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
<figure><img src="image_9.png" /><figcaption>Fig. 8 Errors in the reconstructed thermal data with respect to the reference temperature fields</figcaption></figure>

of complex fluids. Part II: Transient simulation using space-time separated representation. Journal of Non-Newtonian Fluid Mechanics, 144, 98-121, 2007.

4. D. Amsallem, C. Farhat. Interpolation method for adapting reduced-order models and application to aeroelasticity. AIAA Journal, 46, 1803-1813, 2008.

5. D. Amsallem, J. Cortial and C. Farhat. Toward real-time CFD-based aeroelastic computations using a database of reduced-order information. AIAA Journal, 48, 2029-2037, 2010.

6. R.A. Bialecki, A.J. Kassab, A. Fic. Proper orthogonal decomposition and modal analysis for acceleration of transient FEM thermal analysis. Int. J. Numer. Meth. Engrg., 62, 774-797, 2005.

7. T. Bui-Thanh, K. Willcox, O. Ghattas, B. Van Bloemen Waanders. Goal-oriented, model-constrained optimization for reduction of large-scale systems. Journal of Computational Physics, 224/2, 880-896, 2007.

8. V.M. Calo, Y. Efendiev, J. Galvisd, M. Ghommem. Multiscale empirical interpolation for solving nonlinear PDEs. Journal of Computational Physics, 278, 204-220, 2014.

9. F. Chinesta, A. Ammar, F. Lemarchand, P. Beauchene, F. Boust. Alleviating mesh constraints: Model reduction, parallel time integration and high resolution homogenization. Computer Methods in Applied Mechanics and Engineering, 197/5, 400-413, 2008.

10. F. Chinesta, P. Ladeveze, E. Cueto. A short review in model order reduction based on Proper Generalized Decomposition. Archives of Computational Methods in Engineering, 18, 395-404, 2011.

11. F. Chinesta, A. Leygue, F. Bordeu, J.V. Aguado, E. Cueto, D. Gonzalez, I. Alfaro, A. Ammar, A. Huerta, Parametric PGD based computational vademecum for efficient design, optimization and control, Archives of Computational Methods in Engineering, 20/1, 31-59, 2013.

12. F. Chinesta, R. Keunings, A. Leygue, The Proper Generalized Decomposition for advanced numerical simulations. A primer, Springerbriefs, Springer, 2014.

13. M. Cremonesi, D. Neron, P.A. Guidault, P. Ladeveze. PGD-based homogenization technique for the resolution of nonlinear multiscale problems. Computer Methods in Applied Mechanics and Engineering, 267, 275-292, 2013.

14. M.G.D. Geers, V.G. Kouznetsova, W.A.M. Brekelmans. Multi-scale computational homogenization: Trends and challenges. Journal of Computational and Applied Mathematics, 2009, In press, DOI: 10.1016/j.cam.2009.08.077.

15. Ch. Ghnatios, F. Masson, A. Huerta, E. Cueto, A. Leygue, F. Chinesta. Proper Generalized Decomposition Based Dynamic Data-Driven Control of Thermal Processes. Computer Methods in Applied Mechanics and Engineering, 213, 29-41, 2012.
16. M. Girault, E. Videcoq, D. Petit. Estimation of time-varying heat sources through inversion of a low order model built with the Modal Identification Method from in-situ temperature measurements. International Journal of Heat and Mass Transfer, 53, 206-219, 2010.

17. D. Gonzalez, E. Cueto, F. Chinesta. Computational patient avatars for surgery planning. Annals of Biomedical Engineering. In press.

18. F. E. Halabi, D. González, A. Chico, M. Doblaré, FE² multiscale in linear elasticity based on parametrized microscale models using proper generalized decomposition, Computer Methods in Applied Mechanics and Engineering, 2013, 257, 183 - 202.

19. M. Jiang, I. Jasiuk, M. Ostoja-Starzewski. Apparent thermal conductivity of periodic two-dimensional composites. Computational Materials Science, 25/3, 329-338, 2002.

20. T. Kanit, S. Forest, I. Galliet, V. Mounoury, D. Jeulin. Determination of the size of the representative volume element for random composites: statistical and numerical approach. International Journal of Solids and Structures, 40/13-14, 3647-3679, 2003.

21. T. Kanit, F. N'Guyen, S. Forest, D. Jeulin, M. Reed, S. Singleton. Apparent and effective physical properties of heterogeneous materials: Representativity of samples of two materials from food industry. Computer Methods in Applied Mechanics and Engineering, 195/33-36, 3960-3982, 2006.

22. P. Ladevèze, The large time increment method for the analyze of structures with nonlinear constitutive relation described by internal variables, Comptes Rendus Académie des Sciences Paris, 309, 1095-1099, 1989.

23. P. Ladevèze, A. Nouy, On a multiscale computational strategy with time and space homogenization for structural mechanics, Computer Methods In Applied Mechanics and Engineering, 192/28-30, 3061-3087, 2003.

24. P. Ladevèze, D. Néron, P. Gosselet, On a mixed and multiscale domain decomposition method, Computer Methods in Applied Mechanics and Engineering, 96, 1526-1540, 2007.

25. P. Ladevèze, J.-C. Passieux, D. Néron, The latin multiscale computational method and the proper generalized decomposition, Computer Methods In Applied Mechanics and Engineering, 199/21-22, 1287-1296, 2010.

26. H. Lamari, A. Ammar, P. Cartraud, G. Legrain, F. Jacquemin, F. Chinesta. Routes for Efficient Computational Homogenization of Non-Linear Materials Using the Proper Generalized Decomposition. Archives of Computational Methods in Engineering, 17/4, 373-391, 2010.

27. E. Lopez, E. Abisset-Chavanne, F. Lebel, R. Upadhyay, S. Comas-Cardona, C. Bintruy, F. Chinesta. Advanced thermal simulation of processes involving materials exhibiting fine-scale microstructures. International Journal of Material Forming. In press, DOI 10.1007/s12289-015-1222-2

28. Y. Maday, E.M. Ronquist. A reduced-basis element method. C. R. Acad. Sci. Paris, Ser. I, vol. 335, 195-200, 2002.

29. Y. Maday, A.T. Patera, G. Turinici. A priori convergence theory for reduced-basis approximations of single-parametric elliptic partial differential equations. Journal of Scientific Computing, 17/1-4, 437-446, 2002.

30. Y. Maday, E.M. Ronquist. The reduced basis element method: application to a thermal fin problem. SIAM J. Sci. Comput., 26/1, 240-258, 2004.

31. M. Ostoja-Starzewski. Material spatial randomness: From statistical to representative volume element. Probabilistic Engineering Mechanics, 21/2, 112-132, 2006.

32. H.M. Park, D.H. Cho. The use of the Karhunen-Loève decomposition for the modelling of distributed parameter systems. Chem. Engineer. Science, 51, 81-98, 1996.

33. M. Polito, P. Perona. Grouping and dimensionality reduction by Locally Linear Embedding, in Advances in Neural Information Processing Systems 14, MIT Press, 1255-1262, 2001.

34. S.T. Roweis, L.K. Saul. Nonlinear dimensionality reduction by Locally Linear Embedding. Science, 290, 2323-2326, 2000.

35. G. Rozza, D.B.P. Huynh, A.T. Patera. Reduced basis approximation and a posteriori error estimation for affinely parametrized elliptic coercive partial differential equations – application to transport and continuum mechanics. Archives of Computational Methods in Engineering, 15/3, 229-275, 2008.

36. D. Ryckelynck, L. Hermanns, F. Chinesta, E. Alarcon. An efficient a priori model reduction for boundary element models. Engineering Analysis with Boundary Elements, 29, 796-801, 2005.

37. D. Ryckelynck. A priori hyperreduction method: an adaptive approach. Journal of Computational Physics, 202, 346-366, 2005.
38. D. Ryckelynck, F. Chinesta, E. Cueto, A. Ammar. On the a priori model reduction: Overview and recent developments. Archives of Computational Methods in Engineering, State of the Art Reviews, 13/1, 91-128, 2006.

39. K. Sab. On the homogenization and the simulation of random materials. European Journal of Mechanics, A/Solids, 11/5, 585-607, 1992.

40. B. Scholkopf, A. Smola, K.R. Muller. Nonlinear component analysis as a kernel eigenvalue problem. Neural Computation, 10/5, 1299-1319, 1998.

41. B. Scholkopf, A. Smola, K.R. Muller. Kernel principal component analysis. In Advances in Kernel Methods – Support Vector Learning, MIT Press, 327-352, 1999.

42. J.B. Tenenbaum, V. de Silva, J.C. Langford. A global framework for nonlinear dimension-ality reduction. Science, 290, 2319-2323, 2000.

43. E. Videcoq, O. Quemener, M. Lazard, A. Neveu. Heat source identification and on-line temperature control by a Branch Eigenmodes Reduced Model. International Journal of Heat and Mass Transfer, 51 4743-4752, 2008.

44. Q. Wang. Kernel Principal Component Analysis and its Applications in Face Recognition and Active Shape Models. arXiv:1207.3538

45. J. Yvonnet, D. Gonzalez, Q.-C. He, Numerically explicit potentials for the homogenization of nonlinear elastic heterogeneous materials, Computer Methods in Applied Mechanics and Engineering, 198/33-36, 2009, 2723-2737.

46. V.A. Zimmer, K. Lekadir, C. Hoogendoorn, A.F. Frangi, G. Piella. A framework for optimal kernel-based manifold embedding of medical image data. Computerized Medical Imaging and Graphics, 41, 93-107, 2015.

## A On computational homogenization

In this section the simplest procedure related to computational homogenization is revisited. Due to the microscopic heterogeneity, the macroscopic thermal modeling needs a homogenized thermal conductivity which depends on the microscopic details. To compute this homogenized thermal conductivity an appropriate RVE is considered at position $X \in \Omega$, $\omega(X)$ [20, ?, ?] in which the microstructure is perfectly defined at this scale.

In the linear case the local microscopic conductivity $k(x)$ is known at each point $x$ in the microscopic domain $\omega(X)$.

We can define the macroscopic temperature gradient at position $X$, $G(X)$, from:

$$G(X) = \langle g \rangle = \frac{1}{|\omega(X)|} \int_{\omega(X)} g(x) \, dx \quad (11)$$

where the temperature gradient writes $g(x) = \nabla T(x)$.

We also assume the existence of a localization tensor $L(x, X)$ such that

$$g(x) = L(x, X) \cdot G(X) \quad (12)$$

The microscopic heat flux $q$ writes according to Fourier's law

$$q(x) = -k(x) \cdot g(x) \quad (13)$$

and its macroscopic counterpart $Q(X)$ reads:

$$Q(X) = \langle q(x) \rangle = -\langle k(x) \cdot g(x) \rangle = -\langle k(x) \cdot L(x, X) \rangle \cdot G(X) \quad (14)$$

from which the homogenized thermal conductivity can be defined from

$$K(X) = \langle k(x) \cdot L(x, X) \rangle \quad (15)$$

Since $k(x)$ is perfectly known everywhere in the representative volume element $\omega(X)$, the definition of the homogenized thermal conductivity tensor only requires the computation of the localization tensor $L(x, X)$. Several approaches are proposed in the literature to define this tensor, according to the choice of boundary conditions. The objective here is not to discuss this choice. The interested reader can find some details in [39, 19, 20, 14]. For the sake of simplicity, we use essential boundary conditions on $\partial\omega(X)$ corresponding to the assumption of uniform
<figure><img src="image_10.png" /><figcaption>Fig. 9 Homogenization procedure of linear heterogeneous models.</figcaption></figure>

temperature gradient on the RVE $\omega(\mathbf{X})$. We consider the general 3D case that involves the solution of the three boundary-value problems related to the steady state heat transfer model in the microscopic domain $\omega(\mathbf{X})$ for three different boundary conditions on $\partial\omega(\mathbf{X})$:

$$ \left\{ \begin{array}{l} \nabla \cdot (\mathbf{k}(\mathbf{x}) \cdot \nabla \Theta^1(\mathbf{x})) = 0, \\ \Theta^1(\mathbf{x} \in \partial\omega(\mathbf{X})) = x \end{array} \right. , \qquad (16) $$

$$ \left\{ \begin{array}{l} \nabla \cdot (\mathbf{k}(\mathbf{x}) \cdot \nabla \Theta^2(\mathbf{x})) = 0, \\ \Theta^2(\mathbf{x} \in \partial\omega(\mathbf{X})) = y \end{array} \right. , \qquad (17) $$

and

$$ \left\{ \begin{array}{l} \nabla \cdot (\mathbf{k}(\mathbf{x}) \cdot \nabla \Theta^3(\mathbf{x})) = 0 \\ \Theta^3(\mathbf{x} \in \partial\omega(\mathbf{X})) = z \end{array} \right. \qquad (18) $$

It is easy to prove that these three solutions verify

$$ \left\{ \begin{array}{l} \mathbf{G}^1 = \langle \nabla \Theta^1(\mathbf{x}) \rangle = (1, 0, 0)^T \\ \mathbf{G}^2 = \langle \nabla \Theta^2(\mathbf{x}) \rangle = (0, 1, 0)^T \\ \mathbf{G}^3 = \langle \nabla \Theta^3(\mathbf{x}) \rangle = (0, 0, 1)^T \end{array} \right. \qquad (19) $$

where $(\cdot)^T$ denotes the transpose. Thus, the localization tensor results finally:

$$ \mathbf{L}(\mathbf{x}, \mathbf{X}) = (\nabla \Theta^1(\mathbf{x}) \quad \nabla \Theta^2(\mathbf{x}) \quad \nabla \Theta^3(\mathbf{x})) \qquad (20) $$

The resulting non-concurrent homogenization procedure is illustrated in Fig. 9. As soon as tensor $\mathbf{L}(\mathbf{x}, \mathbf{X})$ is known at each position $\mathbf{x}$, the constitutive law relating the macroscopic temperature gradient and the macroscopic heat flux becomes fully defined.