Displacement Finite Element

Theory

This element implements a small-strain displacement-based Finite Element formulation. The incremental strain is computed from the displacement increment via the kinematic matrix B, and internal forces and consistent tangents are evaluated by Gauss quadrature.

\[\Delta \eps = \mathbf{B}\, \Delta \mathbf{u}, \qquad \mathbf{K}_e = \sum_{qp} \mathbf{B}^\mathsf{T}\, \mathbf{C}\, \mathbf{B}\, J_0 w, \qquad \mathbf{P}_e = \sum_{qp} \mathbf{B}^\mathsf{T}\, \sig\, J_0 w.\]

where \(\Delta \eps\) is the linearized strain tensor in Voigt notation, \(\mathbf{B}\) is the strain–displacement matrix, \(\Delta \mathbf{u}\) the nodal displacement increment; \(\mathbf{K}_e\) the element tangent stiffness, \(\mathbf{C}\) the consistent material tangent; \(J_0 = \det \mathbf{J}\) the Jacobian determinant and \(w\) the quadrature weight; \(\mathbf{P}_e\) the element internal force vector, \(\sig\) the Cauchy stress tensor in Voigt notation, and \(\sum_{qp}\) denotes summation over quadrature points.

The mass matrix and body-force vector are obtained from the expressions:

\[\mathbf{M}_e = \sum_{qp} \rho\, \mathbf{N}^\mathsf{T} \mathbf{N}\, J_0 w, \qquad \mathbf{P}_e^{(b)} = \sum_{qp} \mathbf{N}^\mathsf{T} \mathbf{f}\, J_0 w.\]

where \(\mathbf{M}_e\) is the consistent mass matrix, \(\rho\) the mass density, \(\mathbf{N}\) the shape-function matrix of the displacement field, and \(\mathbf{f}\) the body-force vector per unit volume.

Surface tractions and pressures on a boundary face \(\Gamma_e\) are integrated as

\[\mathbf{P}_e^{(t)} = \int_{\Gamma_e} \mathbf{N}^\mathsf{T} \, \mathbf{t} \, \mathrm{d}\Gamma, \qquad \mathbf{P}_e^{(p)} = - \int_{\Gamma_e} p \, \mathbf{N}^\mathsf{T} \, \mathbf{n} \, \mathrm{d}\Gamma.\]

where \(\Gamma_e\) is the element boundary, \(\mathbf{n}\) its outward unit normal, \(\mathbf{t}\) the surface traction vector, \(p\) the pressure magnitude, and \(\mathbf{N}\) the shape-function matrix. The minus sign for pressure follows an outward-normal convention.

For 2D and 1D problems, the integration measure \(J_0 w\) is scaled by thickness \(t\) or cross-sectional area \(A\), respectively: \(J_0 w \leftarrow J_0 w\, t\) (2D), \(J_0 w \leftarrow J_0 w\, A\) (1D). Here, \(t\) denotes thickness and \(A\) denotes cross-sectional area.

Constitutive updates are carried out per quadrature point using a Marmot material model:

• 3D solid

• Plane stress

• Plane strain

• 1D uniaxial stress

Each quadrature point stores stress, strain and a material state vector; the element accumulates the history incrementally and may suggest a reduced time step if required by the material.

Implementation

template<int nDim, int nNodes>
class DisplacementFiniteElement : public MarmotElement, public MarmotGeometryElement<nDim, nNodes>

Inheritance diagram for Marmot::Elements::DisplacementFiniteElement:

Collaboration diagram for Marmot::Elements::DisplacementFiniteElement:

Displacement-based finite element template.

Uses linearized kinematics (small strains) and supports 1D, 2D, and 3D formulations via a section assumption. Holds quadrature-point state and delegates constitutive updates to Marmot materials while assembling element residuals, tangents, and mass matrices.

Template Parameters:

• nDim – Number of spatial dimensions (1, 2, or 3).

• nNodes – Number of element nodes.

Public Types

enum SectionType

Kinematic section assumption used by the element.

Controls which constitutive call is performed at quadrature points.

• UniaxialStress

• PlaneStress

• PlaneStrain

• Solid

Values:

enumerator UniaxialStress

enumerator PlaneStress

enumerator PlaneStrain

enumerator Solid

using ParentGeometryElement = MarmotGeometryElement<nDim, nNodes>

using JacobianSized = typename ParentGeometryElement::JacobianSized

using dNdXiSized = typename ParentGeometryElement::dNdXiSized

using BSized = typename ParentGeometryElement::BSized

using XiSized = typename ParentGeometryElement::XiSized

using RhsSized = Matrix<double, sizeLoadVector, 1>

using KeSizedMatrix = Matrix<double, sizeLoadVector, sizeLoadVector>

using CSized = Matrix<double, ParentGeometryElement::voigtSize, ParentGeometryElement::voigtSize>

using Voigt = Matrix<double, ParentGeometryElement::voigtSize, 1>

Public Functions

DisplacementFiniteElement(int elementID, FiniteElement::Quadrature::IntegrationTypes integrationType, SectionType sectionType)

Construct element with ID, quadrature rule and section assumption.

Parameters:

• elementID – Unique element label.

• integrationType – Integration (quadrature) rule.

• sectionType – Section assumption (1D/2D/3D).

virtual int getNumberOfRequiredStateVars()

Total number of required state variables for this element (sum over all quadrature points).

virtual std::vector<std::vector<std::string>> getNodeFields()

Node-level fields exposed by the element. Returns [“displacement”] for each node.

virtual std::vector<int> getDofIndicesPermutationPattern()

Permutation pattern from local DOF ordering to solver ordering (identity by default).

inline virtual int getNNodes()

Number of nodes of this element type.

inline virtual int getNSpatialDimensions()

Number of spatial dimensions.

inline virtual int getNDofPerElement()

Number of degrees of freedom per element (nNodes * nDim).

inline virtual std::string getElementShape()

Geometric shape of the element (as reported by the parent geometry element).

virtual void assignStateVars(double *stateVars, int nStateVars)

Map the provided element state vector to all quadrature points.

virtual void assignProperty(const ElementProperties &marmotElementProperty)

Assign element properties (e.g., thickness in 2D, area in 1D).

virtual void assignProperty(const MarmotMaterialSection &marmotElementProperty)

Assign material section and instantiate per-quadrature-point materials.

virtual void assignNodeCoordinates(const double *coordinates)

Provide nodal coordinates to the parent geometry element.

virtual void initializeYourself()

Precompute geometry-related quantities at quadrature points (B, detJ, J0xW).

virtual void setInitialConditions(StateTypes state, const double *values)

Initialize state or materials.

Parameters:

• state – MarmotMaterialInitialization, GeostaticStress or MarmotMaterialStateVars.

• values – For GeostaticStress: [sigmaY(z1), y1, sigmaY(z2), y2, kx, kz].

virtual void computeDistributedLoad(MarmotElement::DistributedLoadTypes loadType, double *P, double *K, const int elementFace, const double *load, const double *QTotal, const double *time, double dT)

Assemble distributed surface loads on a boundary face.

Pressure and traction contributions are integrated on the boundary \(\Gamma_e\):

\[\mathbf{P}_e^{(p)} = - \int_{\Gamma_e} p\, \mathbf{N}^\mathsf{T} \mathbf{n}\, \mathrm{d}\Gamma,\qquad \mathbf{P}_e^{(t)} = \int_{\Gamma_e} \mathbf{N}^\mathsf{T} \mathbf{t}\, \mathrm{d}\Gamma. \]

Parameters:

• loadType – Pressure or SurfaceTraction.

• P – Element RHS contribution (accumulated).

• K – Optional stiffness contribution (unused).

• elementFace – Boundary face index.

• load – Pressure magnitude or traction vector (size nDim).

• QTotal – Total DOF vector (unused).

• time – Current time data forwarded to materials.

• dT – Time increment.

virtual void computeBodyForce(double *P, double *K, const double *load, const double *QTotal, const double *time, double dT)

Assemble body force contribution.

Integrates \(\mathbf{P}_e^{(b)} = \int_{\Omega_e} \mathbf{N}^\mathsf{T} \mathbf{f}\, \mathrm{d}\Omega\).

virtual void computeYourself(const double *QTotal, const double *dQ, double *Pe, double *Ke, const double *time, double dT, double &pNewdT)

Compute internal force and consistent tangent stiffness.

Uses the small-strain relation \(\Delta\boldsymbol{\varepsilon}=\mathbf{B}\,\Delta\mathbf{u}\) and integrates

\[\mathbf{K}_e = \sum_{qp} \mathbf{B}^\mathsf{T} \mathbf{C} \mathbf{B}\, J_0 w,\qquad \mathbf{P}_e = \sum_{qp} \mathbf{B}^\mathsf{T} \boldsymbol{\sigma}\, J_0 w. \]

If pNewdT<1, the routine returns early to signal time step reduction.

Parameters:

• QTotal – Total displacement vector.

• dQ – Incremental displacement.

• Pe – Internal force vector (accumulated).

• Ke – Tangent stiffness matrix (accumulated).

• time – Time data forwarded to materials.

• dT – Time increment.

• pNewdT – Suggested scaling of dT by the material; if reduced (<1), the routine returns early.

virtual void computeConsistentInertia(double *M)

Compute consistent mass matrix using material density.

\(\mathbf{M}_e = \sum_{qp} \rho\, \mathbf{N}^\mathsf{T}\mathbf{N}\, J_0 w\).

virtual void computeLumpedInertia(double *M)

Compute lumped mass vector via row-sum of consistent mass.

\(\mathbf{m}_e = \mathrm{rowsum}(\mathbf{M}_e)\).

inline virtual StateView getStateView(const std::string &stateName, int qpNumber)

Access a named state view at a quadrature point.

Note

Using “sdv” returns the raw material state vector and is deprecated.

virtual std::vector<double> getCoordinatesAtCenter()

Get physical coordinates at the element center.

virtual std::vector<std::vector<double>> getCoordinatesAtQuadraturePoints()

Get physical coordinates at each quadrature point.

virtual int getNumberOfQuadraturePoints()

Number of quadrature points of this element.

Public Members

Map<const VectorXd> elementProperties

Element-level properties (e.g., thickness for 2D, area for 1D).

const int elLabel

Element label (ID) used for logging and material creation.

const SectionType sectionType

Section assumption applied by this element instance.

std::vector<QuadraturePoint> qps

Quadrature points owned by the element (one per integration point).

Public Static Attributes

static int sizeLoadVector = nNodes * nDim

static int nCoordinates = nNodes * nDim

struct QuadraturePoint

Data and state associated with a quadrature point.

Holds parent coordinates, integration weight, Jacobian determinant, kinematic (strain-displacement) B-matrix, and a material instance with managed state variables.

Public Functions

inline int getNumberOfRequiredStateVarsQuadraturePointOnly()

inline int getNumberOfRequiredStateVars()

inline void assignStateVars(double *stateVars, int nStateVars)

inline QuadraturePoint(XiSized xi, double weight)

Public Members

const XiSized xi

const double weight

double detJ

double J0xW

BSized B

std::unique_ptr<QPStateVarManager> managedStateVars

std::unique_ptr<MarmotMaterialHypoElastic> material

class QPStateVarManager : public MarmotStateVarVectorManager

Inheritance diagram for Marmot::Elements::DisplacementFiniteElement::QuadraturePoint::QPStateVarManager:

Collaboration diagram for Marmot::Elements::DisplacementFiniteElement::QuadraturePoint::QPStateVarManager:

Manager for per-quadrature-point state variables.

Provides named accessors to stress \(\sig\), strain \(\eps\) and the material state vector. The layout is [stress(6), strain(6), begin of material state(…)] in 3D Voigt notation.

Public Functions

inline QPStateVarManager(double *theStateVarVector, int nStateVars)

Public Members

mVector6d stress

mVector6d strain

Eigen::Map<Eigen::VectorXd> materialStateVars

Public Static Functions

static inline int getNumberOfRequiredStateVarsQuadraturePointOnly()

Private Static Attributes

static const auto layout = makeLayout( {{ .name = "stress", .length = 6 },{ .name = "strain", .length = 6 },{ .name = "begin ofmaterialstate", .length = 0 },} )