Fourier-based schemes for computing the mechanical response of composites with accurate local fields

François Willot

Source-linked AI summary

Fourier-based schemes for computing the mechanical response of composites with accurate local fields

François Willot

arXiv:1412.8398v1Dec 21, 2014math.NAcond-mat.mtrl-sci

TL;DR

The paper addresses discretization choices in Fourier-based elasticity schemes that affect local-field accuracy and convergence. It modifies the Green operator using centered differences on a rotated grid and finds more accurate local fields and faster convergence for stress equilibrium and effective elastic moduli. The approach also improves estimates of effective elastic moduli and reduces dependence on the reference medium.

Problem
The paper examines how discretization affects the accuracy of local stress and strain fields and convergence in Fourier-based elasticity schemes.
Method
The authors modify the Fourier-domain Green operator by deriving it from finite-difference discretizations of continuum mechanics, including centered differences on a rotated grid.
Results
The proposed discretization produces more accurate local fields and superior convergence rates for stress equilibrium and effective properties than the compared methods.
Takeaways & Limitations
The method provides better estimates of effective elastic moduli, and those estimates do not depend on the reference medium.
Takeaways & Limitations
The forward-and-backward discretization breaks axial symmetries, while the centered scheme has undefined Green operators at three frequencies for even grid sizes.

Abstract

from arXiv · show

We modify the Green operator involved in Fourier-based computational schemes in elasticity, in 2D and 3D. The new operator is derived by expressing continuum mechanics in terms of centered differences on a rotated grid. Use of the modified Green operator leads, in all systems investigated, to more accurate strain and stress fields than using the discretizations proposed by Moulinec and Suquet (1994) or Willot and Pellegrini (2008). Moreover, we compared the convergence rates of the "direct" and "accelerated" FFT schemes with the different discretizations. The discretization method proposed in this work allows for much faster FFT schemes with respect to two criteria: stress equilibrium and effective elastic moduli.

1. Introduction

FFT methods efficiently compute local mechanical responses of image-defined composites, including large 2D and 3D microstructures. This work focuses on how discretization affects FFT accuracy and convergence.

FFT methods compute local stress and strain tensors at each pixel or voxel of image-defined composite microstructures.
FFT tools have expanded beyond linear elasticity to problems including viscoplasticity and crack propagation.
They can process images containing billions of voxels and represent multiscale materials such as concrete or mortar.
Recent analysis connected the original Moulinec–Suquet method to a particular approximation space and optimization method, enabling other FFT schemes.
This work examines discretization effects on local fields and convergence rates in direct and accelerated FFT schemes.

2. Microstructure and material elastic response

The paper studies linear elasticity in square or cubic domains containing binary matrix–inclusion media under periodic overall-strain loading. Material behavior is parameterized mainly by the contrast between inclusion and matrix properties.

The elasticity problem is posed in a square or cubic domain Ω=[−1/2,1/2]^d with d=2 or 3.
The local response uses strain, stress, displacement, and elasticity-tensor fields, with isotropic elastic behavior in each phase.
The media are binary, with phase 1 designated as the matrix and phase 2 as the inclusions.
Poisson’s ratios are fixed at ν1=ν2=0.25, giving µα/κα=0.6 in both 2D and 3D.
The property contrast χ parametrizes the material, with χ=0 representing a porous medium and χ=∞ a rigidly reinforced one.
Periodic boundary conditions impose an overall strain, and the effective elastic tensor is computed from spatial averages of the resulting fields.

3. Lippmann-Schwinger equation and FFT methods

The FFT methods derive from Lippmann–Schwinger equations and use a reference elasticity tensor with an associated Green operator. Direct and accelerated schemes differ in iteration strategy, while discretization choices determine the Fourier representation and treatment of high-frequency modes.

Fourier methods follow the Lippmann–Schwinger equations derived from the elasticity and boundary-value formulation.
The Green operator is defined using a homogeneous reference tensor C0 and associated polarization field τ, with its Fourier form specified for nonzero wave vectors.
The direct scheme applies the Lippmann–Schwinger equations iteratively, whereas refined accelerated and augmented-Lagrangian methods address slow convergence in highly contrasted composites.
In the Moulinec–Suquet discretization, convolution is computed algebraically in Fourier space and fields are represented using trigonometric polynomials on an L^d voxel grid.
For even grid sizes, the standard Green operator fails a conjugate-symmetry relation at the highest frequency, producing a nonzero imaginary component in the inverse-transformed strain field.
Modified treatments can enforce zero stress or strain at selected highest-frequency modes, but numerical experiments found little convergence-rate influence except at small resolution.
The accelerated and direct convergence rates depend on the reference tensor C0, with different recommended choices for the two schemes.

4. Discretization and approximation space

The paper replaces the standard Fourier Green operator with discretization-specific operators derived from finite differences, including a new rotated-grid centered scheme. These operators enforce discrete compatibility and equilibrium while addressing high-frequency modes and extending the construction to 3D.

Modified Green operators: The paper derives a general modified Green operator G′ that includes earlier discretizations and introduces a new operator.The resulting operators are denoted GC, GW, and, for the rotated scheme, GR.
Finite-difference formulation: Finite-difference discretizations approximate equilibrium and strain admissibility directly on the voxel grid rather than assuming a continuum trigonometric representation.The displacement, strain, and stress fields are placed at different grid locations according to each scheme.
Two-dimensional schemes: The centered and forward-and-backward schemes use different discrete gradients, with the latter estimating derivatives more locally near interfaces but breaking the symmetry of the two angle bisectors.The rotated scheme instead evaluates displacement at pixel corners and strain and stress at pixel centers using centered differences in a 45°-rotated basis.
High-frequency modes: When the grid size L is even, centered and rotated operators become undefined at modes where their discrete gradient vanishes, so the operators are set to zero there.These treatments enforce admissibility or compatibility and stress equilibrium in the corresponding discrete formulation; GW needs no analogous high-frequency treatment.
Two-dimensional schemes: For the rotated scheme, GR is obtained by substituting the rotated discrete gradient kR into the Green-operator expression.In 2D, the construction uses a rotated basis and places displacement at pixel corners while evaluating strain and stress at centers.
Three-dimensional extension: The 3D construction generalizes the rotated scheme by evaluating strain and stress at voxel centers, displacement at corners, and derivatives across opposite corners.The corresponding operators GC, GW, and GR are derived in 3D, with special zero assignments for GC and GR at even-grid modes where the discrete gradients vanish.
Approximation properties: The finite-difference operators are periodic with high frequencies cut, and GR is expected to provide higher local-field accuracy because it uses centered differences and more appropriate derivative locations.The operators approach the continuum Green operator as q → 0, preserving discretization-independent fields in the fine-resolution limit.

5. Local strain and stress fields accuracy

Across 2D and 3D benchmarks, the modified operator GR produces local stress fields with fewer oscillations and better symmetry than the other discretizations, while effective-modulus estimates converge toward a common limit.

Two-dimensional benchmark: At L = 2048, all methods approach the same local stress field, but G still produces spurious interface oscillations up to 2048^2 pixels.The oscillations persist after local averaging.
Two-dimensional benchmark: GR greatly reduces oscillations relative to G and GC and preserves the problem’s symmetries, unlike GW near inclusion corners.GC produces checkerboard and aligned patterns, while GW remains nonsymmetric along one-pixel-wide lines.
Three-dimensional benchmark: In 3D, G and GC generate strong vertical and horizontal artificial patterns at L = 512, whereas GW and GR have the smallest oscillations.GW remains asymmetric near one corner, while GR produces symmetric fields with almost no oscillations at all tested resolutions.
Periodic array of spheres: The spherical-inclusion effective modulus e C11,11 increases with resolution and approaches 1.208±0.001 for all schemes.The estimate uses resolutions L = 32, 64, 128, 256 and 512, with iterations stopped at a maximum stress variation below 2 10^-10.
Periodic array of spheres: At fixed resolution, G has about twice the prediction error of GR, while GC and GW lie between them.GR provides the best estimate among the methods compared.

6. Convergence rate

The study evaluates direct and accelerated FFT convergence using different Green operators, optimizing the reference modulus and examining stress-equilibrium and effective-modulus criteria. The proposed discretization shows faster convergence in low-contrast cases and substantially improves effective-modulus estimates for quasi-porous media.

6. Convergence rate: Convergence is assessed with an L2-norm criterion based on stress equilibrium, using precision η ≪1 and the Frobenius norm of the mean stress for normalization.All schemes enforce stress equilibrium only at convergence.
6.1. Convergence rate with respect to stress equilibrium: For χ = 10^-5, N(E0) is similar across accelerated schemes below E0 ≲0.03, but rises strongly with E0 for AS while remaining less sensitive for ASC,W,R.The ASC,W,R schemes exhibit local minima; ASR has one near E0 ≈0.09, whereas ASW and ASC have two.
6.1. Convergence rate with respect to stress equilibrium: Poisson-ratio deviations strongly increase N(E0), so the reference Poisson ratio is fixed at ν0 = 0.25 and E0 is optimized by gradient descent.The descent compares iteration counts and, when tied, the achieved precision.
6.2. Convergence rate with respect to the effective elastic moduli: For direct schemes, the optimal reference follows E0 proportional to χ^β with 0.5003 ≤ β ≤ 0.509, while accelerated-scheme optima depend on the Green operator and can approach small constants at low contrast.For AS, E0 is about 0.01 when χ ≤10^-3; ASC and ASR show constants around 0.07 and 0.12, respectively.
6.1. Convergence rate with respect to stress equilibrium: The direct-scheme iteration count scales as χ for χ ≫1 and 1/χ for χ ≪1, whereas AS scales as √χ and 1/√χ, respectively.The accelerated scheme therefore requires fewer iterations, apart from one reported exception.
6.1. Convergence rate with respect to stress equilibrium: For χ < 1, DSC,W,R and ASC,W,R require at most 430 iterations, while ASR requires at most 168 and is the fastest accelerated scheme.Iteration counts are nearly constant over 0 ≤χ ≤10^-5.
6.2. Convergence rate with respect to the effective elastic moduli: After about 7 iterations, ASR estimates the quasi-porous effective modulus to relative precision 10^-2, compared with more than 50 iterations for AS and ASC,W.AS and ASC,W display strong oscillations that ASR substantially reduces.
6.2. Convergence rate with respect to the effective elastic moduli: For quasi-rigid spheres, all schemes require many more iterations to reach precision 10^-2, and ASR is less advantageous than in the quasi-porous case.The slower convergence follows the poorer behavior observed for χ > 1.

7. Conclusion

The proposed discretization improves local-field accuracy, especially near interfaces, and provides better effective-modulus estimates than other methods. Accelerated schemes are compared through elastic-modulus estimates over iterations for a 3D quasi-rigid-sphere model.

Figure 10 tracks the elastic modulus e C11,11 against iteration count for accelerated schemes AS, ASC, ASW, and ASR.The comparison uses a 3D Boolean model of quasi-rigid spheres.
The new discretization produces more accurate local fields than other methods, especially near interfaces.The method is proposed for Fourier-based schemes in both 2D and 3D.
The method provides better estimates for effective elastic moduli.