felupe.constitution#
- class felupe.constitution.LinearElastic(E=None, nu=None)[source]#
Isotropic linear-elastic material formulation.
\[\begin{split}\begin{bmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{12} \\ \sigma_{23} \\ \sigma_{31} \end{bmatrix} = \frac{E}{(1+\nu)(1-2\nu)}\begin{bmatrix} 1-\nu & \nu & \nu & 0 & 0 & 0\\ \nu & 1-\nu & \nu & 0 & 0 & 0\\ \nu & \nu & 1-\nu & 0 & 0 & 0\\ 0 & 0 & 0 & \frac{1-2\nu}{2} & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1-2\nu}{2} & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{1-2\nu}{2} \end{bmatrix} \cdot \begin{bmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \varepsilon_{33} \\ 2 \varepsilon_{12} \\ 2 \varepsilon_{23} \\ 2 \varepsilon_{31} \end{bmatrix}\end{split}\]with the strain tensor
\[\boldsymbol{\varepsilon} = \frac{1}{2} \left( \frac{\partial \boldsymbol{u}}{\partial \boldsymbol{X}} + \left( \frac{\partial \boldsymbol{u}}{\partial \boldsymbol{X}} \right)^T \right)\]- Parameters
E (float) – Young’s modulus.
nu (float) – Poisson ratio.
- gradient(extract, E=None, nu=None)[source]#
Evaluate the stress tensor (as a function of the deformation gradient).
- Parameters
extract (list of ndarray) – List with Deformation gradient
F
(3x3) as first itemE (float, optional) – Young’s modulus (default is None)
nu (float, optional) – Poisson ratio (default is None)
- Returns
Stress tensor (3x3)
- Return type
ndarray
- hessian(extract=None, E=None, nu=None, shape=(1, 1), region=None)[source]#
Evaluate the elasticity tensor. The Deformation gradient is only used for the shape of the trailing axes.
- Parameters
extract (list of ndarray, optional) – List with Deformation gradient
F
(3x3) as first item (default is None)E (float, optional) – Young’s modulus (default is None)
nu (float, optional) – Poisson ratio (default is None)
shape ((int, int), optional (default is (1, 1))) – Tuple with shape of the trailing axes (default is None)
region (Region, optional) – A numeric region for shape of the trailing axes (default is None)
- Returns
elasticity tensor (3x3x3x3)
- Return type
ndarray
- class felupe.constitution.LinearElasticTensorNotation(E=None, nu=None, parallel=False)[source]#
Isotropic linear-elastic material formulation.
\[ \begin{align}\begin{aligned}\boldsymbol{\sigma} &= 2 \mu \ \boldsymbol{\varepsilon} + \gamma \ \text{tr}(\boldsymbol{\varepsilon}) \ \boldsymbol{I}\\\frac{\boldsymbol{\partial \sigma}}{\partial \boldsymbol{\varepsilon}} &= 2 \mu \ \boldsymbol{I} \odot \boldsymbol{I} + \gamma \ \boldsymbol{I} \otimes \boldsymbol{I}\end{aligned}\end{align} \]with the strain tensor
\[\boldsymbol{\varepsilon} = \frac{1}{2} \left( \frac{\partial \boldsymbol{u}}{\partial \boldsymbol{X}} + \left( \frac{\partial \boldsymbol{u}}{\partial \boldsymbol{X}} \right)^T \right)\]- Parameters
E (float) – Young’s modulus.
nu (float) – Poisson ratio.
- gradient(extract=None, E=None, nu=None)[source]#
Evaluate the stress tensor (as a function of the deformation gradient).
- Parameters
extract (list of ndarray) – List with Deformation gradient
F
(3x3) as first itemE (float, optional) – Young’s modulus (default is None)
nu (float, optional) – Poisson ratio (default is None)
- Returns
Stress tensor (3x3)
- Return type
ndarray
- hessian(extract=None, E=None, nu=None, shape=(1, 1), region=None)[source]#
Evaluate the elasticity tensor. The Deformation gradient is only used for the shape of the trailing axes.
- Parameters
extract (list of ndarray. optional) – List with Deformation gradient
F
(3x3) as first item (default is None)E (float, optional) – Young’s modulus (default is None)
nu (float, optional) – Poisson ratio (default is None)
shape ((int, int), optional (default is (1, 1))) – Tuple with shape of the trailing axes (default is None)
region (Region, optional) – A numeric region for shape of the trailing axes (default is None)
- Returns
elasticity tensor (3x3x3x3)
- Return type
ndarray
- class felupe.constitution.LinearElasticPlaneStress(E, nu)[source]#
Plane-stress isotropic linear-elastic material formulation.
- Parameters
E (float) – Young’s modulus.
nu (float) – Poisson ratio.
- gradient(extract, E=None, nu=None)[source]#
Evaluate the 2d-stress tensor from the deformation gradient.
- Parameters
extract (list of ndarray) – List with In-plane components (2x2) of the Deformation gradient ``F``as first item
E (float, optional) – Young’s modulus (default is None)
nu (float, optional) – Poisson ratio (default is None)
- Returns
In-plane components of stress tensor (2x2)
- Return type
ndarray
- hessian(extract=None, E=None, nu=None, shape=(1, 1), region=None)[source]#
Evaluate the elasticity tensor from the deformation gradient.
- Parameters
extract (list of ndarray, optional) – List with In-plane components (2x2) of the Deformation gradient ``F``as first item (default is None)
E (float, optional) – Young’s modulus (default is None)
nu (float, optional) – Poisson ratio (default is None)
shape ((int, int), optional (default is (1, 1))) – Tuple with shape of the trailing axes (default is None)
region (Region, optional) – A numeric region for shape of the trailing axes (default is None)
- Returns
In-plane components of elasticity tensor (2x2x2x2)
- Return type
ndarray
- strain(extract, E=None, nu=None)[source]#
Evaluate the strain tensor from the deformation gradient.
- Parameters
extract (list of ndarray) – List with In-plane components (2x2) of the Deformation gradient ``F``as first item
E (float, optional) – Young’s modulus (default is None)
nu (float, optional) – Poisson ratio (default is None)
- Returns
e – Strain tensor (3x3)
- Return type
ndarray
- stress(extract, E=None, nu=None)[source]#
“Evaluate the 3d-stress tensor from the deformation gradient.
- Parameters
extract (list of ndarray) – List with In-plane components (2x2) of the Deformation gradient ``F``as first item
E (float, optional) – Young’s modulus (default is None)
nu (float, optional) – Poisson ratio (default is None)
- Returns
Stress tensor (3x3)
- Return type
ndarray
- class felupe.constitution.LinearElasticPlaneStrain(E, nu)[source]#
Plane-strain isotropic linear-elastic material formulation.
- Parameters
E (float) – Young’s modulus.
nu (float) – Poisson ratio.
- gradient(extract, E=None, nu=None)[source]#
Evaluate the 2d-stress tensor from the deformation gradient.
- Parameters
extract (list of ndarray) – List with In-plane components (2x2) of the Deformation gradient ``F``as first item
E (float, optional) – Young’s modulus (default is None)
nu (float, optional) – Poisson ratio (default is None)
- Returns
In-plane components of stress tensor (2x2)
- Return type
ndarray
- hessian(extract, E=None, nu=None)[source]#
Evaluate the 2d-elasticity tensor from the deformation gradient.
- Parameters
extract (list of ndarray) – List with In-plane components (2x2) of the Deformation gradient ``F``as first item
E (float, optional) – Young’s modulus (default is None)
nu (float, optional) – Poisson ratio (default is None)
- Returns
In-plane components of elasticity tensor (2x2x2x2)
- Return type
ndarray
- strain(extract, E=None, nu=None)[source]#
Evaluate the strain tensor from the deformation gradient.
- Parameters
extract (list of ndarray) – List with In-plane components (2x2) of the Deformation gradient ``F``as first item
E (float, optional) – Young’s modulus (default is None)
nu (float, optional) – Poisson ratio (default is None)
- Returns
e – Strain tensor (3x3)
- Return type
ndarray
- stress(extract, E=None, nu=None)[source]#
“Evaluate the 3d-stress tensor from the deformation gradient.
- Parameters
extract (list of ndarray) – List with In-plane components (2x2) of the Deformation gradient ``F``as first item
E (float, optional) – Young’s modulus (default is None)
nu (float, optional) – Poisson ratio (default is None)
- Returns
Stress tensor (3x3)
- Return type
ndarray
- class felupe.constitution.NeoHooke(mu=None, bulk=None, parallel=False)[source]#
Nearly-incompressible isotropic hyperelastic Neo-Hooke material formulation. The strain energy density function of the Neo-Hookean material formulation is a linear function of the trace of the isochoric part of the right Cauchy-Green deformation tensor.
In a nearly-incompressible constitutive framework the strain energy density is an additive composition of an isochoric and a volumetric part. While the isochoric part is defined on the distortional part of the deformation gradient, the volumetric part of the strain energy function is defined on the determinant of the deformation gradient.
\[ \begin{align}\begin{aligned}\psi &= \hat{\psi}(\hat{\boldsymbol{C}}) + U(J)\\\hat\psi(\hat{\boldsymbol{C}}) &= \frac{\mu}{2} (\text{tr}(\hat{\boldsymbol{C}}) - 3)\end{aligned}\end{align} \]with
\[ \begin{align}\begin{aligned}J &= \text{det}(\boldsymbol{F})\\\hat{\boldsymbol{F}} &= J^{-1/3} \boldsymbol{F}\\\hat{\boldsymbol{C}} &= \hat{\boldsymbol{F}}^T \hat{\boldsymbol{F}}\end{aligned}\end{align} \]The volumetric part of the strain energy density function is a function the volume ratio.
\[U(J) = \frac{K}{2} (J - 1)^2\]The first Piola-Kirchhoff stress tensor is evaluated as the gradient of the strain energy density function. The hessian of the strain energy density function enables the corresponding elasticity tensor.
\[ \begin{align}\begin{aligned}\boldsymbol{P} &= \frac{\partial \psi}{\partial \boldsymbol{F}}\\\mathbb{A} &= \frac{\partial^2 \psi}{\partial \boldsymbol{F}\ \partial \boldsymbol{F}}\end{aligned}\end{align} \]A chain rule application leads to the following expression for the stress tensor. It is formulated as a sum of the physical-deviatoric (not the mathematical deviator!) and the physical-hydrostatic stress tensors.
\[ \begin{align}\begin{aligned}\boldsymbol{P} &= \boldsymbol{P}' + \boldsymbol{P}_U\\\boldsymbol{P}' &= \frac{\partial \hat{\psi}}{\partial \hat{\boldsymbol{F}}} : \frac{\partial \hat{\boldsymbol{F}}}{\partial \boldsymbol{F}} = \bar{\boldsymbol{P}} - \frac{1}{3} (\bar{\boldsymbol{P}} : \boldsymbol{F}) \boldsymbol{F}^{-T}\\\boldsymbol{P}_U &= \frac{\partial U(J)}{\partial J} \frac{\partial J}{\partial \boldsymbol{F}} = U'(J) J \boldsymbol{F}^{-T}\end{aligned}\end{align} \]with
\[ \begin{align}\begin{aligned}\frac{\partial \hat{\boldsymbol{F}}}{\partial \boldsymbol{F}} &= J^{-1/3} \left( \boldsymbol{I} \overset{ik}{\otimes} \boldsymbol{I} - \frac{1}{3} \boldsymbol{F} \otimes \boldsymbol{F}^{-T} \right)\\\frac{\partial J}{\partial \boldsymbol{F}} &= J \boldsymbol{F}^{-T}\\\bar{\boldsymbol{P}} &= J^{-1/3} \frac{\partial \hat{\psi}}{\partial \hat{\boldsymbol{F}}}\end{aligned}\end{align} \]With the above partial derivatives the first Piola-Kirchhoff stress tensor of the Neo-Hookean material model takes the following form.
\[\boldsymbol{P} = \mu J^{-2/3} \left( \boldsymbol{F} - \frac{1}{3} (\boldsymbol{F} : \boldsymbol{F}) \boldsymbol{F}^{-T} \right) + K (J - 1) J \boldsymbol{F}^{-T}\]Again, a chain rule application leads to an expression for the elasticity tensor.
\[ \begin{align}\begin{aligned}\mathbb{A} &= \mathbb{A}' + \mathbb{A}_{U}\\\mathbb{A}' &= \bar{\mathbb{A}} - \frac{1}{3} \left( (\bar{\mathbb{A}} : \boldsymbol{F}) \otimes \boldsymbol{F}^{-T} + \boldsymbol{F}^{-T} \otimes (\boldsymbol{F} : \bar{\mathbb{A}}) \right ) + \frac{1}{9} (\boldsymbol{F} : \bar{\mathbb{A}} : \boldsymbol{F}) \boldsymbol{F}^{-T} \otimes \boldsymbol{F}^{-T}\\\mathbb{A}_{U} &= (U''(J) J + U'(J)) J \boldsymbol{F}^{-T} \otimes \boldsymbol{F}^{-T} - U'(J) J \boldsymbol{F}^{-T} \overset{il}{\otimes} \boldsymbol{F}^{-T}\end{aligned}\end{align} \]with
\[\bar{\mathbb{A}} = J^{-1/3} \frac{\partial^2 \hat\psi}{\partial \hat{\boldsymbol{F}}\ \partial \hat{\boldsymbol{F}}} J^{-1/3}\]With the above partial derivatives the (physical-deviatoric and -hydrostatic) parts of the elasticity tensor associated to the first Piola-Kirchhoff stress tensor of the Neo-Hookean material model takes the following form.
\[ \begin{align}\begin{aligned}\mathbb{A} &= \mathbb{A}' + \mathbb{A}_{U}\\\mathbb{A}' &= J^{-2/3} \left(\boldsymbol{I} \overset{ik}{\otimes} \boldsymbol{I} - \frac{1}{3} \left( \boldsymbol{F} \otimes \boldsymbol{F}^{-T} + \boldsymbol{F}^{-T} \otimes \boldsymbol{F} \right ) + \frac{1}{9} (\boldsymbol{F} : \boldsymbol{F}) \boldsymbol{F}^{-T} \otimes \boldsymbol{F}^{-T} \right)\\\mathbb{A}_{U} &= K J \left( (2J - 1) \boldsymbol{F}^{-T} \otimes \boldsymbol{F}^{-T} - (J - 1) \boldsymbol{F}^{-T} \overset{il}{\otimes} \boldsymbol{F}^{-T} \right)\end{aligned}\end{align} \]- Parameters
mu (float) – Shear modulus
bulk (float) – Bulk modulus
- function(extract, mu=None, bulk=None)[source]#
Strain energy density function per unit undeformed volume of the Neo-Hookean material formulation.
- Parameters
extract (list of ndarray) – List with the Deformation gradient
F
(3x3) as first itemmu (float, optional) – Shear modulus (default is None)
bulk (float, optional) – Bulk modulus (default is None)
- gradient(extract, mu=None, bulk=None)[source]#
Gradient of the strain energy density function per unit undeformed volume of the Neo-Hookean material formulation.
- Parameters
extract (list of ndarray) – List with the Deformation gradient
F
(3x3) as first itemmu (float, optional) – Shear modulus (default is None)
bulk (float, optional) – Bulk modulus (default is None)
- hessian(extract, mu=None, bulk=None)[source]#
Hessian of the strain energy density function per unit undeformed volume of the Neo-Hookean material formulation.
- Parameters
extract (list of ndarray) – List with the Deformation gradient
F
(3x3) as first itemmu (float, optional) – Shear modulus (default is None)
bulk (float, optional) – Bulk modulus (default is None)
- class felupe.constitution.ThreeFieldVariation(material, parallel=False)[source]#
Hu-Washizu hydrostatic-volumetric selective \((\boldsymbol{u},p,J)\) - three-field variation for nearly- incompressible material formulations. The total potential energy for nearly-incompressible hyperelasticity is formulated with a determinant-modified deformation gradient. Pressure and volume ratio fields should be kept one order lower than the interpolation order of the displacement field, e.g. linear displacement fields should be paired with element-constant (mean) values of pressure and volume ratio.
The total potential energy of internal forces is defined with a strain energy density function in terms of a determinant-modified deformation gradient and an additional control equation.
\[ \begin{align}\begin{aligned}\Pi &= \Pi_{int} + \Pi_{ext}\\\Pi_{int} &= \int_V \psi(\boldsymbol{F}) \ dV \qquad \rightarrow \qquad \Pi_{int}(\boldsymbol{u},p,J) = \int_V \psi(\overline{\boldsymbol{F}}) \ dV + \int_V p (J-\overline{J}) \ dV\\\overline{\boldsymbol{F}} &= \left(\frac{\overline{J}}{J}\right)^{1/3} \boldsymbol{F}\end{aligned}\end{align} \]The variations of the total potential energy w.r.t. \((\boldsymbol{u},p,J)\) lead to the following expressions. We denote first partial derivatives as \(\boldsymbol{f}_{(\bullet)}\) and second partial derivatives as \(\boldsymbol{A}_{(\bullet,\bullet)}\).
\[ \begin{align}\begin{aligned}\delta_{\boldsymbol{u}} \Pi_{int} &= \int_V \boldsymbol{f}_{\boldsymbol{u}} : \delta \boldsymbol{F} \ dV = \int_V \left( \frac{\partial \psi}{\partial \overline{\boldsymbol{F}}} : \frac{\partial \overline{\boldsymbol{F}}}{\partial \boldsymbol{F}} + p J \boldsymbol{F}^{-T} \right) : \delta \boldsymbol{F} \ dV\\\delta_{p} \Pi_{int} &= \int_V f_{p} \ \delta p \ dV = \int_V (J - \overline{J}) \ \delta p \ dV\\\delta_{\overline{J}} \Pi_{int} &= \int_V f_{\overline{J}} \ \delta \overline{J} \ dV = \int_V \left( \frac{\partial \psi}{\partial \overline{\boldsymbol{F}}} : \frac{\partial \overline{\boldsymbol{F}}}{\partial \overline{J}} - p \right) : \delta \overline{J} \ dV\end{aligned}\end{align} \]The projection tensors from the variations lead the following results.
\[ \begin{align}\begin{aligned}\frac{\partial \overline{\boldsymbol{F}}}{\partial \boldsymbol{F}} &= \left(\frac{\overline{J}}{J}\right)^{1/3} \left( \boldsymbol{I} \overset{ik}{\odot} \boldsymbol{I} - \frac{1}{3} \boldsymbol{F} \otimes \boldsymbol{F}^{-T} \right)\\\frac{\partial \overline{\boldsymbol{F}}}{\partial \overline{J}} &= \frac{1}{3 \overline{J}} \overline{\boldsymbol{F}}\end{aligned}\end{align} \]The double-dot products from the variations are now evaluated.
\[ \begin{align}\begin{aligned}\overline{\boldsymbol{P}} &= \frac{\partial \psi}{\partial \overline{\boldsymbol{F}}} = \overline{\overline{\boldsymbol{P}}} - \frac{1}{3} \left( \overline{\overline{\boldsymbol{P}}} : \boldsymbol{F} \right) \boldsymbol{F}^{-T} \qquad \text{with} \qquad \overline{\overline{\boldsymbol{P}}} = \left(\frac{\overline{J}}{J}\right)^{1/3} \frac{\partial \psi}{\partial \overline{\boldsymbol{F}}}\\\frac{\partial \psi}{\partial \overline{\boldsymbol{F}}} : \frac{1}{3 \overline{J}} \overline{\boldsymbol{F}} &= \frac{1}{3 \overline{J}} \overline{\overline{\boldsymbol{P}}} : \boldsymbol{F}\end{aligned}\end{align} \]We now have three formulas; one for the first Piola Kirchhoff stress and two additional control equations.
\[ \begin{align}\begin{aligned}\boldsymbol{f}_{\boldsymbol{u}} (= \boldsymbol{P}) &= \overline{\overline{\boldsymbol{P}}} - \frac{1}{3} \left( \overline{\overline{\boldsymbol{P}}} : \boldsymbol{F} \right) \boldsymbol{F}^{-T}\\f_p &= J - \overline{J}\\f_{\overline{J}} &= \frac{1}{3 \overline{J}} \left( \overline{\overline{\boldsymbol{P}}} : \boldsymbol{F} \right) - p\end{aligned}\end{align} \]A linearization of the above formulas gives six equations (only results are given here).
\[ \begin{align}\begin{aligned}\mathbb{A}_{\boldsymbol{u},\boldsymbol{u}} &= \overline{\overline{\mathbb{A}}} + \frac{1}{9} \left( \boldsymbol{F} : \overline{\overline{\mathbb{A}}} : \boldsymbol{F} \right) \boldsymbol{F}^{-T} \otimes \boldsymbol{F}^{-T} - \frac{1}{3} \left( \boldsymbol{F}^{-T} \otimes \left( \overline{\overline{\boldsymbol{P}}} + \boldsymbol{F} : \overline{\overline{\mathbb{A}}} \right) + \left( \overline{\overline{\boldsymbol{P}}} + \overline{\overline{\mathbb{A}}} : \boldsymbol{F} \right) \otimes \boldsymbol{F}^{-T} \right)\\&+\left( p J + \frac{1}{9} \overline{\overline{\boldsymbol{P}}} : \boldsymbol{F} \right) \boldsymbol{F}^{-T} \otimes \boldsymbol{F}^{-T} - \left( p J - \frac{1}{3} \overline{\overline{\boldsymbol{P}}} : \boldsymbol{F} \right) \boldsymbol{F}^{-T} \overset{il}{\odot} \boldsymbol{F}^{-T}\\A_{p,p} &= 0\\A_{\overline{J},\overline{J}} &= \frac{1}{9 \overline{J}^2} \left( \boldsymbol{F} : \overline{\overline{\mathbb{A}}} : \boldsymbol{F} \right) - 2 \left( \overline{\overline{\boldsymbol{P}}} : \boldsymbol{F} \right)\\\boldsymbol{A}_{\boldsymbol{u},p} &= \boldsymbol{A}_{p, \boldsymbol{u}} = J \boldsymbol{F}^{-T}\\\boldsymbol{A}_{\boldsymbol{u},\overline{J}} &= \boldsymbol{A}_{\overline{J}, \boldsymbol{u}} = \frac{1}{3 \overline{J}} \left( \boldsymbol{P}' + \boldsymbol{F} : \overline{\overline{\mathbb{A}}} - \frac{1}{3} \left( \boldsymbol{F} : \overline{\overline{\mathbb{A}}} : \boldsymbol{F} \right) \boldsymbol{F}^{-T} \right)\\A_{p,\overline{J}} &= A_{\overline{J}, p} = -1\end{aligned}\end{align} \]with
\[\overline{\overline{\mathbb{A}}} = \left(\frac{\overline{J}}{J}\right)^{1/3} \frac{\partial^2 \psi}{\partial \overline{\boldsymbol{F}} \partial \overline{\boldsymbol{F}}} \left(\frac{\overline{J}}{J}\right)^{1/3}\]as well as
\[\boldsymbol{P}' = \boldsymbol{P} - p J \boldsymbol{F}^{-T}\]- Parameters
material (Material) – A material definition with
gradient
andhessian
methods.
- fun_P#
Method for gradient evaluation
- Type
function
- fun_A#
Method for hessian evaluation
- Type
function
- detF#
Determinant of deformation gradient
- Type
ndarray
- iFT#
Transpose of inverse of the deformation gradient
- Type
ndarray
- Fb#
Determinant-modified deformation gradient
- Type
ndarray
- Pb#
First Piola-Kirchhoff stress tensor (in determinant-modified framework)
- Type
ndarray
- Pbb#
Determinant-modification multiplied by
Pb
- Type
ndarray
- PbbF#
Double-dot product of
Pb
and the deformation gradient- Type
ndarray
- gradient(extract)[source]#
List of variations of total potential energy w.r.t displacements, pressure and volume ratio.
δ_u(Π_int) = ∫_V (∂ψ/∂F + p cof(F)) : δF dV δ_p(Π_int) = ∫_V (det(F) - J) δp dV δ_J(Π_int) = ∫_V (∂U/∂J - p) δJ dV
- Parameters
extract (list of ndarray) – List of extracted field values with Deformation gradient
F
as first, the hydrostatic pressurep
as second and the volume ratioJ
as third item.- Returns
List of gradients w.r.t. the input variables F, p and J
- Return type
list of ndarrays
- hessian(extract)[source]#
List of linearized variations of total potential energy w.r.t displacements, pressure and volume ratio (these expressions are symmetric;
A_up = A_pu
if derived from a total potential energy formulation). List entries have to be arranged as a flattened list from the upper triangle blocks:Δ_u(δ_u(Π_int)) = ∫_V δF : (∂²ψ/(∂F∂F) + p ∂cof(F)/∂F) : ΔF dV Δ_p(δ_u(Π_int)) = ∫_V δF : J cof(F) Δp dV Δ_J(δ_u(Π_int)) = ∫_V δF : ∂²ψ/(∂F∂J) ΔJ dV Δ_p(δ_p(Π_int)) = ∫_V δp 0 Δp dV Δ_J(δ_p(Π_int)) = ∫_V δp (-1) ΔJ dV Δ_J(δ_J(Π_int)) = ∫_V δJ ∂²ψ/(∂J∂J) ΔJ dV [[0 1 2], [ 3 4], [ 5]] --> [0 1 2 3 4 5]
- Parameters
extract (list of ndarray) – List of extracted field values with Deformation gradient
F
as first, the hydrostatic pressurep
as second and the volume ratioJ
as third item.- Returns
List of hessians in upper triangle order
- Return type
list of ndarrays
- class felupe.constitution.LineChange(parallel=False)[source]#
Line Change.
\[d\boldsymbol{x} = \boldsymbol{F} d\boldsymbol{X}\]Gradient:
\[\frac{\partial \boldsymbol{F}}{\partial \boldsymbol{F}} = \boldsymbol{I} \overset{ik}{\otimes} \boldsymbol{I}\]
- class felupe.constitution.AreaChange(parallel=False)[source]#
Area Change.
\[d\boldsymbol{a} = J \boldsymbol{F}^{-T} d\boldsymbol{A}\]Gradient:
\[\frac{\partial J \boldsymbol{F}^{-T}}{\partial \boldsymbol{F}} = J \left( \boldsymbol{F}^{-T} \otimes \boldsymbol{F}^{-T} - \boldsymbol{F}^{-T} \overset{il}{\otimes} \boldsymbol{F}^{-T} \right)\]- function(extract, N=None, parallel=None)[source]#
Area change.
- Parameters
extract (list of ndarray) – List of extracted field values with Deformation gradient as first item.
N (ndarray or None, optional) – Area normal vector (default is None)
- Returns
Cofactor matrix of the deformation gradient
- Return type
ndarray
- gradient(extract, N=None, parallel=None)[source]#
Gradient of area change.
- Parameters
extract (list of ndarray) – List of extracted field values with Deformation gradient as first item.
N (ndarray or None, optional) – Area normal vector (default is None)
- Returns
Gradient of cofactor matrix of the deformation gradient
- Return type
ndarray
- class felupe.constitution.VolumeChange(parallel=False)[source]#
Volume Change.
\[d\boldsymbol{v} = \text{det}(\boldsymbol{F}) d\boldsymbol{V}\]Gradient and hessian (equivalent to gradient of area change):
\[ \begin{align}\begin{aligned}\frac{\partial J}{\partial \boldsymbol{F}} &= J \boldsymbol{F}^{-T}\\\frac{\partial^2 J}{\partial \boldsymbol{F}\ \partial \boldsymbol{F}} &= J \left( \boldsymbol{F}^{-T} \otimes \boldsymbol{F}^{-T} - \boldsymbol{F}^{-T} \overset{il}{\otimes} \boldsymbol{F}^{-T} \right)\end{aligned}\end{align} \]- function(extract)[source]#
Gradient of volume change.
- Parameters
extract (list of ndarray) – List of extracted field values with Deformation gradient as first item.
- Returns
J – Determinant of the deformation gradient
- Return type
ndarray