Reduced FFT-Based Simulation of a Mechanically Loaded Clustered Microstructure Using an Adaptive Set of Fourier Modes

Processes, such as deep rolling or induction hardening, have a remarkable influence on the material properties within the surface layer of a work piece. Our overall goal is to develop efficient two-scale methods, which are able to show the microstructural evolution of the machined material. The calculation of a spatially resolved microstructure comes along with a high computational effort. To reduce the computational costs, we combine a clustered description of the structure [1] with a model order reduction technique for the performed fast Fourier transformations (FFT) [2]. We choose a reduced set of Fourier modes, which is adapted to the underlying microstructure and thus based on the occurring strain field [3]. By that, we analyze the influence of a mechanical impact on an elasto-plastically deforming material.


Introduction
Our overall goal is a two-scale simulation of different thermo-mechanical machining processes in which we use the finite element (FE) method for the macroscopic boundary value problem and fast Fourier transforms (FFT) for the analysis at the micro scale, see for instance [4,5]. Nevertheless, in this article, we want to focus on model order reduction techniques applied to the microstructural computations.
Compared to the more popular FE 2 methods, the FE-FFT concept enables a promising reduction of the computational effort, see for instance [6,7,8]. The FFT approach is introduced by Moulinec and Suquet [2,9] who used a basic fixed-point scheme for the simulation of a microscopic composite structure. Further investigations by different research groups deal with the modifications of the solution scheme, which show an improved computational performance, such as the implementation of Newton-Krylov solvers [10,11,12], augmented Lagrangians [13,14,15] and polarization-based algorithms [16]. An overview of FFT-based computations and their applications can be found in [17].
Although the FFT-based computation already shows an improved performance compared to the classical FE simulation, the calculation of the highly resolved microstructure still needs high computational costs. Due to that, model order reduction techniques are a condition for efficient complex two-scale process simulations. Thus, in the current work, we want to present a coupled approach which combines the usage of a reduced set of Fourier modes with a clustered microstructure, e.g. [18,19].
After this short introduction, the next section gives an overview of the used FFT method as well as of the applied model order reduction technique, which is based on the idea of using a reduced set of Fourier modes. Afterwards the description of the applied clustering technique and of the implemented algorithm completes the presentation of the computational structure. Numerical results show the performance of our combined approach, before the last section finally gives a summary and an outlook on future investigations.

Notation.
Here are some remarks for the better understanding of the following notations: we are using Einstein's summation convention. This means that if one index occurs exactly twice, one can leave out the summation symbols: ∑ = . In addition, the double contraction denoted by (•): (•) of two tensors reads : = . Quantities that refer to the macroscopic scale are marked by (•) ��� , whereby variables related to the microscale waive a specific notation.

Microstructural Analysis Using Fast Fourier Transformations
General Concept. Within a two-scale simulation, every macroscopic point � corresponds to an underlying microstructure. Fig. 1 shows an exemplary inhomogeneous microstructure, in which the stiffness as well as inner variables, such as a plastic strain or a phase concentration, depend on the microstructural position . The total microscopic strain (�, ) accounts for two contributions: one of the macroscale �(�) and another one due to the fluctuations within the microstructure �(�, ). One can write For a better overview, we neglect the dependency on the macroscopic position in the following. Besides Eq. 1, the balance of linear momentum div ( ) = as well as a constitutive relation for the stress ( ) = ( , ( ), ( )), which might depend on a set of inner variables ( ), complete the heterogeneous boundary value problem. As presented in Eq. 2, the microscale is free of external forces, which already influence the macroscopic boundary value problem. Based on [20], we now introduce the so-called polarization stress which describes the stress deviation compared to a homogeneous reference material with stiffness ℂ 0 . Inserting the relation Eq. 3 into the balance equation Eq. 2 enables and hence, the handling of the divergence of the polarization stress as an external force. Based on the works by [21,22], the influence of on ( ) results from the convolution integral Herein, the so-called Green's operator Γ 0 ( − ′ ) maps the polarization stress at a point ′ on the strain at a point . The rather complex solution of the convolution integral motivates the calculations performed in Fourier's space based on [2,9]. In the following, all the quantities which are related to Fourier's space are given by (⋅) � . Additionally, corresponding to a coordinate in the real space, the quantities depend on Fourier modes . Accordingly, Eq. 5 simplifies to

Achievements and Trends in Material Forming
wherein Green's operator and function in Fourier's space read The back transformation of �( ) into the real space gives the total microscopic strain ( ). The calculations using Eq. 3 and Eq. 6 as well as the Fourier transformations and back transformations are repeated until convergence.
Reduced Set of Fourier Modes. The basic idea of the applied model order reduction technique is that not all the Fourier modes are necessary for the computations in Fourier space. Instead of using the firstly introduced fixed sampling pattern [23], we use a sampling pattern ℛ that is adapted to the underlying microstructure. We introduced two methods to identify the choice of Fourier modes [3,24].
The first method considers the underlying geometry of a material that consists of two different phases, namely a matrix material and inclusions. By the introduction of a characteristic function ( ) that is zero in the matrix Ω and one in the inclusions Ω and its transformation into Fourier's space, one can identify the modes with the highest amplitudes for this specific material. The choice of ℛ % of the modes with the largest amplitudes then defines the related sampling pattern, see Fig. 2. More details can be found in [24].
The geometrically adapted sampling pattern serves as a good choice for linear elastic materials. For a more complex material behavior, we introduced a strain-based sampling pattern in [3]. The idea is similar to the procedure in the geometrically adapted sampling pattern: we also choose the ℛ % modes with the highest amplitudes. Instead of the function ( ), the norm of the converged microscopic strain field of the previous load step leads to the updated and thus strain-based sampling pattern.
In the current work, we choose the second approach. Only for the first load step, we will apply the geometrically adapted reduced set of Fourier modes.

Clustering Approach and Coupled Algorithm
K-Means Clustering. In addition to the strain-based model order reduction technique for the calculations performed in Fourier space, we want to use a clustered microstructure based on k-means clustering to reduce the computational effort in real space. The works of [18,19] give an overview of the method. The starting point is the construction of a strain-concentration tensor ( , �) which describes the microscopic strain ( , �) resulting from a macroscopic strain �(�):
The main idea of the clustering technique is now to assign the grid points of the microstructure to the clusters in such a way that the sum of the differences of the stress concentration tensors of the individual grid points compared to the related averaged tensor for the clusters � is minimized. The minimization problem reads wherein contains the optimized sets of grid points distributed to each cluster .
Coupled Algorithm. Our current work combines the use of a strain-based sampling pattern for the FFT-based calculations and the clustering approach presented in the former sections. The calculation of the clusters as well as of the geometrically adapted sampling pattern for the first time step takes place in an initialization step. For each time step, we applied the algorithm presented in Fig. 3, which is based on a fixed-point iteration scheme by [2]. The transformation of the polarization stress into Fourier's space enables the calculation of the strain in Fourier space based on Eq. 6 and its back transformation in the real space. Afterwards, the strain is homogenized for each cluster by calculating the average for the grid points. For each cluster , we then calculate the stress and inner variables. We repeat this procedure until convergence of the strain field. Finally, the Fourier transformation of the norm of the converged strain field leads to the updated sampling pattern for the next load step.

Numerical Results
To show the performance of the coupled algorithm, we present the numerical result for a microstructure consisting of two different materials as presented in Fig. 4 (left). The elastic parameters for the inclusions (black) read = 2.0 GPa and = 2.0 GPa. The matrix (gray) with = 1.0 GPa and = 1.0 GPa is also able to deform plastically with a yield limit Thereby, = 256 × 256 grid points and an unreduced set of Fourier modes is used without clustering. In the following, this result serves as the reference solution. We performed calculations for different numbers of clusters for the matrix and inclusions, and , as well as for a varied percentage of used Fourier modes ℛ. Fig. 5 shows the norm of the difference of the stress component compared to the reference solution and given by The related strain-based sampling pattern of the last load step is presented in Fig. 6. For a more detailed analysis of the model order reduction technique, we kindly refer to [3].
The results in Fig. 5 show that as expected the increase of the number of clusters as well as the percentage of used Fourier modes lead to a decreased error in the stress distribution. But, already for a low number of clusters and Fourier modes the results show a good accordance to the reference fullfield solution. Fig. 7 allows a more detailed analysis of the computations. It shows the microscopic error (left) depending on the percentage of Fourier modes ℛ and calculated by In Fig. 7, the black lines refer to the method using a reduced set of modes but without clustering. The results underline the expectations and findings from the contour plot in Fig. 5: Also, the microscopic error decreases with a higher number of used clusters as well as of used Fourier modes. The error, which is only slightly higher than that one without clustering, shows a converging behavior. The simulations also show an improvement of the speed up factor, see Fig. 7 (right), which relates the

2290
Achievements and Trends in Material Forming computation time used for the individual simulation to the time needed for the reference solution.
Obviously, the decrease of used clusters not necessarily leads to a further reduction of computation times. The reason for that lies in the worse convergence of the strain field due to the clustered stress field within the microstructure. This motivates future investigations on that topic to improve the correspondence between the spectral solver and the clustered microstructure.

Summary
We presented a coupled approach to analyze microstructural evolutions in a highly resolved manner. The ansatz accounts for two different model order reductions techniques: a strain-based sampling pattern for the used fast Fourier transformations and a clustered microstructure for the computations in real space. The results show an improvement of the speed up factor that comes along with an acceptable increase of the microscopic error. In future works, we will further investigate the convergence behavior of the calculations also using more robust solvers like Newton-Krylov, see for instance [10,11,12]. The motivation for that is the high number of iterations that is necessary when a lower number of clusters is used. Additionally, the implementation within a finite element software will enable two-scale process simulations with a reduced computation time.