Skip to content

Aging linear Volterra behavior (ALV) ​

Before this page

The Laplace-Carson route, which treats the non-aging case contrasted below, and the elastic chain transposed on this page: The Eshelby inclusion problem, Concentration and contribution tensors and Homogenization schemes.

If the material does not age, there is a cheaper route

Everything on this page treats a kernel in which and enter independently, which is what forces the Volterra operators to be discretized directly. When the kernel depends on t - t' alone the constitutive law is a convolution, a transform turns it into a product, and the whole problem reduces to an elastic one at each Carson variable — see the Laplace-Carson route.

The two algebras are worth contrasting explicitly. Creep and relaxation are related by     here, with the relaxation kernel and the creep kernel — a convolution, discretized below into the inverse of a block-triangular matrix — and by the pointwise reciprocal    there. That single difference is the whole cost and the whole benefit of each route.

ALV stands for aging linear Volterra: a linear behavior whose response at time depends on the whole history through a two-time kernel, the kernel being free to age. In mechanics it is aging linear viscoelasticity, the subject of most of this page. The same algebra holds at order 2 for a hereditary conductivity, diffusivity or permittivity, whose kernel maps the history of the gradient to the flux analog    of Conventions; every statement below transposes to it by the dictionary of that page.

Everything in the elastic part of this documentation — the Eshelby problem, the Hill polarization tensor, the concentration and contribution tensors, the schemes built on them — carries over to aging linear viscoelasticity essentially unchanged, provided two substitutions are made systematically:

ElasticAging viscoelastic
tensor two-time kernel
double contraction  Volterra product  
tensor inverse Volterra inverse
identity  , the Heaviside function

This page follows that substitution from the constitutive law to the homogenization schemes. The derivations are those of Barthélémy et al. (2016) (the Eshelby problem and the Hill kernel) and Barthélémy et al. (2019) (the schemes); the time discretization is that of Sanahuja (2013).

Notation used on this page

A kernel is a function of two times — observation time and loading time — and is written in the same typeface as its elastic counterpart. denotes the Volterra product defined in the next section, the Volterra inverse — never the ordinary tensor inverse. In MeanFieldHomogenization these two operations are volterra_product and volterra_inverse.

1. The aging linear viscoelastic behavior ​

The strain and stress histories are related by a Stieltjes integral (Barthélémy et al., 2016):

with the creep compliance kernel and the relaxation kernel. Causality imposes   for  .

The non-aging case is the special one where the kernels depend on and only through their difference  . There, and only there, the Laplace–Carson correspondence principle applies and the problem reduces to an elastic one with complex moduli — the route taken in Viscoelastic composites and cross-checked against the present one in Frequency or time?. When a phase ages — its properties evolve with its own maturity, as a hydrating cement paste does — and enter independently, there is no convolution to transform, and the time domain is the only available route.

Volterra algebra ​

The relation above is written compactly   , extending to tensors the scalar Volterra operator (Barthélémy et al., 2016). Between two kernels the same symbol denotes

This product is associative and distributive over addition, but not commutative: commutativity holds only for non-aging kernels (Barthélémy et al., 2016; after Maghous and Creus). Every formula below therefore keeps its factors in order, including the apparently scalar ones.

The Heaviside function acts as the identity,   , so the identity elements are   for 4-tensor kernels and   for 2-tensor kernels. The Volterra inverse is defined by      ; relaxation and creep kernels are Volterra inverses of one another,  .

Discretization: kernels become block matrices ​

On a time grid     the Stieltjes integral is approximated by the trapezoidal rule of Sanahuja (2013). Strain and stress histories become block column vectors, and each kernel a lower block-triangular matrix of size     whose blocks are

A block is not a kernel value

is a weighted combination of kernel values at neighboring times, not . Reading a relaxation function off a column of the matrix therefore gives the wrong answer; the physical extraction is to apply the matrix to a unit strain step, which for a step at amounts to summing each block row. The lower-triangular structure encodes causality.

The pay-off of this representation is that the Volterra product becomes an ordinary matrix product and the Volterra inverse an ordinary matrix inverse (Sanahuja, 2013), as in the appendix of Barthélémy et al. (2019). In MeanFieldHomogenization the discretization is trapezoidal_matrix and the inverse volterra_inverse.

2. The Eshelby problem in ALV ​

Consider an ellipsoid of shape tensor embedded in an infinite medium of relaxation kernel , carrying a uniform polarization history :

with the characteristic function of . Momentum balance and the decay condition at infinity give (Barthélémy et al., 2016)

where is the surface Dirac distribution on the boundary of and its outward normal. This is literally the elastic problem with replaced by .

Solving it through the ALV Green kernel yields the central result: the strain is uniform inside the ellipsoid, exactly as in elasticity, and

3. The Hill polarization kernel ​

The kernel appearing above is the ALV counterpart of the elastic Hill polarization tensor (Barthélémy et al., 2016):

It differs from the elastic Hill tensor in exactly one place: the inverse of the acoustic tensor    is a Volterra inverse. Like its elastic counterpart it depends only on the shape and orientation of the ellipsoid and on the reference kernel.

The Eshelby kernel follows by the same definition as in elasticity,   , and relates the uniform strain inside to a uniform eigenstrain history :

Isotropic matrix: time and space decouple ​

If the reference kernel is isotropic,     , the acoustic tensor is diagonal in the ``(\underline{\xi}\otimes\underline{\xi}, MarkdownAST.LineBreak()

\boldsymbol{1}-\underline{\xi}\otimes\underline{\xi})`` decomposition and can be inverted in the Volterra sense analytically. The Hill kernel then factorizes (Barthélémy et al., 2016; Barthélémy et al., 2019):

where and are the purely geometric tensors of the elastic theory (Hill polarization tensors):

This is what makes ALV tractable: the geometry sits entirely in and — the elastic ones, computed once per shape — and the time dependence in two scalar Volterra inverses, so a   tensor-kernel inversion becomes two   scalar ones. hill_kernel discretizes the matrix law, extracts , inverts   and in the Volterra sense, and assembles against tens_UA and tens_VA. An anisotropic reference kernel has no such fast path.

Why the longitudinal modulus, and not the bulk one. The Hill tensor inverts the acoustic tensor of the reference. For an isotropic one, ``\underline{\xi}\cdot\mathbb{C}_0\cdot\underline{\xi} = \bigl(k + \tfrac{4}{3}\mu\bigr),\underline{\xi}\otimes\underline{\xi}

  • \mu,\bigl(\boldsymbol{1} - \underline{\xi}\otimes\underline{\xi}\bigr)``

for a unit : its two eigenvalues are the longitudinal (P-wave) modulus   and the shear modulus . These two, not the bulk modulus alone, are the scalars the formula above inverts, which is why the code works with the pair

with   the coefficients of and , and with their Volterra inverses   and  , so that ``\widetilde{\mathbb{P}} = J_{\mathrm{long}},\mathbb{U}^{\boldsymbol{A}}

  • J_{\mathrm{shear}},\bigl(\mathbb{V}^{\boldsymbol{A}} - \mathbb{U}^{\boldsymbol{A}}\bigr)``,

block by block in time.

Going back to bulk and shear is linear on the relaxation kernels:

It is not linear on their inverses: is the Volterra inverse of the longitudinal kernel, not a combination of the bulk and shear creep functions. Those are and , one Volterra inversion each of the kernels above; in the elastic limit    .

The example below follows these steps on a Maxwell law discretized on five instants. trapezoidal_matrix gives the   block matrix of the relaxation kernel, iso_params_from_blocks its   kernels and , and volterra_inverse the discrete Volterra inverse of a scalar kernel (block_size = 1).

julia
using MeanFieldHomogenization, TensND

law = maxwell_iso(10.0, 4.0, 1.0, 0.5)   # k, μ and their relaxation times
times = collect(range(0.0, 2.0; length = 5))
α, β = iso_params_from_blocks(trapezoidal_matrix(law, times))   # 3k and 2μ

# longitudinal and shear kernels, the pair the Hill kernel inverts
M_long, M_shear = (α .+ 2 .* β) ./ 3, β ./ 2

# back to bulk and shear: linear on the relaxation kernels
k, μ = M_long .- (4 / 3) .* M_shear, M_shear
(k ≈ α ./ 3, μ ≈ β ./ 2)
(true, true)

The inverses do not convert this way: the longitudinal creep kernel differs from the bulk one, already at the first instant, where both reduce to the elastic values   and .

julia
J_long = volterra_inverse(M_long; block_size = 1)
J_shear = volterra_inverse(M_shear; block_size = 1)
J_bulk = volterra_inverse(k; block_size = 1)
(J_long[1, 1], J_bulk[1, 1], 1 / (10.0 + 4 / 3 * 4.0), 1 / 10.0)
(0.06521739130434782, 0.1, 0.06521739130434782, 0.1)

Assembling ``J_{\mathrm{long}},\mathbb{U}^{\boldsymbol{A}}

  • J_{\mathrm{shear}},(\mathbb{V}^{\boldsymbol{A}} - \mathbb{U}^{\boldsymbol{A}})``

by hand, in Kelvin–Mandel components, gives the Hill kernel that hill_kernel returns:

julia
ell = Ellipsoid(3.0, 1.0, 1.0)
U = TensND.mandel(TensND.canonical_components(tens_UA(ell)); symmetrize = true)
D = TensND.mandel(TensND.canonical_components(tens_VA(ell)); symmetrize = true) - U
kron(J_long, U) + kron(J_shear, D) ≈ hill_kernel(ell, law, times)
true

For a sphere the geometric tensors are    and  , giving the closed form (Barthélémy et al., 2016)

An anisotropic reference kernel does not enjoy this decoupling; as in elasticity, that is the case where the surface integral must be evaluated numerically.

4. From the inclusion to the inhomogeneity ​

Replace the polarization by a genuine inhomogeneity: the ellipsoid now has its own relaxation kernel , and the medium is loaded by a remote strain history . Superposing the remote field and the response to the fictitious polarization    gives (Barthélémy et al., 2016)

with the dilute strain concentration kernel

The strain remains uniform inside : it depends on time alone. The associated contribution kernel is

the second form following from the identity        , which holds in any associative algebra and so survives the loss of commutativity.

5. Schemes ​

With concentration kernels in hand, every matrix-based scheme transposes term by term (Barthélémy et al., 2019). Writing for the volume fraction of phase , for its relaxation kernel and for the reference kernel, the general form is   , or equivalently

and the schemes differ only in how is estimated. Each one is the elastic scheme of the chapter on schemes, with every product a Volterra product and every inverse a Volterra inverse. With and the contribution and dilute concentration kernels of phase ,   , the creep kernel of the matrix and the Heaviside kernel, the identity of the Volterra algebra:

SchemeEffective kernel
Voigt, Reuss  and  
Dilute 
DiluteDual  , with      
Mori–Tanaka      
Maxwell, PCW     
Self-consistent    , the kernels computed in the estimate itself
Asymmetric self-consistent    , the kernels computed in the estimate itself
Differentialthe incorporation ODE, integrated on the kernels

in the Maxwell row is the Hill kernel of the distribution shape , not of an inclusion. When is invertible the same estimate reads    (Barthélémy et al., 2019), by the identity of §4: the far-field argument of the elastic derivation carries over unchanged, since it uses nothing but the linearity of the integral equation. The two self-consistent schemes and the differential one take the running estimate as reference, so that their Hill kernels, built for an isotropic reference, require every phase to keep that estimate isotropic.

The one new difficulty is bookkeeping: does not commute, so the order of the factors is prescribed — even in the isotropic case, where every factor looks scalar. Two kernels commute only on a uniform time grid with non-aging laws, or when one of them is elastic; an aging law, or a logarithmic grid, makes the order matter at the percent level. In the Mori–Tanaka row the inverse stands on the right of  , because it maps the average strain to that of the matrix, on which the contributions act. All ten schemes are implemented by homogenize_alv at order 4; at order 2 all of them but the asymmetric self-consistent one.

6. The n-layer composite sphere ​

The -layer sphere of Hervé and Zaoui (1993) (see Layered spheres) transposes by the same rule. Its elastic construction propagates a state vector across the shells by a product of transfer matrices —   for the bulk () harmonic,   for the shear () one — whose entries are rational expressions in the scalar moduli of each layer and in the layer radii.

In the ALV setting every scalar modulus becomes its   trapezoidal Volterra matrix, every product a matrix product and every reciprocal a Volterra inverse; the transfer matrices become block matrices of size   and   respectively. Both harmonics carry over, together with the per-layer concentration kernels and the imperfect-interface transfers. The composite sphere then enters the schemes exactly as in elasticity — through its volume-averaged concentration kernel, having no Hill tensor of its own.

In conduction the temperature of layer is    and, with the flux analog   , the radial component is    , the amplitudes , being Volterra operators that act on the history of the remote gradient. A jump is the outer value minus the inner one. A Kapitza interface of resistance and a surface-conductive one of conductance give, for this harmonic,

the minus sign of the Kapitza law being that of the dictionary, and the amplitudes beyond the interface follow in closed form,

so that only is ever inverted and an impermeable core needs no inverse. As in elasticity, the concentration of the whole sphere adds to the material averages   the temperature jumps  , and its flux the surface current    of the surface-conductive interfaces, the outer interface counting unless it is left to the matrix.

7. Symmetry classes and structured storage ​

ALV operators inherit the symmetry classes of their elastic counterparts, and those classes are closed under Volterra product and inverse. That closure is what makes compact storage possible:

ClassStored componentsFull  Storage costClosure operation
ISO (18×)scalar Volterra products / inverses
TI (6×)  block-Volterra + 2 scalars
ORTHO (3×)  block-Volterra + 3 scalars
Genericfull    block-LU

with ISO ⊂ TI ⊂ ORTHO ⊂ generic. The types ALVKernelISO, ALVKernelTI and ALVKernelOrtho wrap these compact representations as AbstractMatrix, so they flow through generic Julia matrix code while preserving both the storage saving and the algebraic closure.

Where to go next ​

The theory closes with two appendices, starting with The finite Eshelby cell, the numerical device behind the finite-element inclusions. The laminate of Periodic multilayer — the laminate cell is transposed to aging viscoelasticity in its §10, by the substitution of this page. The schemes of §5 are compared on one composite in the tutorial Aging viscoelastic schemes side by side, and the manual page Viscoelastic homogenization gives the syntax of homogenize_alv.