Skip to content

API — Elasticity ​

MeanFieldHomogenization.Elasticity Module
julia
MeanFieldHomogenization.Elasticity

Hill polarization tensors for ellipsoidal inclusions (2D / 3D, isotropic / anisotropic matrix). Public entry point: hill_tensor.

MeanFieldHomogenization.Elasticity.EllipsoidShape Type
julia
EllipsoidShape

Abstract supertype for the shape classification of an Ellipsoid.

Concrete subtypes (3-D): Spherical, Prolate, Oblate, Triaxial.

Concrete subtypes (2-D): Circular, Elliptic.

MeanFieldHomogenization.Elasticity.Spherical Type

3-D sphere: a = b = c.

MeanFieldHomogenization.Elasticity.Prolate Type

3-D prolate spheroid: a > b = c (axis of revolution = e₁).

MeanFieldHomogenization.Elasticity.Oblate Type

3-D oblate spheroid: a = b > c (axis of revolution = e₃).

MeanFieldHomogenization.Elasticity.Triaxial Type

3-D triaxial ellipsoid: a > b > c.

MeanFieldHomogenization.Elasticity.Circular Type

2-D circle: a = b.

MeanFieldHomogenization.Elasticity.Elliptic Type

2-D ellipse: a > b.

MeanFieldHomogenization.Elasticity.Ellipsoid Type
julia
Ellipsoid{dim, S<:EllipsoidShape, T<:Number, B<:AbstractBasis}

Ellipsoidal inclusion of the Eshelby problem (Eshelby, 1957). In the Echoes convention, an ellipsoid is described by an invertible second-order shape tensor through

semi_axes stores the eigenvalues (sorted in decreasing order for real-valued types) and basis stores the orthonormal frame relative to the canonical frame.

The shape S is determined at construction time:

  • 3-D: Spherical, Prolate, Oblate, or Triaxial

  • 2-D: Circular or Elliptic

T can be any Number subtype (Float64, ForwardDiff.Dual, SymPyPythonCall.Sym, Symbolics.Num, …).

MeanFieldHomogenization.Elasticity.Spheroid Function
julia
Spheroid(ω; euler_angles = ())

Convenience constructor for an axisymmetric ellipsoid (spheroid). The aspect ratio is the ratio of the polar semi-axis (along the axis of revolution) to the equatorial semi-axis:

  •   ⇒ oblate (disc-like, axis of revolution is the short axis)

  •   ⇒ prolate (needle-like, axis of revolution is the long axis)

  •   ⇒ sphere

The two equatorial semi-axes are fixed to . Eshelby/Hill computations are scale-invariant, so the absolute size of the inclusion has no effect on the concentration tensor and only matters.

Optional euler_angles (tuple of length 0–3, ZYZ convention, radians) rotate the symmetry axis from . Without angles, the axis of revolution is (oblate) or (prolate) — i.e. the result is strictly equivalent to Ellipsoid(1, 1, ω; euler_angles), with internal sorting and basis permutation handled by the regular Ellipsoid constructor.

Examples

julia
Spheroid(0.2)                                      # oblate (axis ‖ ez)
Spheroid(5.0)                                      # prolate (axis ‖ ez after sort)
Spheroid(0.2; euler_angles = (π/4, π/3, 0))        # tilted axis
MeanFieldHomogenization.Elasticity.Cylinder Type
julia
Cylinder{S<:CylindricalShape, T<:Number, B<:AbstractBasis} <:
    AbstractEllipsoidalInclusion{3, T}

Infinite cylindrical inclusion with transverse semi-axes (with    when T <: Real, or in the caller-provided order when T is symbolic). The cylinder axis is the first column of the local basis — consistent with the Prolate convention where the axis of revolution is also .

S encodes the cross-section shape:

T can be any Number subtype (Float64, ForwardDiff.Dual, SymPyPythonCall.Sym, Symbolics.Num, …) — all analytical paths are type-generic.

MeanFieldHomogenization.Elasticity.CylindricalShape Type
julia
CylindricalShape

Abstract supertype for the shape classification of a Cylinder. Sub-traits:

All cylinders are infinite along in the inclusion's local basis.

MeanFieldHomogenization.Elasticity.CircularCylindrical Type

3-D infinite cylinder with circular base: b = c, a → ∞, axis = e₁.

MeanFieldHomogenization.Elasticity.EllipticCylindrical Type

3-D infinite cylinder with elliptic base: b > c, a → ∞, axis = e₁.

MeanFieldHomogenization.Elasticity.tens_IA Function
julia
tens_IA(ell::Ellipsoid{3}) -> AbstractTens{2,3}
tens_IA(ell::Ellipsoid{2}) -> AbstractTens{2,2}

Newton-potential geometric tensor of an ellipsoid, defined by

is symmetric with the same eigenvectors as and with diagonal components satisfying   (Kellogg, 1929; Eshelby, 1957; Parnell, 2016).

MeanFieldHomogenization.Elasticity.tens_UA Function
julia
tens_UA(ell::Ellipsoid{3}) -> AbstractTens{4,3}
tens_UA(ell::Ellipsoid{2}) -> AbstractTens{4,2}

4th-order Newton-potential geometric tensor :

In the principal frame, the non-zero Kelvin–Mandel components are

with the and coefficients given in the Echoes appendix (tables / ).

In two dimensions the weight is    on the unit circle, the one of tens_IA and tens_VA, and an ellipse of semi-axes has, with   ,

so that the Hill tensor of an isotropic medium is     in both dimensions. These closed forms hold at  .

MeanFieldHomogenization.Elasticity.tens_VA Function
julia
tens_VA(ell::Ellipsoid{3}) -> AbstractTens{4,3}
tens_VA(ell::Ellipsoid{2}) -> AbstractTens{4,2}

4th-order Newton-potential geometric tensor :

In the principal frame,   and     for  .

MeanFieldHomogenization.Elasticity.hill_tensor Function
julia
hill_tensor(ell, C₀; method=:auto, abstol=1e-8, reltol=1e-6, maxiters=1_000_000)
    → AbstractTens

Hill polarization tensor ( in conduction) for an ellipsoidal inclusion ell embedded in a reference medium C₀. C₀ can be a 4th-order stiffness (elasticity) or a 2nd-order conductivity tensor — dispatch selects the appropriate formulation automatically.

The general expression of the elastic polarization tensor is (Willis, 1977; Mura, 1987):

The isotropic case (C₀::TensISO) is evaluated analytically; the anisotropic case uses the Cauchy-residue reduction of (Masson, 2008) (trait Residue) or the DECUHR adaptive cubature of Espelid and Genz (1994) (trait DECUHR). See the Hill polarization tensors theory page for the full dispatch table and return types.

julia
hill_tensor(incl::NeuralHillInclusion, P₀; kw...) -> AbstractTens

Hill tensor evaluated by the surrogate: features off the geometry and P₀, one forward pass, then an exact decode into the structured tensor of the surrogate's class, in the global frame.

Keyword arguments of the analytic entry points (method, tolerances) are accepted and ignored — there is no quadrature to steer — so a scheme that forwards them works unchanged.

MeanFieldHomogenization.Core.HillCache Type
julia
HillCache(; maxentries = 1024)

A store of Hill tensors for with_hill_cache: each numerical Hill tensor (anisotropic matrix: residues, nested Gauss–Kronrod, DECUHR) is computed once per shape, frame, matrix, method and tolerances, and read back afterwards. When it holds maxentries tensors it is emptied before storing the next one. hits and misses count the reads and the computations. Thread-safe.

MeanFieldHomogenization.Core.with_hill_cache Function
julia
with_hill_cache(f, cache::HillCache)

Run f() with the numerical Hill tensors read from cache when it holds them, and stored in it otherwise. For a series of calls that keeps the shapes and the anisotropic matrix:

julia
cache = HillCache()
estimates = with_hill_cache(cache) do
    [homogenize(cell(f), MoriTanaka(), :C) for f in 0.01:0.01:0.3]
end

Each tensor is keyed by the exact values that determine it (the shape, its semi-axes and frame, the canonical components of the matrix, the method and the tolerances), and only Float64 ones are stored: a dual or a symbolic input is computed as usual, so that derivatives are unaffected. The results are the bits the calls give without the cache.

MeanFieldHomogenization.Core.eshelby_tensor Function
julia
eshelby_tensor(incl, C₀; method=:auto, abstol, reltol, maxiters) -> AbstractTens

Eshelby tensor of the inclusion incl embedded in a matrix of stiffness / conductivity C₀, derived from the Hill polarization tensor (or ) by the relations

in elasticity (order 4), and

in conductivity / diffusion (order 2).

The appropriate method is selected by dispatch on the order of C₀: an AbstractTens{4, 3} (elasticity) triggers the double contraction  , while an AbstractTens{2, 3} (conductivity) triggers the simple contraction  .

All keyword arguments (method, abstol, reltol, maxiters) are forwarded verbatim to hill_tensor; see its docstring for the set of admissible algorithm traits.

See also hill_tensor.

MeanFieldHomogenization.Core.second_hill_tensor Function
julia
second_hill_tensor(incl, C₀; method=:auto, abstol, reltol, maxiters) -> AbstractTens

Second Hill tensor of the inclusion incl in a matrix of stiffness or conductivity C₀:

in elasticity (order 4), and

in conductivity / diffusion (order 2), with () the Hill tensor. It gives the stress (flux) inside the inclusion under an eigenstrain (eigengradient),    , and, unlike , has a finite limit when the ellipsoid flattens into a crack, whose compliance is  ; see The Eshelby inclusion problem.

The order of C₀ selects the formula, and the structure of the operands is kept: an isotropic or a coaxial transversely isotropic pair gives a structured . The keyword arguments are forwarded to hill_tensor. For an aging-viscoelastic matrix, see second_hill_kernel and second_hill_kernel_order2.

See also eshelby_tensor.

Nanoinclusion interfaces ​

MeanFieldHomogenization.Elasticity.surface_stiffness Function
julia
surface_stiffness(spheroid, κs, μs) -> TensTI{4}

Average surface stiffness of a spheroidal nanoinclusion carrying a Gurtin-Murdoch interface of surface bulk modulus (κs) and surface shear modulus (μs) (Dormieux et al., 2016).

spheroid is an Ellipsoid of revolution with in-plane semi-axis and polar semi-axis ; the aspect ratio is  , oblate for   and prolate for  . The result is transversely isotropic about the symmetry axis, with

and four companion components given in the same reference. κs and μs have the dimension of a stiffness times a length, so scales as : the stiffening it produces is a size effect, controlled by the smallest dimension of the particle.

Adding it to the bulk stiffness gives the equivalent particle of the paper, equivalent_particle, which any scheme of the package accepts as an ordinary inclusion property.

Limiting cases are handled exactly: the spherical case   — where the closed form has a removable singularity — is evaluated from its Taylor series, and reproduces the isotropic tensor      of the same reference.

MeanFieldHomogenization.Elasticity.equivalent_particle Function
julia
equivalent_particle(C_I, spheroid, κs, μs) -> AbstractTens

Stiffness of the equivalent particle — the nanoinclusion C_I together with its interface — of Dormieux et al. (2016):

The paper's central result is that the strain concentration rule and the homogenized stiffness keep their classical form once this substitution is made, so a nanocomposite estimate needs no new scheme:

julia
C_eq = equivalent_particle(C_i, Spheroid(0.1), κs, μs)
add_phase!(rve, :nano, Spheroid(0.1), Dict(:C => C_eq); fraction = f)
homogenize(rve, MoriTanaka(), :C)

Note that the strain average rule is unaffected by the interface, but the stress average rule is not: the interface contributes   to the macroscopic stress, which is precisely what using in the scheme accounts for.

Parameter conversions ​

The isotropic moduli are TensND's, re-exported by the package unchanged: iso_stiffness(k, μ), k_mu(ℂ), E_nu(ℂ) and iso_stiffness_E_nu(E, ν) (TensND, isotropic elastic moduli). The Hoenig parameters of a transversely isotropic stiffness are read and built here, on TensND's ti_constants and ti_from_constants with form = :hoenig, which act on the compliance.

MeanFieldHomogenization.Elasticity.hoenig_params Function
julia
hoenig_params(t::TensND.TensTI{4}) -> (E1, h, nu1, nu2, gamma)
hoenig_params(t::TensND.AbstractTens{4,3}, axis) -> (E1, h, nu1, nu2, gamma)

Hoenig (1978) parametrization of a transversely isotropic stiffness tensor: E1 the in-plane Young's modulus, h the axial/in-plane modulus ratio, nu1 the in-plane Poisson's ratio, nu2 the out-of-plane Poisson's ratio, and gamma the shear anisotropy ratio: TensND's TensND.ti_constants of the compliance inv(t).

For a tensor that is not already TensTI, the 2-argument form projects onto the TI span about axis first (TensND.best_fit_ti). For a compliance tensor, use hoenig_params(inv(S)) / hoenig_params(inv(S), axis).

MeanFieldHomogenization.Elasticity.hoenig_stiffness Function
julia
hoenig_stiffness(E1, h, nu1, nu2, gamma, axis) -> TensND.TensTI{4,T,5}

Build the TI stiffness tensor from its Hoenig (1978) parameters and symmetry axis — the reciprocal of hoenig_params: the inverse of TensND's TensND.ti_from_constants compliance.

The cubic symmetry class ​

The class itself lives in TensND: TensCubic stores the three constants and a cube frame, its products and inverses are componentwise, cubic_from_constants / cubic_constants convert to and from , cubic_anisotropy is the Zener-type departure from isotropy, and project(t, Val(:CUBIC), frame) — or best_fit_cubic — projects onto the class. None of that is reimplemented here.

What this package adds is a name for the projection residual, which is the quantity a solver reads, and the symmetry trait. The trait is honest only now that a storage type exists: material_symmetry answers from the container, and before TensND had TensCubic a CubicSym would have let dispatch claim a structure the object did not carry.

Read cubic_residual and cubic_anisotropy together. The first is the distance to a class the answer belongs to by group theory, so it is discretization error and nothing else; the second lives inside the class. An artifact breaks the symmetry, a real morphological anisotropy does not.

MeanFieldHomogenization.Elasticity.cubic_residual Function
julia
cubic_residual(t, frame = CanonicalBasis{3,Float64}()) -> Real

Relative Frobenius distance from t to its projection onto the cubic class with cube axes frame — the third value returned by project(t, Val(:CUBIC), frame), named because it is the one that gets read.

This is a free error bar. Where the morphology and the medium leave the octahedral group invariant, the answer belongs to the class by group theory, so this distance is bounded by the discretization error and by nothing else: no reference solution appears in it. What lives inside the class, and is therefore not measured here, is TensND.cubic_anisotropy. Read together they separate a real morphological anisotropy from a numerical artifact — an artifact breaks the symmetry, a real anisotropy does not.

The default frame is the canonical one. Pass the inclusion's own basis when the cube is not aligned with it: a rotated cubic tensor is not cubic about the canonical axes, and the residual would say so, correctly but uselessly.

MeanFieldHomogenization.Elasticity.CubicSym Type

Cubic material (TensND TensCubic) — three independent constants.