Homogenization

This module provides mean-field homogenization methods for predicting effective properties of composite materials.

Eshelby Tensor

The foundation of mean-field homogenization is the Eshelby tensor, which describes the strain field inside an ellipsoidal inclusion embedded in an infinite matrix.

Homogenization Schemes

Functions for computing effective properties using Mori-Tanaka, self-consistent, and other mean-field schemes.

class cylinder_multi : public phase_multi
#include <cylinder_multi.hpp>

Cylindrical inclusion class for fiber-reinforced composite homogenization.

This class implements the mechanics of cylindrical (fiber) inclusions for mean-field homogenization of fiber-reinforced composites. Cylindrical inclusions are a limiting case of ellipsoids with one infinite axis.

Applications:

  • Unidirectional fiber composites

  • Long fiber reinforced polymers

  • Carbon fiber composites

Concentration Tensors:

The strain concentration relates macroscopic to local strain:

\[\boldsymbol{\varepsilon}^{fiber} = \mathbf{A} : \bar{\boldsymbol{\varepsilon}} \]

Local and global concentration tensors are related by rotation:

\[\mathbf{A} = \mathbf{R}^T : \mathbf{A}_{loc} : \mathbf{R} \]

See also

phase_multi for base class documentation

See also

ellipsoid_multi for general ellipsoidal inclusions

Public Functions

cylinder_multi()

Default constructor.

cylinder_multi(const arma::mat&, const arma::mat&, const arma::mat&, const arma::mat&, const arma::vec&, const arma::mat&, const arma::mat&, const arma::mat&, const arma::mat&)

Full constructor with all parameters.

Parameters:
  • A – Strain concentration tensor (global)

  • A_start – Initial strain concentration tensor

  • B – Stress concentration tensor (global)

  • B_start – Initial stress concentration tensor

  • A_in – Inelastic concentration vector

  • T_loc – Local interaction tensor

  • T – Global interaction tensor

  • A_loc – Local strain concentration tensor

  • B_loc – Local stress concentration tensor

cylinder_multi(const cylinder_multi&)

Copy constructor.

Parameters:

cm – cylinder_multi object to copy

~cylinder_multi()

Destructor.

virtual cylinder_multi &operator=(const cylinder_multi&)

Assignment operator.

Parameters:

cm – cylinder_multi object to assign from

Returns:

Reference to this object

Public Members

arma::mat T_loc

Interaction tensor in local coordinates (6×6).

Computed for cylinder geometry in the local frame where the cylinder axis is aligned with x₃.

arma::mat T

Interaction tensor in global coordinates (6×6).

Rotated from local to global using cylinder orientation

arma::mat A_loc

Strain concentration tensor in local coordinates (6×6).

Local frame concentration tensor before rotation to global frame

arma::mat B_loc

Stress concentration tensor in local coordinates (6×6).

Local frame concentration tensor before rotation to global frame

Friends

friend std::ostream &operator<<(std::ostream&, const cylinder_multi&)

Output stream operator for debugging.

Parameters:
Returns:

Reference to output stream

class eshelby_memo
#include <ellipsoid_multi.hpp>

Ellipsoidal inclusion class for Eshelby-based homogenization.

This class implements the mechanics of ellipsoidal inclusions embedded in an infinite matrix, based on Eshelby’s classical solution (1957). It is used in mean-field homogenization methods such as:

  • Mori-Tanaka scheme

  • Self-consistent scheme

  • Differential scheme

Key Tensors:

The Eshelby tensor \( \mathbf{S} \) relates the constrained strain in an inclusion to its eigenstrain:

\[\boldsymbol{\varepsilon}^c = \mathbf{S} : \boldsymbol{\varepsilon}^* \]

Hill’s polarization tensor \( \mathbf{P} \) is related to the Eshelby tensor by:

\[\mathbf{P} = \mathbf{S} : \mathbf{L}_0^{-1} \]

The interaction tensor \( \mathbf{T} \) (Wu’s tensor) relates inclusion stress to the strain difference:

\[\mathbf{T} = \mathbf{L}_0 : (\mathbf{I} - \mathbf{S}) \]

Memo of Eshelby tensors within one homogenization pass.

See also

phase_multi for base class documentation

See also

Eshelby, J.D. (1957). The determination of the elastic field of an ellipsoidal inclusion, and related problems. Proc. R. Soc. Lond. A 241, 376-396.

The numerical Eshelby tensor \( \mathbf{S} \) of an inclusion depends only on its semi-axes and on the medium stiffness written in the inclusion frame. In a self-consistent or Mori-Tanaka pass every inclusion of the same geometry embedded in the same medium therefore shares one integration (the most expensive step of the scheme, about 0.35 ms per 50x50 quadrature). A memo is created by the scheme for one pass and handed to fillT(); the key is the exact (bitwise) pair (rotated medium stiffness, semi-axes), so a hit returns exactly the tensor a fresh integration would. It is a local object, never static: no state is shared between calls or threads.

Note

All tensors are computed in the local coordinate system of the ellipsoid and can be rotated to global coordinates using the ellipsoid orientation.

Public Functions

const arma::mat *find(const arma::mat &Lt_local, const double &a1, const double &a2, const double &a3) const

The memoised tensor for this key, or nullptr.

const arma::mat &store(const arma::mat &Lt_local, const double &a1, const double &a2, const double &a3, const arma::mat &S)

Store and return the tensor of this key.

void clear()
inline std::size_t size() const
class ellipsoid_multi : public phase_multi
#include <ellipsoid_multi.hpp>

Public Functions

ellipsoid_multi()

Default constructor.

ellipsoid_multi(const arma::mat&, const arma::mat&, const arma::mat&, const arma::mat&, const arma::vec&, const arma::mat&, const arma::mat&, const arma::mat&, const arma::mat&, const arma::mat&, const arma::mat&)

Full constructor with all parameters.

Parameters:
  • A – Strain concentration tensor

  • A_start – Initial strain concentration tensor

  • B – Stress concentration tensor

  • B_start – Initial stress concentration tensor

  • A_in – Inelastic concentration vector

  • S_loc – Local Eshelby tensor

  • P_loc – Local polarization tensor

  • T_loc – Local interaction tensor

  • T – Global interaction tensor

  • T_in_loc – Local inelastic interaction tensor

  • T_in – Global inelastic interaction tensor

ellipsoid_multi(const ellipsoid_multi&)

Copy constructor.

Parameters:

em – ellipsoid_multi object to copy

~ellipsoid_multi()

Destructor.

virtual void fillS_loc(const arma::mat&, const ellipsoid&)

Compute the Eshelby tensor in local coordinates.

Computes \( \mathbf{S} \) based on the matrix stiffness and ellipsoid geometry using numerical integration over the unit sphere.

Parameters:
  • L0 – Matrix (reference medium) stiffness tensor (6×6)

  • ell – Ellipsoid geometry (aspect ratios and orientation)

virtual void fillP_loc(const arma::mat&, const ellipsoid&)

Compute Hill’s polarization tensor in local coordinates.

Computes \( \mathbf{P} = \mathbf{S} : \mathbf{L}_0^{-1} \)

Parameters:
  • L0 – Matrix stiffness tensor (6×6)

  • ell – Ellipsoid geometry

virtual void fillT(const arma::mat&, const arma::mat&, const ellipsoid&, eshelby_memo* = nullptr)

Compute the interaction tensor.

Computes the interaction tensor \( \mathbf{T} \) that relates the stress in the inclusion to the strain mismatch with the matrix.

Parameters:
  • L0 – Matrix stiffness tensor (6×6)

  • Lt – Inclusion tangent stiffness tensor (6×6)

  • ell – Ellipsoid geometry

  • memo – Optional eshelby_memo of the current homogenization pass: the Eshelby tensor is taken from it when the same inclusion geometry has already been met in the same (rotated) medium, and stored in it otherwise.

virtual void fillT_iso(const arma::mat&, const arma::mat&, const ellipsoid&, eshelby_memo* = nullptr)

Compute the interaction tensor for isotropic matrix.

Optimized computation when the matrix is isotropic, using analytical expressions for the Eshelby tensor.

Parameters:
  • L0 – Isotropic matrix stiffness tensor (6×6)

  • Lt – Inclusion tangent stiffness tensor (6×6)

  • ell – Ellipsoid geometry

virtual void fillT_mec_in(const arma::mat&, const arma::mat&, const ellipsoid&, eshelby_memo* = nullptr)

Compute the inelastic interaction tensor.

Computes the tensor relating eigenstrain/transformation strain to the resulting stress perturbation.

Parameters:
  • L0 – Matrix stiffness tensor (6×6)

  • Lt – Inclusion tangent stiffness tensor (6×6)

  • ell – Ellipsoid geometry

arma::mat eshelby_memoized(eshelby_memo *memo, const arma::mat &Lt_local, const ellipsoid &ell)

Eshelby tensor of ell in the medium Lt_local (inclusion frame): reused from this inclusion’s last integration when the medium is unchanged, else from the pass memo, else integrated (and recorded in both).

virtual ellipsoid_multi &operator=(const ellipsoid_multi&)

Assignment operator.

Parameters:

em – ellipsoid_multi object to assign from

Returns:

Reference to this object

Public Members

arma::mat S_loc

Eshelby tensor in local coordinates (6×6 in Voigt notation).

The Eshelby tensor depends only on the ellipsoid aspect ratios and the matrix Poisson ratio (for isotropic matrix).

arma::mat P_loc

Hill’s polarization tensor in local coordinates (6×6 in Voigt notation).

Related to Eshelby tensor: \( \mathbf{P} = \mathbf{S} : \mathbf{L}_0^{-1} \)

arma::mat T_loc

Interaction tensor (Wu’s tensor) in local coordinates (6×6).

\( \mathbf{T}_{loc} = \mathbf{L}_0 : (\mathbf{I} - \mathbf{S}) \)

arma::mat T

Interaction tensor in global coordinates (6×6).

Rotated from local to global using ellipsoid orientation

arma::mat T_in_loc

Inelastic interaction tensor in local coordinates (6×6).

Used for eigenstrain/transformation strain contributions

arma::mat T_in

Inelastic interaction tensor in global coordinates (6×6).

eshelby_memo S_loc_memo

Key of the memoised S_loc: the (rotated) medium stiffness and semi-axes it was integrated for. Lets an inclusion embedded in an unchanged medium (elastic matrix, successive increments and Newton iterations) reuse its own Eshelby tensor.

Public Static Attributes

static int mp

Number of integration points in polar direction for numerical integration.

static int np

Number of integration points in azimuthal direction.

static arma::vec &x

Gauss points for polar integration.

The four Gauss vectors are computed once (costly) and shared by every ellipsoid. They are references to heap vectors that are never destroyed: a static armadillo object’s destructor would run at DLL unload, which is unsafe on Windows.

static arma::vec &wx

Gauss weights for polar integration.

static arma::vec &y

Gauss points for azimuthal integration.

static arma::vec &wy

Gauss weights for azimuthal integration.

Friends

friend std::ostream &operator<<(std::ostream&, const ellipsoid_multi&)

Output stream operator for debugging.

Parameters:
Returns:

Reference to output stream

class layer_multi : public phase_multi
#include <layer_multi.hpp>

Layer class for laminated composite homogenization.

This class implements the mechanics of layered composites for mean-field homogenization. It is used for materials with a stratified microstructure where layers are stacked in a specific direction.

Applications:

  • Laminated composites (e.g., carbon fiber laminates)

  • Layered geological media

  • Multilayer coatings

  • Sandwich structures

Coordinate System:

  • The stacking direction is aligned with the normal direction (n)

  • The in-plane directions are tangential (t)

Stress/Strain Decomposition:

For layered media, it is convenient to decompose stress and strain into normal (n) and tangential (t) components:

  • Normal stress components: \( \sigma_{nn}, \sigma_{nt} \)

  • Tangential strain components: \( \varepsilon_{tt} \)

The continuity conditions at layer interfaces impose:

  • Continuous normal stress: \( \sigma_n^{(1)} = \sigma_n^{(2)} \)

  • Continuous tangential strain: \( \varepsilon_t^{(1)} = \varepsilon_t^{(2)} \)

See also

phase_multi for base class documentation

Public Functions

layer_multi()

Default constructor.

layer_multi(const arma::mat&, const arma::mat&, const arma::mat&, const arma::mat&, const arma::vec&, const arma::mat&, const arma::mat&, const arma::mat&, const arma::mat&, const arma::vec&, const arma::vec&)

Full constructor with all parameters.

Parameters:
  • A – Strain concentration tensor

  • A_start – Initial strain concentration tensor

  • B – Stress concentration tensor

  • B_start – Initial stress concentration tensor

  • A_in – Inelastic concentration vector

  • Dnn – Normal-normal stiffness block

  • Dnt – Normal-tangential stiffness block

  • dXn – Deformation gradient derivative (normal)

  • dXt – Deformation gradient derivative (tangential)

  • sigma_hat – Interface stress variable

  • dzdx1 – Field gradient through thickness

layer_multi(const layer_multi&)

Copy constructor.

Parameters:

lm – layer_multi object to copy

~layer_multi()

Destructor.

virtual layer_multi &operator=(const layer_multi&)

Assignment operator.

Parameters:

lm – layer_multi object to assign from

Returns:

Reference to this object

Public Members

arma::mat Dnn

Normal-normal block of the tangent modulus (3×3).

Submatrix of the stiffness tensor relating normal stress to normal strain: \( \mathbf{D}_{nn} = \frac{\partial \boldsymbol{\sigma}_n}{\partial \boldsymbol{\varepsilon}_n} \)

arma::mat Dnt

Normal-tangential block of the tangent modulus (3×3).

Submatrix coupling normal stress to tangential strain: \( \mathbf{D}_{nt} = \frac{\partial \boldsymbol{\sigma}_n}{\partial \boldsymbol{\varepsilon}_t} \)

arma::mat dXn

Derivative of deformation gradient with respect to normal position (3×3).

Used in incremental formulations for geometric nonlinearity

arma::mat dXt

Derivative of deformation gradient with respect to tangential position (3×3).

Used in incremental formulations for geometric nonlinearity

arma::vec sigma_hat

Stress-like quantity for interface conditions (6×1).

Intermediate variable used in the layer homogenization algorithm

arma::vec dzdx1

Derivative of local field with respect to layer position (6×1).

Used for computing gradients through the layer thickness

Friends

friend std::ostream &operator<<(std::ostream&, const layer_multi&)

Output stream operator for debugging.

Parameters:
  • os – Output stream

  • lm – layer_multi object to output

Returns:

Reference to output stream

class phase_multi
#include <phase_multi.hpp>

Base class for phases in multiphase homogenization schemes.

The phase_multi class serves as the base class for all phase types used in mean-field homogenization methods. It stores the concentration tensors that relate macroscopic strain/stress to local phase-level quantities.

Concentration Tensor Relationships:

The strain concentration tensor \( \mathbf{A} \) relates macroscopic strain to local strain:

\[\boldsymbol{\varepsilon}^{(r)} = \mathbf{A}^{(r)} : \bar{\boldsymbol{\varepsilon}} \]

The stress concentration tensor \( \mathbf{B} \) relates macroscopic stress to local stress:

\[\boldsymbol{\sigma}^{(r)} = \mathbf{B}^{(r)} : \bar{\boldsymbol{\sigma}} \]

For inelastic behavior, an additional concentration vector handles eigenstrain contributions:

\[\boldsymbol{\varepsilon}^{(r)} = \mathbf{A}^{(r)} : \bar{\boldsymbol{\varepsilon}} + \mathbf{A}_{in}^{(r)} \]

Derived Classes:

Subclassed by cylinder_multi, ellipsoid_multi, layer_multi

Public Functions

phase_multi()

Default constructor.

Initializes concentration tensors to zero matrices/vectors of appropriate size

phase_multi(const arma::mat &A_init, const arma::mat &A_start_init, const arma::mat &B_init, const arma::mat &B_start_init, const arma::vec &A_in_init)

Constructor with parameters.

Parameters:
  • A_init – Initial strain concentration tensor (6×6)

  • A_start_init – Initial strain concentration tensor for start state (6×6)

  • B_init – Initial stress concentration tensor (6×6)

  • B_start_init – Initial stress concentration tensor for start state (6×6)

  • A_in_init – Initial inelastic concentration vector (6×1)

phase_multi(const phase_multi &pm)

Copy constructor.

Parameters:

pm – phase_multi object to copy

virtual ~phase_multi()

Virtual destructor.

void to_start()

Reset current tensors to start-of-increment values.

Sets A = A_start and B = B_start for iterative schemes

void set_start()

Store current tensors as start-of-increment values.

Sets A_start = A and B_start = B after convergence

virtual phase_multi &operator=(const phase_multi &pm)

Assignment operator.

Parameters:

pm – phase_multi object to assign from

Returns:

Reference to this object

Public Members

arma::mat A

Strain concentration tensor (6×6 matrix in Voigt notation).

Relates macroscopic strain to local phase strain: \( \boldsymbol{\varepsilon}^{local} = \mathbf{A} : \bar{\boldsymbol{\varepsilon}} \)

arma::mat A_start

Strain concentration tensor at the start of increment.

Stored for incremental schemes and convergence checks

arma::mat B

Stress concentration tensor (6×6 matrix in Voigt notation).

Relates macroscopic stress to local phase stress: \( \boldsymbol{\sigma}^{local} = \mathbf{B} : \bar{\boldsymbol{\sigma}} \)

arma::mat B_start

Stress concentration tensor at the start of increment.

Stored for incremental schemes and convergence checks

arma::vec A_in

Inelastic strain concentration vector (6×1 in Voigt notation).

Accounts for eigenstrain/transformation strain contributions: \( \boldsymbol{\varepsilon}^{local} = \mathbf{A} : \bar{\boldsymbol{\varepsilon}} + \mathbf{A}_{in} \)

Friends

friend std::ostream &operator<<(std::ostream &os, const phase_multi &pm)

Output stream operator for debugging.

Parameters:
  • os – Output stream

  • pm – phase_multi object to output

Returns:

Reference to output stream