Homogenization schemes

MeanFieldHom.Schemes computes the effective property tensor of a multi-phase medium from (i) the phase geometries, (ii) the phase properties and (iii) the phase volume fractions or crack densities.

Notation

The Representative Volume Element (RVE) consists of:

  • a matrix phase of property tensor $\mathbb C_0$ (or $\mathbf K_0$ for the 2nd-order conductivity problem),
  • one or more inclusion phases of property tensors $\mathbb C_i$, geometries $\mathcal G_i$, and amounts $f_i$ (volume fraction) or $\varepsilon_i$ (crack density).

For each inclusion the dilute strain concentration tensor $\mathbb A_\mathrm{dil}^{(i)}$ and the size-independent stiffness contribution $\mathbb N_i = (\mathbb C_i - \mathbb C_0): \mathbb A_\mathrm{dil}^{(i)}$ are the natural building blocks (Kachanov & Sevostianov 2018). The dual compliance contribution $\mathbb H_i = (\mathbb S_i - \mathbb S_0): \mathbb A_\sigma^{(i)}$ is more natural for cracks (whose stiffness contribution is the rank-1 limit of a divergent eigenvalue).

Bounds

SchemeFormula
Voigt$\langle \mathbb C \rangle = \sum_i f_i \mathbb C_i$ (upper bound, Hill 1965)
Reuss$\langle \mathbb S \rangle^{-1}$ (lower bound)

Cracks are ignored in both bounds: their volume contribution vanishes in the penny limit (c → 0) while their density stays finite.

One-shot schemes (require a matrix)

SchemeEffective stiffness
Dilute$\mathbb C_0 + \sum_i f_i \mathbb N_i$ (first order in f)
DiluteDual$\big(\mathbb S_0 + \sum_i f_i \mathbb H_i\big)^{-1}$
Mori-Tanaka$\mathbb C_0 + \big(\sum_i f_i \mathbb N_i\big) : \big(f_m\,\mathbb I + \sum_i f_i \mathbb A_\mathrm{dil}^{(i)}\big)^{-1}$ (Mori-Tanaka 1973, Christensen 1990)
Maxwell$\mathbb C_0 + \Sigma : (\mathbb I - \mathbb P_d : \Sigma)^{-1}$ with P_d the Hill tensor of the outer distribution shape
PCWidentical algebraic form, distribution-shape-aware ensemble interpretation (Ponte-Castañeda & Willis 1995)

The distribution shape is stored at the RVE level (default: unit sphere ⇒ Mori-Tanaka limit). Any AbstractInclusion can be used; the hierarchy AbstractDistributionShape leaves room for a future PairwiseDistribution extension following Willis 1982.

Iterative schemes

SchemeIteration
SelfConsistent (McLaughlin 1977)$\mathbb C^{(n+1)} = \big(\sum_i f_i \mathbb C_i : \mathbb A_\mathrm{dil}^{(i)}(\mathbb C^{(n)})\big) : \big(\sum_i f_i \mathbb A_\mathrm{dil}^{(i)}(\mathbb C^{(n)})\big)^{-1}$
AsymmetricSelfConsistentswitches between stiffness- and compliance-form iteration based on the matrix-vs-Voigt-bound contrast

The default solver is a damped Picard fixed point (Anderson with memory 1, Dual-safe). Loading NonlinearSolve.jl activates the MeanFieldHomNonlinearSolveExt extension, which accepts every SciML non-linear algorithm (NewtonRaphson(), TrustRegion(), Anderson(), …) via the algorithm keyword of SelfConsistent.

Differential scheme

The DifferentialScheme integrates the multi-phase Norris ODE (Norris 1985) on a fictitious incorporation time $\tau \in [0, 1]$,

\[\frac{\mathrm d \mathbb C^{hom}}{\mathrm d \tau} = \sum_i \dot\varphi_i \, \mathbb N_i(\mathbb C^{hom}) , \qquad \dot\varphi_i = \dot f_i + \frac{f_i}{f_0} \sum_j \dot f_j ,\]

the increments $\dot\varphi_i$ following from the volume balance by Sherman-Morrison, along a user-selectable trajectory (Proportional, Sequential, CustomPath, Path). The dual form on the compliance is available through formulation = :compliance, and cracks — which have no volume but a finite density — enter with a balance of their own.

The trajectories agree in the dilute limit (f → 0) and diverge like f at finite fractions — a physical feature of the scheme.

The full derivation, the crack case, the closed form of the homothetic trajectory and the SciML resolution are in The differential scheme.

Number-type compatibility

Every scheme is mandated to support:

  • Float64 — default;
  • ForwardDiff.Dual — sensitivity analysis through fractions, moduli, geometric parameters;
  • Complex{Float64} — frequency-domain viscoelasticity (parity with the C++ ECHOES library, templated on T = double | complex<double>);
  • SymPy.Sym, Symbolics.Num, BigFloat — best-effort, with explicit documentation of any limitation (the iterative SC solvers are not symbolic-friendly because the linear-system Jacobian must be numeric).