Internal Functions

These are additional functions that do not clearly fit into another category.

typedef PetscErrorCode RatelTSMonitorFunction(TS, PetscInt, PetscReal, Vec, void*)

Monitor routine signature.

PetscClassId RATEL_CLASSID

PetscLogStage RATEL_Setup

PetscLogEvent RATEL_DMSetupByOrder

PetscLogEvent RATEL_DMSetupSolver

PetscLogEvent RATEL_Diagnostics

PetscLogEvent RATEL_Diagnostics_CeedOp

PetscLogEvent RATEL_Residual

PetscLogEvent RATEL_Residual_CeedOp

PetscLogEvent RATEL_Jacobian

PetscLogEvent RATEL_Jacobian_CeedOp

static inline CeedMemType MemTypePetscToCeed(PetscMemType mem_type)

Translate PetscMemType to CeedMemType.

Parameters:

• mem_type – [in] PetscMemType

Returns:

Equivalent CeedMemType

static inline CeedInt GetCeedQuadratureSize(CeedEvalMode eval_mode, CeedInt dim, CeedInt num_components)

Return the quadrature size from the eval mode, dimension, and number of components.

Parameters:

• eval_mode – [in] The basis evaluation mode

• dim – [in] The basis dimension

• num_components – [in] The basis number of components

Returns:

The maximum of the two integers

int MMSError_CEED_ScalarBPs(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute error from true solution for CEED scalar BPs MMS.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSCEEDBPsParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - quadrature point coordinates

  • 2 - computed solution

• out – [out] Output array

  • 0 - error from true solution

Returns:

An error code: 0 - success, otherwise - failure

int MMSError_CEED_VectorBPs(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute error from true solution for CEED vector BPs MMS.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSCEEDBPsParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - quadrature point coordinates

  • 2 - computed solution

• out – [out] Output array

  • 0 - error from true solution

Returns:

An error code: 0 - success, otherwise - failure

int MMSError_Linear(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute error from true solution for linear elasticity MMS.

Uses an MMS adapted from doi:10.1016/j.compstruc.2019.106175

Parameters:

• ctx – [in] QFunction context, unused

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - quadrature point coordinates

  • 2 - computed solution

• out – [out] Output array

  • 0 - error from true solution

Returns:

An error code: 0 - success, otherwise - failure

int BodyForce(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute forcing term with user-specified body force.

Parameters:

• ctx – [in] QFunction context, holding RatelForcingBodyParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

• out – [out] Output array

  • 0 - forcing term

Returns:

An error code: 0 - success, otherwise - failure

int BodyForceEnergy(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute external energy cause by forcing term with user-specified body force.

Parameters:

• ctx – [in] QFunction context, holding RatelForcingBodyParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - displacement u

• out – [out] Output array

  • 0 - energy due to body force

Returns:

An error code: 0 - success, otherwise - failure

int MMSForce_CEED_BP1(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute forcing term for CEED scalar BP1 MMS.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSCEEDBPsParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - quadrature point coordinates

• out – [out] Output array

  • 0 - forcing term

Returns:

An error code: 0 - success, otherwise - failure

int MMSForce_CEED_BP2(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute forcing term for CEED vector BP2 MMS.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSCEEDBPsParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - quadrature point coordinates

• out – [out] Output array

  • 0 - forcing term

Returns:

An error code: 0 - success, otherwise - failure

int MMSForce_CEED_BP3(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute forcing term for CEED scalar BP3 MMS.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSCEEDBPsParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - quadrature point coordinates

• out – [out] Output array

  • 0 - forcing term

Returns:

An error code: 0 - success, otherwise - failure

int MMSForce_CEED_BP4(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute forcing term for CEED vector BP4 MMS.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSCEEDBPsParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - quadrature point coordinates

• out – [out] Output array

  • 0 - forcing term

Returns:

An error code: 0 - success, otherwise - failure

int MMSForce_Linear(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute forcing term for linear elasticity MMS.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSLinearElasticityParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - quadrature point coordinates

• out – [out] Output array

  • 0 - forcing term

Returns:

An error code: 0 - success, otherwise - failure

int MMSForceEnergy_Linear(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute external energy cause by forcing term for linear elasticity MMS.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSLinearElasticityParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - quadrature point coordinates

  • 2 - displacement u

• out – [out] Output array

  • 0 - energy due to mms force

Returns:

An error code: 0 - success, otherwise - failure

int MMSForce_MixedLinear(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute forcing term for mixed linear elasticity MMS.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSMixedLinearElasticityParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - quadrature point coordinates

• out – [out] Output array

  • 0 - forcing term, displacement field

Returns:

An error code: 0 - success, otherwise - failure

int MMSForceEnergy_MixedLinear(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute external energy cause by forcing term for mixed linear elasticity MMS.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSMixedLinearElasticityParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - quadrature point coordinates

  • 2 - displacement u

• out – [out] Output array

  • 0 - energy due to mms force

Returns:

An error code: 0 - success, otherwise - failure

int MPMBodyForce(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute forcing term with user-specified body force and density defined at points.

Parameters:

• ctx – [in] QFunction context, holding RatelForcingBodyParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - state data

  • 2 - density

• out – [out] Output array

  • 0 - forcing term

Returns:

An error code: 0 - success, otherwise - failure

int MMSForce_PoroElasticityLinear(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute forcing term for linear poroelasticity MMS.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSPoroElasticityLinearParams

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - qdata

  • 1 - quadrature point coordinates

• out – [out] Output array

  • 0 - forcing term, displacement field

  • 1 - forcing term, pressure field

Returns:

An error code: 0 - success, otherwise - failure

int Mass(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Apply mass operator.

Parameters:

• ctx – [in] QFunction context, holding the number of components

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - volumetric qdata

  • 1 - u

• out – [out] Output array

  • 0 - v

Returns:

An error code: 0 - success, otherwise - failure

int RatelScalarBPsMMSTrueSolution(void *ctx, const CeedScalar coords[3], CeedScalar *true_solution)

Compute true solution for scalar BPs problems.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSCEEDBPsParams

• coords – [in] Coordinate array

• true_solution – [out] True solution

Returns:

An error code: 0 - success, otherwise - failure

int RatelScalarBPsPoissonForcing(void *ctx, const CeedScalar coords[3], CeedScalar *poisson_forcing)

Compute forcing term for Poisson problems.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSCEEDBPsParams

• coords – [in] Coordinate array

• poisson_forcing – [out] Forcing term needed for Poisson problem

Returns:

An error code: 0 - success, otherwise - failure

int RatelVectorBPsMMSTrueSolution(void *ctx, const CeedScalar coords[3], CeedScalar true_solution[3])

Compute true solution for vector BPs problems.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSCEEDBPsParams

• coords – [in] Coordinate array

• true_solution – [out] True solution

Returns:

An error code: 0 - success, otherwise - failure

int RatelVectorBPsPoissonForcing(void *ctx, const CeedScalar coords[3], CeedScalar poisson_forcing[3])

Compute forcing term Poisson problems.

Parameters:

• ctx – [in] QFunction context, holding RatelMMSCEEDBPsParams

• coords – [in] Coordinate array

• poisson_forcing – [out] Forcing term needed for Poisson problem

Returns:

An error code: 0 - success, otherwise - failure

int RatelLinearElasticityMMSTrueSolution(void *ctx, const CeedScalar coords[3], CeedScalar true_solution[3])

Compute true solution for linear elasticity MMS.

Uses an MMS adapted from doi:10.1016/j.compstruc.2019.106175

Parameters:

• ctx – [in] QFunction context, holding RatelMMSLinearElasticityParams

• coords – [in] Coordinate array

• true_solution – [out] True solution

Returns:

An error code: 0 - success, otherwise - failure

int RatelLinearElasticityMMSForcing(void *ctx, const CeedScalar coords[3], CeedScalar mms_forcing[3])

Compute forcing term for linear elasticity MMS.

Uses an MMS adapted from doi:10.1016/j.compstruc.2019.106175

Parameters:

• ctx – [in] QFunction context, holding RatelMMSLinearElasticityParams

• coords – [in] Coordinate array

• mms_forcing – [out] Forcing term

Returns:

An error code: 0 - success, otherwise - failure

int RatelMixedLinearElasticityMMSTrueSolution(void *ctx, const CeedScalar coords[3], CeedScalar true_solution_u[3], CeedScalar *true_solution_p)

Compute true solution for mixed linear elasticity MMS.

Uses an MMS adapted from doi:10.1016/j.compstruc.2019.106175

Parameters:

• ctx – [in] QFunction context, holding RatelMMSMixedLinearElasticityParams

• coords – [in] Coordinate array

• true_solution_u – [out] True solution, displacement

• true_solution_p – [out] True solution, pressure

Returns:

An error code: 0 - success, otherwise - failure

int RatelMixedLinearElasticityMMSForcing(void *ctx, const CeedScalar coords[3], CeedScalar mms_forcing_u[3])

Compute forcing term for mixed linear elasticity MMS.

Uses an MMS adapted from doi:10.1016/j.compstruc.2019.106175

Parameters:

• ctx – [in] QFunction context, holding RatelMMSMixedLinearElasticityParams

• coords – [in] Coordinate array

• mms_forcing_u – [out] Forcing term for displacement

Returns:

An error code: 0 - success, otherwise - failure

int RatelPoroElasticityMMSTrueSolution(void *ctx, const CeedScalar coords[3], CeedScalar true_solution_u[3], CeedScalar *true_solution_p)

Compute true solution for linear poroelasticity MMS.

Uses an MMS adapted from doi:10.1016/j.compstruc.2019.106175

Parameters:

• ctx – [in] QFunction context, holding RatelMMSPoroElasticityLinearParams

• coords – [in] Coordinate array

• true_solution_u – [out] True solution, displacement

• true_solution_p – [out] True solution, pressure

Returns:

An error code: 0 - success, otherwise - failure

int RatelPoroElasticityMMSForcing(void *ctx, const CeedScalar coords[3], CeedScalar mms_forcing_u[3])

Compute displacement forcing term for linear poroelasticity MMS.

Uses an MMS adapted from doi:10.1016/j.compstruc.2019.106175

Parameters:

• ctx – [in] QFunction context, holding RatelMMSPoroElasticityLinearParams

• coords – [in] Coordinate array

• mms_forcing_u – [out] Forcing term for displacement

Returns:

An error code: 0 - success, otherwise - failure

int RatelPoroElasticityMMSGamma(void *ctx, const CeedScalar coords[3], CeedScalar *mms_forcing_p)

Compute pressure forcing term for linear poroelasticity MMS.

Uses an MMS adapted from doi:10.1016/j.compstruc.2019.106175

Parameters:

• ctx – [in] QFunction context, holding RatelMMSPoroElasticityLinearParams

• coords – [in] Coordinate array

• mms_forcing_p – [out] Forcing term for pressure

Returns:

An error code: 0 - success, otherwise - failure

int RatelMigratePoints(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Apply mass operator.

Parameters:

• ctx – [in] QFunction context, not used

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - delta u

• out – [out] Output array

  • 0 - point positions (x)

Returns:

An error code: 0 - success, otherwise - failure

int ScaleLumpedDualDiagnosticTerms(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Apply inverse of lumped mass operator to nodal values.

Parameters:

• ctx – [in] QFunction context, holding the number of components

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - element restriction multiplicity

  • 1 - composite operator nodal multiplicity

  • 2 - dual diagnostics: nodal volume, other values

• out – [out] Output array

  • 0 - lumped dual diagnostics: nodal volume, other values divided by nodal volume

Returns:

An error code: 0 - success, otherwise - failure

int ScaledMass(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Apply scaled mass operator.

Parameters:

• ctx – [in] QFunction context, holding RatelScaledMassParams struct

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - volumetric qdata

  • 1 - u

• out – [out] Output array

  • 0 - v

Returns:

An error code: 0 - success, otherwise - failure

int ComputeCentroid(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)

Compute face centroid.

Parameters:

• ctx – [in] QFunction context, unused

• Q – [in] Number of quadrature points

• in – [in] Input arrays

  • 0 - volumetric qdata

  • 1 - quadrature point coordinates

• out – [out] Output array

  • 0 - centroid

Returns:

An error code: 0 - success, otherwise - failure

bool RatelCoordinatesInBoundingBox(const CeedScalar x[3], const CeedScalar min[3], const CeedScalar max[3])

Check if a point is in a bounding box.

Parameters:

• x – [in] Point to check

• min – [in] Minimum values of the bounding box

• max – [in] Maximum values of the bounding box

CeedScalar RatelLog1pSeries(CeedScalar x)

Series approximation of log1p()

log1p() is not vectorized in libc. The series expansion up to the sixth term is accurate to 1e-10 relative error in the range 2/3 < J < 3/2. This is accurate enough for many hyperelastic applications, but perhaps not for materials like foams that may encounter volumetric strains on the order of 1/5 < J < 5.

Parameters:

• x – [in] Scalar value

Returns:

A scalar value: log1p(x)

CeedScalar RatelExpm1Series(CeedScalar x)

Series approximation of expm1()

Parameters:

• x – [in] Scalar value

Returns:

A scalar value: expm1(x)

CeedScalar RatelAtanSeries(CeedScalar x)

Approximation of atan(x) with max error 1.95815585968262e-14 on [-1,1].

Computed via Remez algorithm with 31 degree polynomial using github.com/DKenefake/OptimalPoly

Parameters:

• x – [in] Input value

Returns:

An error code: 0 - success, otherwise - failure

CeedInt RatelSign(CeedScalar x)

Compute sign of scalar value.

Parameters:

• x – [in] Scalar value

Returns:

An integer: -1 if x is negative and 1 if x is positive

int RatelQdataPack(CeedInt Q, CeedInt i, CeedScalar wdetJ, const CeedScalar local[3][3], CeedScalar stored[10][CEED_Q_VLA])

Unpack qdata at quadrature point.

Parameters:

• Q – [in] Number of quadrature points

• i – [in] Current quadrature point

• wdetJ – [in] Quadrature weight

• local – [in] Qdata for quadrature point i

• stored – [out] All qdata

Returns:

An error code: 0 - success, otherwise - failure

int RatelQdataUnpack(CeedInt Q, CeedInt i, const CeedScalar stored[10][CEED_Q_VLA], CeedScalar local[3][3])

Unpack qdata at quadrature point.

Parameters:

• Q – [in] Number of quadrature points

• i – [in] Current quadrature point

• stored – [in] All qdata

• local – [out] Qdata for quadrature point i

Returns:

An error code: 0 - success, otherwise - failure

int RatelStoredValuesPack(CeedInt Q, CeedInt i, CeedInt start, CeedInt num_comp, const CeedScalar *local, CeedScalar *stored)

Pack stored values at quadrature point.

Parameters:

• Q – [in] Number of quadrature points

• i – [in] Current quadrature point

• start – [in] Starting index to store components

• num_comp – [in] Number of components to store

• local – [in] Local values for quadrature point i

• stored – [out] Stored values

Returns:

An error code: 0 - success, otherwise - failure

int RatelStoredValuesUnpack(CeedInt Q, CeedInt i, CeedInt start, CeedInt num_comp, const CeedScalar *stored, CeedScalar *local)

Unpack stored values at quadrature point.

Parameters:

• Q – [in] Number of quadrature points

• i – [in] Current quadrature point

• start – [in] Starting index to store components

• num_comp – [in] Number of components to store

• stored – [in] Stored values

• local – [out] Local values for quadrature point i

Returns:

An error code: 0 - success, otherwise - failure

int RatelGradUnpack(CeedInt Q, CeedInt i, const CeedScalar grad[3][3][CEED_Q_VLA], CeedScalar local[3][3])

Unpack gradient at quadrature point.

Parameters:

• Q – [in] Number of quadrature points

• i – [in] Current quadrature point

• grad – [in] Gradient of u, x (or du in jacobian) w.r.t X, X is ref coordinate [-1, 1]^3

• local – [out] dudX, dxdX (or ddudX in jacobian) for quadrature point i

Returns:

An error code: 0 - success, otherwise - failure

CeedScalar RatelDot3(const CeedScalar a[3], const CeedScalar b[3])

Compute dot product of two length 3 vectors.

Parameters:

• a – [in] First input vector

• b – [in] Second input vector

Returns:

A scalar value: a . b

CeedScalar RatelNorm3(const CeedScalar a[3])

Compute norm of a length 3 vector.

Parameters:

• a – [in] Input vector

Returns:

A scalar value: |a|

int RatelScalarVecMult(CeedScalar alpha, const CeedScalar a[3], CeedScalar b[3])

Compute b = alpha a for a length 3 vector.

Parameters:

• alpha – [in] Scaling factor

• a – [in] Input vector

• b – [out] Output vector

Returns:

An error code: 0 - success, otherwise - failure

int RatelVecVecSubtract(const CeedScalar a[3], const CeedScalar b[3], CeedScalar c[3])

Compute c = a - b for two length 3 vectors.

Parameters:

• a – [in] First input vector

• b – [in] Second input vector

• c – [out] Output vector

Returns:

An error code: 0 - success, otherwise - failure

int RatelVecVecCross(const CeedScalar a[3], const CeedScalar b[3], CeedScalar c[3])

Compute c = a x b for two length 3 vectors.

Parameters:

• a – [in] First input vector

• b – [in] Second input vector

• c – [out] Output vector

Returns:

An error code: 0 - success, otherwise - failure

int RatelVecOuterMult(const CeedScalar a[3], CeedScalar A_sym[6])

Compute outer product of a vector length 3 with it self results in symmetric matrix A_sym = a*a^T

Parameters:

• a – [in] Input vector

• A_sym – [out] Symmetric matrix

Returns:

An error code: 0 - success, otherwise - failure

int RatelVecVecOuterMult(const CeedScalar a[3], const CeedScalar b[3], CeedScalar A[3][3])

Compute outer product of two vectors of length 3 results in a matrix A = a*b^T

Parameters:

• a – [in] Input vector

• b – [in] Input vector

• A – [out] Matrix

Returns:

An error code: 0 - success, otherwise - failure

int RatelScalarMatMultSymmetric(CeedScalar alpha, const CeedScalar A_sym[6], CeedScalar B_sym[6])

Compute B_sym = alpha A_sym

Parameters:

• alpha – [in] Scaling factor

• A_sym – [in] Input 3x3 matrix, in symmetric representation

• B_sym – [out] Output 3x3 matrix, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatVecMult(CeedScalar alpha, const CeedScalar A[3][3], const CeedScalar x[3], CeedScalar b[3])

Compute b = alpha A x for a length 3 vector and a 3x3 matrix.

Parameters:

• alpha – [in] Scaling factor

• A – [in] Input matrix

• x – [in] Input vector

• b – [out] Output vector

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatTransposeVecMult(CeedScalar alpha, const CeedScalar A[3][3], const CeedScalar x[3], CeedScalar b[3])

Compute b = alpha A^T x for a length 3 vector and a 3x3 matrix.

Parameters:

• alpha – [in] Scaling factor

• A – [in] Input matrix

• x – [in] Input vector

• b – [out] Output vector

Returns:

An error code: 0 - success, otherwise - failure

int RatelSymmetricMatPack(const CeedScalar full[3][3], CeedScalar sym[6])

Pack symmetric 3x3 matrix from full to symmetric representation.

Matrix is packed as

\( \begin{vmatrix} s_0 & s_5 & s_4 \\ s_5 & s_1 & s_3 \\ s_4 & s_3 & s_2 \end{vmatrix} \)

where \(s_0, s_1, \dots \) are the indices of the vector in symmetric representation.

Parameters:

• full – [in] Input full 3x3 matrix

• sym – [out] Output matrix as length 6 vector, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelSymmetricMatUnpack(const CeedScalar sym[6], CeedScalar full[3][3])

Unpack symmetric 3x3 matrix from symmetric to full representation.

Matrix is packed as

\( \begin{vmatrix} s_0 & s_5 & s_4 \\ s_5 & s_1 & s_3 \\ s_4 & s_3 & s_2 \end{vmatrix} \)

where \(s_0, s_1, \dots \) are the indices of the vector in symmetric representation.

Parameters:

• sym – [in] Input matrix as length 6 vector, in symmetric representation

• full – [out] Output full 3x3 matrix

Returns:

An error code: 0 - success, otherwise - failure

void RatelKelvinMandelMatUnpack(const CeedScalar sym[6], CeedScalar full[3][3])

Unpack symmetric tensor into full tensor via Kelvin-Mandel normalization.

Matrix is packed as

\( \begin{vmatrix} s_0 & s_5/ \sqrt{2} & s_4/ \sqrt{2} \\ s_5/ \sqrt{2} & s_1 & s_3/ \sqrt{2} \\ s_4/ \sqrt{2} & s_3/ \sqrt{2} & s_2/ \sqrt{2} \end{vmatrix} \)

where \(s_0, s_1, \dots \) are the indices of the vector in symmetric representation.

Parameters:

• sym – [in] Input matrix as length 6 vector, in symmetric representation

• full – [out] Output full 3x3 matrix

Returns:

An error code: 0 - success, otherwise - failure

void RatelKelvinMandelMatPack(const CeedScalar full[3][3], CeedScalar sym[6])

Pack full tensor into symmetric tensor of length 6 via Kelvin-Mandel normalization.

Matrix is packed as

\( \begin{vmatrix} s_0 & s_5/ \sqrt{2} & s_4/ \sqrt{2} \\ s_5/ \sqrt{2} & s_1 & s_3/ \sqrt{2} \\ s_4/ \sqrt{2} & s_3/ \sqrt{2} & s_2/ \sqrt{2} \end{vmatrix} \)

where \(s_0, s_1, \dots \) are the indices of the vector in symmetric representation.

Parameters:

• full – [in] Input full 3x3 matrix

• sym – [out] Output matrix as length 6 vector, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

CeedScalar RatelMatTrace(const CeedScalar A[3][3])

Compute the trace of a 3x3 matrix.

Parameters:

• A – [in] Input 3x3 matrix

Returns:

A scalar value: trace(A)

CeedScalar RatelMatTraceSymmetric(const CeedScalar A_sym[6])

Compute the trace of a 3x3 matrix.

Parameters:

• A_sym – [in] Input 3x3 matrix, in symmetric representation

Returns:

A scalar value: trace(A)

int RatelMatMatAdd(CeedScalar alpha, const CeedScalar A[3][3], CeedScalar beta, const CeedScalar B[3][3], CeedScalar C[3][3])

Compute C = alpha A + beta B for 3x3 matrices.

Parameters:

• alpha – [in] First scaling factor

• A – [in] First input matrix

• beta – [in] Second scaling factor

• B – [in] Second input matrix

• C – [out] Output matrix

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatCopy(const CeedScalar A[3][3], CeedScalar B[3][3])

Copy B = A for 3x3 matrices.

Parameters:

• A – [in] Input matrix

• B – [out] Output matrix

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatMatAddSymmetric(CeedScalar alpha, const CeedScalar A_sym[6], CeedScalar beta, const CeedScalar B_sym[6], CeedScalar C_sym[6])

Compute C = alpha A + beta B for 3x3 matrices, in symmetric representation.

Parameters:

• alpha – [in] First scaling factor

• A_sym – [in] First input matrix, in symmetric representation

• beta – [in] Second scaling factor

• B_sym – [in] Second input matrix, in symmetric representation

• C_sym – [out] Output matrix, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatMatMatAddSymmetric(CeedScalar alpha, const CeedScalar A_sym[6], CeedScalar beta, const CeedScalar B_sym[6], CeedScalar gamma, const CeedScalar C_sym[6], CeedScalar D_sym[6])

Compute D = alpha A + beta B + gamma C for 3x3 matrices, in symmetric representation.

Parameters:

• alpha – [in] First scaling factor

• A_sym – [in] First input matrix, in symmetric representation

• beta – [in] Second scaling factor

• B_sym – [in] Second input matrix, in symmetric representation

• gamma – [in] Third scaling factor

• C_sym – [in] Third input matrix, in symmetric representation

• D_sym – [out] Output matrix, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatMatMatAdd(CeedScalar alpha, const CeedScalar A[3][3], CeedScalar beta, const CeedScalar B[3][3], CeedScalar gamma, const CeedScalar C[3][3], CeedScalar D[3][3])

Compute D = alpha A + beta B + gamma C for 3x3 matrices.

Parameters:

• alpha – [in] First scaling factor

• A – [in] First input matrix

• beta – [in] Second scaling factor

• B – [in] Second input matrix

• gamma – [in] Third scaling factor

• C – [in] Third input matrix

• D – [out] Output matrix

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatDeviatoricSymmetric(CeedScalar trace_A, const CeedScalar A_sym[6], CeedScalar A_sym_dev[6])

Compute deviatoric part of symmetric matrix; A_sym_dev = A_sym - (1/3) trace(A) I for 3x3 matrix, in symmetric representation.

Parameters:

• trace_A – [in] Trace of input matrix

• A_sym – [in] Input matrix, in symmetric representation

• A_sym_dev – [out] Deviatoric part of input matrix, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatMatMult(CeedScalar alpha, const CeedScalar A[3][3], const CeedScalar B[3][3], CeedScalar C[3][3])

Compute C = alpha A * B for 3x3 matrices.

Parameters:

• alpha – [in] Scaling factor

• A – [in] First input matrix

• B – [in] Second input matrix

• C – [out] Output matrix

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatTransposeMatMult(CeedScalar alpha, const CeedScalar A[3][3], const CeedScalar B[3][3], CeedScalar C[3][3])

Compute C = alpha A^T * B for 3x3 matrices.

Parameters:

• alpha – [in] Scaling factor

• A – [in] First input matrix

• B – [in] Second input matrix

• C – [out] Output matrix

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatMatTransposeMult(CeedScalar alpha, const CeedScalar A[3][3], const CeedScalar B[3][3], CeedScalar C[3][3])

Compute C = alpha A * B^T for 3x3 matrices.

Parameters:

• alpha – [in] Scaling factor

• A – [in] First input matrix

• B – [in] Second input matrix

• C – [out] Output matrix

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatMatMultAtQuadraturePoint(CeedInt Q, CeedInt i, CeedScalar alpha, const CeedScalar A[3][3], const CeedScalar B[3][3], CeedScalar C[3][3][CEED_Q_VLA])

Compute C = alpha A * B for 3x3 matrices at a quadrature point.

Parameters:

• Q – [in] Number of quadrature points

• i – [in] Current quadrature point

• alpha – [in] Scaling factor

• A – [in] First input matrix

• B – [in] Second input matrix

• C – [out] Output matrix, stored for quadrature point i

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatMatTransposeMultAtQuadraturePoint(CeedInt Q, CeedInt i, CeedScalar alpha, const CeedScalar A[3][3], const CeedScalar B[3][3], CeedScalar C[3][3][CEED_Q_VLA])

Compute C = alpha A * B^T for 3x3 matrices at a quadrature point.

Parameters:

• Q – [in] Number of quadrature points

• i – [in] Current quadrature point

• alpha – [in] Scaling factor

• A – [in] First input matrix

• B – [in] Second input matrix

• C – [out] Output matrix, stored for quadrature point i

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatVecMultAtQuadraturePoint(CeedInt Q, CeedInt i, CeedScalar alpha, const CeedScalar A[3][3], const CeedScalar x[3], CeedScalar b[3][CEED_Q_VLA])

Compute b = alpha A x for a length 3 vector and a 3x3 matrix at a quadrature point.

Parameters:

• Q – [in] Number of quadrature points

• i – [in] Current quadrature point

• alpha – [in] Scaling factor

• A – [in] First input matrix

• x – [in] First input vector

• b – [out] Output vector, stored for quadrature point i

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatMatMultPlusMatMatMult(const CeedScalar A[3][3], const CeedScalar B[3][3], const CeedScalar C[3][3], const CeedScalar D[3][3], CeedScalar E[3][3])

Compute E = A * B + C * D for 3x3 matrices.

Parameters:

• A – [in] First input matrix

• B – [in] Second input matrix

• C – [in] Third input matrix

• D – [in] Fourth input matrix

• E – [out] Output matrix

Returns:

An error code: 0 - success, otherwise - failure

CeedScalar RatelMatMatContract(CeedScalar alpha, const CeedScalar A[3][3], const CeedScalar B[3][3])

Compute c = alpha A : B for 3x3 matrices.

Parameters:

• alpha – [in] Scaling factor

• A – [in] First input matrix

• B – [in] Second input matrix

Returns:

A scalar value: alpha A : B

CeedScalar RatelMatMatContractSymmetric(CeedScalar alpha, const CeedScalar A_sym[6], const CeedScalar B_sym[6])

Compute c = alpha A_sym : B_sym for 3x3 matrices, in symmetric representation.

Parameters:

• alpha – [in] Scaling factor

• A_sym – [in] First input matrix, in symmetric representation

• B_sym – [in] Second input matrix, in symmetric representation

Returns:

A scalar value: alpha A_sym : B_sym

CeedScalar RatelMatDetA(const CeedScalar A[3][3])

Compute det(A) for a 3x3 matrix.

Parameters:

• A – [in] Input matrix

Returns:

A scalar value: det(A)

CeedScalar RatelMatDetAM1(const CeedScalar A[3][3])

Compute det(A) - 1 for a 3x3 matrix.

Parameters:

• A – [in] Input matrix

Returns:

A scalar value: det(A) - 1

CeedScalar RatelMatDetAM1Symmetric(const CeedScalar A_sym[6])

Compute det(A_sym) - 1 for a 3x3 matrix, in symmetric representation.

Parameters:

• A_sym – [in] Input matrix, in symmetric representation

Returns:

A scalar value: det(A_sym) - 1

int RatelMatMatMultSymmetric(CeedScalar alpha, const CeedScalar A_sym[6], const CeedScalar B_sym[6], CeedScalar C_sym[6])

Compute C_sym = alpha A_sym * B_sym for 3x3 matrices, in symmetric representation.

Note C_sym will not be symmetric in general. Only use this function when A_sym and B_sym have same eigenvectors.

Parameters:

• alpha – [in] Scaling factor

• A_sym – [in] First input matrix, in symmetric representation

• B_sym – [in] Second input matrix, in symmetric representation

• C_sym – [out] Output matrix, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

CeedScalar RatelMatNorm(const CeedScalar A[3][3])

Compute the Frobenius norm ||A|| = sqrt(sum_ij |a_ij|^2) for a symmetric 3x3 matrix.

Parameters:

• A – [in] Input matrix

Returns:

A scalar value: ||A||

int RatelMatMatMatMultSymmetric(CeedScalar alpha, const CeedScalar A_sym[6], const CeedScalar B_sym[6], CeedScalar C_sym[6])

Compute C_sym = A_sym * B_sym * A_sym for 3x3 matrices, in symmetric representation.

Parameters:

• alpha – [in] Scaling factor

• A_sym – [in] First input matrix, in symmetric representation

• B_sym – [in] Second input matrix, in symmetric representation

• C_sym – [out] Output matrix, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatInverse(const CeedScalar A[3][3], CeedScalar det_A, CeedScalar A_inv[3][3])

Compute A^-1 for a 3x3 matrix.

Parameters:

• A – [in] Input matrix

• det_A – [in] Determinant of A

• A_inv – [out] Output matrix inverse

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatInverseSymmetric(const CeedScalar A[3][3], CeedScalar det_A, CeedScalar A_inv[6])

Compute A^-1 for a symmetric 3x3 matrix.

Parameters:

• A – [in] Input matrix

• det_A – [in] Determinant of A

• A_inv – [out] Output matrix inverse, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelScaledMassApplyAtQuadraturePoint(CeedInt Q, CeedInt i, CeedInt num_comp, CeedScalar alpha, const CeedScalar *u, CeedScalar *v)

Apply mass matrix for all components at a quadrature point.

Parameters:

• Q – [in] Number of quadrature points

• i – [in] Current quadrature point

• num_comp – [in] Number of components in u

• alpha – [in] Scaling factor

• u – [in] Input array

• v – [out] Output array

Returns:

An error code: 0 - success, otherwise - failure

int RatelCInverse(const CeedScalar E_sym[6], CeedScalar Jm1, CeedScalar C_inv_sym[6])

Compute inverse of right Cauchy-Green tensor.

Parameters:

• E_sym – [in] Strain tensor, in symmetric representation

• Jm1 – [in] Determinant of deformation gradient(F) - 1; J = det(F)

• C_inv_sym – [out] Inverse of right Cauchy-Green tensor, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelCInverse_fwd(const CeedScalar C_inv_sym[6], const CeedScalar dE_sym[6], CeedScalar dC_inv_sym[6])

Compute dC_inv_sym = -2 C_inv_sym * dE_sym * C_inv_sym

Parameters:

• C_inv_sym – [in] Inverse of right Cauchy-Green tensor, in symmetric representation

• dE_sym – [in] Incremental strain energy tensor, in symmetric representation

• dC_inv_sym – [out] Inverse of incremental right Cauchy-Green tensor, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelLinearStrain(const CeedScalar grad_u[3][3], CeedScalar e_sym[6])

Compute linear strain e = (grad_u + grad_u^T) / 2

Parameters:

• grad_u – [in] Gradient of u

• e_sym – [out] Linear strain, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelLinearStrain_fwd(const CeedScalar grad_du[3][3], CeedScalar de_sym[6])

Compute incremental linear strain de = (grad_du + grad_du^T) / 2

Parameters:

• grad_du – [in] Gradient of incremental change to u

• de_sym – [out] Incremental linear strain, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelGreenEulerStrain(const CeedScalar grad_u[3][3], CeedScalar e_sym[6])

Compute Green Euler strain e = (b - I) / 2 = (grad_u + grad_u^T + grad_u * grad_u^T) / 2 for current configuration.

Parameters:

• grad_u – [in] Gradient of u

• e_sym – [out] Green Euler strain, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelGreenEulerStrain_fwd(const CeedScalar grad_du[3][3], const CeedScalar b[3][3], CeedScalar de_sym[6])

Compute incremental Green Euler strain de = db / 2 = (grad_du b + b grad_du^T) / 2 for current configuration.

grad_du = d(du)/dx where x is physical coordinate in current configuration.

Parameters:

• grad_du – [in] Gradient of incremental change to u

• b – [in] Left Cauchy-Green deformation tensor

• de_sym – [out] Incremental Green Euler strain, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelGreenLagrangeStrain(const CeedScalar grad_u[3][3], CeedScalar E_sym[6])

Compute Green Lagrange strain E = (C - I)/2 = (grad_u + grad_u^T + grad_u^T * grad_u)/2 for initial configuration.

Parameters:

• grad_u – [in] Gradient of u

• E_sym – [out] Green Lagrange strain, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelGreenLagrangeStrain_fwd(const CeedScalar grad_du[3][3], const CeedScalar F[3][3], CeedScalar dE_sym[6])

Compute incremental Green Lagrange strain dE = (grad_du^T * F + F^T * grad_du) / 2 for initial configuration.

Parameters:

• grad_du – [in] Gradient of incremental change to u

• F – [in] Deformation gradient

• dE_sym – [out] Incremental Green Lagrange strain, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelLinearStress(CeedScalar two_mu, CeedScalar lambda, CeedScalar trace_e, const CeedScalar e_sym[6], CeedScalar sigma_sym[6])

Compute linear stress sigma = lambda*trace(e)I + 2*mu*e

Parameters:

• two_mu – [in] Shear modulus multiplied by 2

• lambda – [in] Lame parameter

• trace_e – [in] Trace of linear strain

• e_sym – [in] Linear strain, in symmetric representation

• sigma_sym – [out] Linear stress, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelLinearStress_fwd(CeedScalar two_mu, CeedScalar lambda, CeedScalar trace_de, const CeedScalar de_sym[6], CeedScalar dsigma_sym[6])

Compute incremental linear stress dsigma = lambda*trace(de)I + 2*mu*de

Parameters:

• two_mu – [in] Shear modulus multiplied by 2

• lambda – [in] Lame parameter

• trace_de – [in] Trace of incremental linear strain

• de_sym – [in] Linear strain incremental, in symmetric representation

• dsigma_sym – [out] Linear stress incremental, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelMixedLinearStress(CeedScalar two_mu, CeedScalar bulk_primal, CeedScalar trace_e, CeedScalar p, const CeedScalar e_dev_sym[6], CeedScalar sigma_sym[6])

Compute mixed linear stress sigma = (bulk_primal * tr(e) - p) * I + 2 * mu * e_dev

Parameters:

• two_mu – [in] Shear modulus multiplied by 2

• bulk_primal – [in] Primal bulk modulus

• trace_e – [in] Trace of linear strain

• p – [in] Pressure

• e_dev_sym – [in] Deviatoric linear strain, in symmetric representation

• sigma_sym – [out] Mixed linear stress, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelMixedLinearStress_fwd(CeedScalar two_mu, CeedScalar bulk_primal, CeedScalar trace_de, CeedScalar dp, const CeedScalar de_dev_sym[6], CeedScalar dsigma_sym[6])

Compute incremental mixed linear stress dsigma = (bulk_primal * tr(de) - dp) * I + 2 * mu * de_dev

Parameters:

• two_mu – [in] Shear modulus multiplied by 2

• bulk_primal – [in] Primal bulk modulus

• trace_de – [in] Trace of incremental linear strain

• dp – [in] Pressure incremental

• de_dev_sym – [in] Deviatoric linear strain incremental, in symmetric representation

• dsigma_sym – [out] Mixed linear stress incremental, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelPoroElasticityLinearStress(CeedScalar two_mu_d, CeedScalar lambda_d, CeedScalar alpha, CeedScalar trace_e, CeedScalar p, const CeedScalar e_sym[6], CeedScalar sigma_sym[6])

Compute poroelastic linear stress sigma = ((lambda_d * trace_e * I) - (alpha * p * I) + (2_mu_d * e_sym)

Parameters:

• two_mu_d – [in] Shear modulus (drained) multiplied by 2

• lambda_d – [in] First Lame parameter

• alpha – [in] Biot’s coefficient

• trace_e – [in] Trace of linear strain

• p – [in] Pressure

• e_sym – [in] Poroelastic linear strain, in symmetric representation

• sigma_sym – [out] Poroelastic linear stress, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelPoroElasticityLinearStress_fwd(CeedScalar two_mu_d, CeedScalar lambda_d, CeedScalar alpha, CeedScalar trace_de, CeedScalar dp, const CeedScalar de_sym[6], CeedScalar dsigma_sym[6])

Compute incremental poroelastic linear stress dsigma = (lambda_d trace_de I) - (alpha dp I) + (2_mu_d de_sym)

Parameters:

• two_mu_d – [in] Shear modulus (drained) multiplied by 2

• lambda_d – [in] First Lame parameter

• alpha – [in] Biot’s coefficient

• trace_de – [in] Trace of incremental linear strain

• dp – [in] Pressure incremental

• de_sym – [in] linear strain incremental, in symmetric representation //Deviatoric ?

• dsigma_sym – [out] Mixed linear stress incremental, in symmetric representation

Returns:

An error code: 0 - success, otherwise - failure

int RatelOrthogonalComplement(const CeedScalar w[3], CeedScalar u[3], CeedScalar v[3])

Compute right handed orthonormal set {w, u, v} for length 3 unit vectors.

Ref https://www.geometrictools.com/GTE/Mathematics/SymmetricEigensolver3x3.h

Parameters:

• w – [in] Input vector, unit length

• u – [out] First output vector

• v – [out] Second output vector

Returns:

An error code: 0 - success, otherwise - failure

int RatelComputeEigenvector0(const CeedScalar A_sym[6], CeedScalar e_val_0, CeedScalar e_vec_0[3])

Compute unit length eigenvector for the first eigenvalue of 3x3 symmetric matrix A_sym.

Ref https://www.geometrictools.com/GTE/Mathematics/SymmetricEigensolver3x3.h

Parameters:

• A_sym – [in] Symmetric matrix

• e_val_0 – [in] First eigenvalue

• e_vec_0 – [out] First eigenvector

Returns:

An error code: 0 - success, otherwise - failure

int RatelComputeEigenvector1(const CeedScalar A_sym[6], const CeedScalar e_vec_0[3], CeedScalar e_val_1, CeedScalar e_vec_1[3])

Compute unit length eigenvector for the second eigenvalue of 3x3 symmetric matrix A_sym.

Ref https://www.geometrictools.com/GTE/Mathematics/SymmetricEigensolver3x3.h

Parameters:

• A_sym – [in] Symmetric matrix

• e_vec_0 – [in] First eigenvector

• e_val_1 – [in] Second eigenvalue

• e_vec_1 – [out] Second eigenvector

Returns:

An error code: 0 - success, otherwise - failure

int RatelMatComputeEigensystemSymmetric(const CeedScalar A_sym[6], CeedScalar e_vals[3], CeedScalar e_vecs[3][3])

Compute eigenvalues/vectors of 3x3 symmetric matrix A_sym.

Ref https://onlinelibrary.wiley.com/doi/10.1002/nme.7153?af=R

The rows of e_vecs are the eigenvectors for each eigenvalue.

Parameters:

• A_sym – [in] Symmetric matrix

• e_vals – [out] Array of eigenvalues

• e_vecs – [out] Array of eigenvectors

Returns:

An error code: 0 - success, otherwise - failure

int RatelEigenValue_fwd(const CeedScalar dA_sym[6], const CeedScalar e_vecs[3][3], CeedScalar de_vals[3])

Compute eigenvalue’s increment de_vals[i] = (dA*e_vecs[i], e_vecs[i]) for Newton linearization.

Parameters:

• dA_sym – [in] Increment of symmetric matrix A

• e_vecs – [in] Eigenvectors of A_sym

• de_vals – [out] Increment eigenvalues

Returns:

An error code: 0 - success, otherwise - failure

int RatelEigenVector_fwd(const CeedScalar dA_sym[6], const CeedScalar e_vals[3], const CeedScalar e_vecs[3][3], CeedInt i, CeedScalar de_vec[3])

Compute eigenvector’s increment de_vec[i] = sum_{j!=i} (dA*e_vecs[i], e_vecs[j]) * e_vecs[j] / (e_vals[i] - e_vals[j]) for Newton linearization.

Parameters:

• dA_sym – [in] Increment of symmetric matrix A

• e_vals – [in] Eigenvalues of A_sym

• e_vecs – [in] Eigenvectors of A_sym

• i – [in] ith de_vec that we want to return

• de_vec – [out] Increment of ith eigenvector

Returns:

An error code: 0 - success, otherwise - failure

int RatelPrincipalStretch(const CeedScalar e_vals[3], CeedScalar pr_str[3])

Compute principal stretch pr[i] = sqrt(1 + e_vals[i]) where e_vals are eigenvalue of Green-Lagrange strain tensor E (in initial configuration) or Green-Euler strain tensor e (in current configuration)

Parameters:

• e_vals – [in] Array of eigenvalues of E or e

• pr_str – [out] Array of principal stretch

Returns:

An error code: 0 - success, otherwise - failure

int RatelPrincipalStretch_fwd(const CeedScalar dA_sym[6], const CeedScalar e_vecs[3][3], const CeedScalar pr_str[3], CeedScalar dpr_str[3])

Compute principal stretch increment dpr_str[i] = (dA*e_vecs[i], e_vecs[i]) / pr[i] for Newton linearization.

Parameters:

• dA_sym – [in] Increment of Green-Lagrange dE (in initial configuration) or Green-Euler de (in current configuration)

• e_vecs – [in] Array of eigenvectors of A_sym

• pr_str – [in] Array of principal stretch

• dpr_str – [out] Array of increment principal stretch

Returns:

An error code: 0 - success, otherwise - failure

int RatelEigenVectorOuterMult(const CeedScalar e_vals[3], const CeedScalar e_vecs[3][3], CeedScalar outer_mult[3][6])

Compute outer product of eigenvector outer_mult[i][6] = (N[i] * N[i]^T) for i=0:2.

Parameters:

• e_vals – [in] Input eigenvalues

• e_vecs – [in] Input eigenvectors 3x3

• outer_mult – [out] Matrix of size 3x6

Returns:

An error code: 0 - success, otherwise - failure

int RatelEigenVectorOuterMult_fwd(const CeedScalar dA_sym[6], const CeedScalar e_vals[3], const CeedScalar e_vecs[3][3], CeedScalar outer_mult_fwd[3][6])

Compute outer product increment of eigenvector Outer_mult_fwd[i][6] = d(N[i] * N[i]^T) for Newton linearization.

Parameters:

• dA_sym – [in] Increment of symmetric matrix A

• e_vecs – [in] Eigenvectors of A

• e_vals – [in] Eigenvalues of A

• outer_mult_fwd – [out] Matrix of size 3x6

Returns:

An error code: 0 - success, otherwise - failure

int RatelVecSignum(const CeedScalar x[3], CeedScalar y[3], CeedScalar *norm)

Compute signum function of vector, defined by sgn(x) = ||x|| > 0 ? x / ||x|| : 0

Not collective across MPI processes.

Parameters:

• x – [in] Input vector

• y – [out] Output vector, y = sgn(x)

• norm – [out] Norm of input vector

Returns:

An error code: 0 - success, otherwise - failure

int RatelVecSignum_fwd(const CeedScalar x[3], const CeedScalar dx[3], CeedScalar y[3], CeedScalar dy[3], CeedScalar *norm)

Compute increment of signum of vector d(sgn(x))/dx = ||x|| > 0 ? (I - (x/||x||)⊗(x/||x||)) / ||x|| * dx : 0

Parameters:

• x – [in] Input vector

• dx – [in] Increment of input vector

• y – [out] Output vector, d = sgn(x)

• dy – [out] Increment of output vector, dy = d(sgn(x))/dx

• norm – [out] Norm of input vector

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelProcessCommandLineOptions(Ratel ratel)

Process general command line options.

Collective across MPI processes.

Parameters:

• ratel – [inout] Ratel context object

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelIncrementMeshRemapState(Ratel ratel)

Increment the state for mesh remapping to invalidate cached data.

Not collective across MPI processes.

Parameters:

• ratel – [inout] Ratel context

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelGetMeshRemapState(Ratel ratel, PetscInt *state)

Get the state for mesh remapping.

Not collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• state – [out] Current mesh remapping state

Returns:

An error code: 0 - success, otherwise - failure

void RatelRemapScale_inner(const PetscInt cdim, const PetscScalar x_prev[], const PetscReal t, const PetscInt num_constants, const PetscScalar constants[], PetscScalar x[])

Remap coordinates by scaling in a direction.

Not collective across MPI processes.

Parameters:

• cdim – [in] Coordinate dimension

• x_prev – [in] Old coordinates

• t – [in] Time

• num_constants – [in] Number of constants, should be 6

• constants – [in] Array of constants: [direction, z0, h0, hf, h_prev, h]

• x – [out] New coordinates

Returns:

none

void RatelRemapScale(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt num_constants, const PetscScalar constants[], PetscScalar f[])

Callback function to remap coordinates by scaling in a direction.

Not collective across MPI processes.

Note

This function signature is required by PETSc. SeeRatelRemapScale_inner() for the actual remapping logic

Parameters:

• dim – [in] Solution dimension

• Nf – [in] Number of fields, 1

• NfAux – [in] Number of auxiliary fields, 0

• uOff – [in] Offset of this point into local solution vector, [0, cdim]

• uOff_x – [in] Unused

• u – [in] Previous coordinates at point

• u_t – [in] Unused

• u_x – [in] Unused

• aOff – [in] Unused

• aOff_x – [in] Unused

• a – [in] Unused

• a_t – [in] Unused

• a_x – [in] Unused

• t – [in] Current time

• x – [in] Unused

• num_constants – [in] Number of constants, should be 6

• constants – [in] Array of constants

• f – [out] Remapped coordinates at points

Returns:

none

PetscErrorCode RatelSetRemapScaleParametersFromOptions(Ratel ratel, DM dm)

Set remap scale parameters from command-line options.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [inout] Solution mesh DM

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelRemapScaleCoordinates(Ratel ratel, PetscReal t, DM dm)

Remap coordinates by scaling in a direction.

Not collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• t – [in] Current time

• dm – [inout] Solution mesh DM, may have coordinates remapped by scaling

Returns:

An error code: 0 - success, otherwise - failure

static CeedScalar RatelTransition(CeedScalar a, CeedScalar b, CeedScalar x)

Transition from a value of “a” for x=0, to a value of “b” for x=1.

Collective across MPI processes.

Parameters:

• a – [in] Scalar value

• b – [in] Scalar value

• x – [in] Scalar value

Returns:

a + (b - a) * (x)

static CeedScalar RatelTransformRightBoundary_1D(CeedScalar eps, CeedScalar x)

1D transformation at the right boundary.

Collective across MPI processes.

Parameters:

• eps – [in] Scalar value in (0, 1]

• x – [in] Scalar value

Returns:

1D transformation at right boundary

static CeedScalar RatelTransformLeftBoundary_1D(CeedScalar eps, CeedScalar x)

1D transformation at the left boundary.

Collective across MPI processes.

Parameters:

• eps – [in] Scalar value in (0, 1]

• x – [in] Scalar value

Returns:

1D transformation at left boundary

PetscErrorCode RatelKershaw(DM dm, PetscScalar eps)

Apply 3D Kershaw mesh transformation.

Collective across MPI processes.

Parameters:

• dm – [inout] DMPlex holding mesh

• eps – [in] Scalar value in (0, 1]. Uniform mesh is recovered for eps=1

Returns:

An error code: 0 - success, otherwise - failure

static PetscErrorCode RatelDMPlexCeedElemRestrictionDestroy(void **ctx)

Destroy ElemRestriction in a PetscContainer.

Not collective across MPI processes.

Parameters:

• ctx – [inout] CeedElemRestriction

Returns:

An error code: 0 - success, otherwise - failure

static PetscErrorCode RatelDMPlexCeedBasisDestroy(void **ctx)

Destroy CeedBasis in a PetscContainer.

Not collective across MPI processes.

Parameters:

• ctx – [inout] CeedBasis

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMPlexCreateFromOptions(Ratel ratel, VecType vec_type, DM *dm)

Read mesh and setup DMPlex from options.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• vec_type – [in] Default vector type (can be changed at run-time using -dm_vec_type)

• dm – [out] New distributed DMPlex

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMSwarmCreateFromOptions(Ratel ratel, DM dm_mesh, DM *dm)

Setup DMSwarm from options.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm_mesh – [out] Existing cellDMPlex

• dm – [out] New distributed DMSwarm

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMSwarmInitalizePointLocations(Ratel ratel, RatelMPMOptions mpm_options, DM dm_swarm)

Initialize point locations for MPM using RatelMPMOptions context.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• mpm_options – [in] RatelMPMOptions context

• dm_swarm – [out] DMSwarm to initialize point locations for

Returns:

An error code: 0 - success, otherwise - failure

static PetscErrorCode RatelDMGetDomainISLocal(Ratel ratel, DM dm, DMLabel domain_label, PetscInt domain_value, PetscInt height, IS *is)

Get the index set for the domain of a given height.

Not collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex to get domain index set from

• domain_label – [in] DMLabel for domain stratum

• domain_value – [in] Value in domain stratum

• height – [in] Topological height, e.g., 0 for cells, 1 for faces

• is – [out] Output index set for domain and height

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMSwarmCreateReferenceCoordinates(Ratel ratel, DM dm_swarm, DMLabel domain_label, PetscInt domain_value, IS *is_points, Vec *X_points_ref)

DMSwarmCreateReferenceCoordinates - Compute the cell reference coordinates for local DMSwarm points.

Collective across MPI processes.

The index set contains ranges of indices for each local cell. This range contains the indices of every point in the cell.

total_num_cells
cell_0_start_index
cell_1_start_index
cell_2_start_index
...
cell_n_start_index
cell_n_stop_index
cell_0_point_0
cell_0_point_0
...
cell_n_point_m

Parameters:

• ratel – [in] Ratel context

• dm_swarm – [in] DMSwarm to compute reference coordinates for

• domain_label – [in] Label for domain

• domain_value – [in] Value for domain label

• is_points – [out] Index set for each local cell

• X_points_ref – [out] Reference coordinates for each local point

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMSwarmCeedElemRestrictionPointsCreate(Ratel ratel, DM dm_swarm, DMLabel domain_label, PetscInt domain_value, CeedVector *x_ref_points, CeedElemRestriction *restriction_x_points)

Create a CeedElemRestriction for the reference coordinates of the DMSwarm points.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm_swarm – [in] DMSwarm to create reference coordinates for

• domain_label – [in] DMLabel for domain in mesh DM

• domain_value – [in] Value for domain label

• x_ref_points – [out] Output CeedVector for reference coordinates

• restriction_x_points – [out] Output CeedElemRestriction for reference coordinates

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMSwarmViewFromOptions(Ratel ratel, DM dm_swarm, const char option[])

View DMSwarm fields from options.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm_swarm – [in] DMSwarm to view fields from

• option – [in] View option string

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelGetSolutionMeshDM(Ratel ratel, DM *dm)

Get the solution mesh DMPlex.

Caller is responsible for destroying the returned DMPlex.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [out] DMPlex holding solution mesh

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelGetEnergyDM(Ratel ratel, DM *dm)

Get the energy mesh DMPlex.

Caller is responsible for destroying the returned DMPlex.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [out] DMPlex holding energy mesh

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMPlexCreateFaceLabel(Ratel ratel, DM dm, RatelMaterial material, PetscInt dm_face, PetscInt *num_label_values, char **face_label_name)

Create a label with separate values for the FE orientations for all cells with this material for the DMPlex face.

Collective across MPI processes.

Parameters:

• ratel – [inout] Ratel context

• dm – [inout] DMPlex holding mesh

• material – [in] RatelMaterial to setup operators for

• dm_face – [in] Face number on DMPlex

• num_label_values – [out] Number of label values for newly created label

• face_label_name – [out] Label name for newly created label, or NULL if no elements match the material and dm_face

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode PetscFECreateLagrangeFromOptions(MPI_Comm comm, PetscInt dim, PetscInt num_comp, PetscBool is_simplex, PetscInt order, PetscInt q_order, const char prefix[], PetscFE *fem)

Setup DM with FE space of appropriate degree.

Parameters:

• comm – [in] MPI communicator

• dim – [in] Spatial dimension

• num_comp – [in] Number of components

• is_simplex – [in] Flax for simplex or tensor product element

• order – [in] Polynomial order of space

• q_order – [in] Quadrature order

• prefix – [in] The options prefix, or NULL

• fem – [out] PetscFE object

Returns:

An error code: 0 - success, otherwise - failure

static PetscErrorCode RatelDMFieldToDSField(Ratel ratel, DM dm, DMLabel domain_label, PetscInt dm_field, PetscInt *ds_field)

Convert DM field to DS field.

Not collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DM holding mesh

• domain_label – [in] Label for DM domain

• dm_field – [in] Index of DM field

• ds_field – [out] Index of DS field

Returns:

An error code: 0 - success, otherwise - failure

static inline PetscErrorCode RatelGetClosurePermutationAndFieldOffsetAtDepth(Ratel ratel, DM dm, PetscInt depth, PetscInt field, IS *permutation, PetscInt *field_offset)

Get the permutation and field offset for a given depth.

Not collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex holding mesh

• depth – [in] Depth of DMPlex topology

• field – [in] Index of DMPlex field

• permutation – [out] Permutation for DMPlex field

• field_offset – [out] Offset to DMPlex field

Returns:

An error code: 0 - success, otherwise - failure

static inline PetscErrorCode RatelGetQuadratureDataP2C(Ratel ratel, PetscFE fe, PetscQuadrature quadrature, PetscInt *q_dim, PetscInt *num_comp, PetscInt *Q, CeedScalar **q_points, CeedScalar **q_weights)

Get quadrature data from PetscQuadrature for use with libCEED.

Not collective across MPI processes.

Note

q_weights and q_points are allocated using PetscCalloc1 and must be freed by the user.

Parameters:

• ratel – [in] Ratel context

• fe – [in] PetscFE object

• quadrature – [in] PETSc basis quadrature

• q_dim – [out] Quadrature dimension, or NULL if not needed

• num_comp – [out] Number of components, or NULL if not needed

• Q – [out] Number of quadrature points, or NULL if not needed

• q_points – [out] Quadrature points in libCEED orientation, or NULL if not needed

• q_weights – [out] Quadrature weights in libCEED orientation, or NULL if not needed

Returns:

An error code: 0 - success, otherwise - failure

static inline PetscErrorCode RatelCreate1DTabulation_Tensor(Ratel ratel, PetscFE fe, PetscTabulation *tabulation, PetscReal **q_points_1d, CeedScalar **q_weights_1d)

Create 1D tabulation from PetscFE.

Collective across MPI processes.

Note

q_weights and q_points are allocated using PetscCalloc1 and must be freed by the user.

Parameters:

• ratel – [in] Ratel context

• fe – [in] PETSc basis FE object

• tabulation – [out] PETSc basis tabulation

• q_points_1d – [out] Quadrature points in libCEED orientation

• q_weights_1d – [out] Quadrature weights in libCEED orientation

Returns:

An error code: 0 - success, otherwise - failure

static inline PetscErrorCode RatelComputeFieldTabulationP2C(Ratel ratel, DM dm, PetscInt field, PetscInt face, PetscTabulation tabulation, CeedScalar **interp, CeedScalar **grad)

Compute field tabulation from PetscTabulation.

Not collective across MPI processes.

Note

interp and grad are allocated using PetscCalloc1 and must be freed by the user.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex holding mesh

• field – [in] Index of DMPlex field

• face – [in] Index of basis face, or 0

• tabulation – [in] PETSc basis tabulation

• interp – [out] Interpolation matrix in libCEED orientation

• grad – [out] Gradient matrix in libCEED orientation

Returns:

An error code: 0 - success, otherwise - failure

static inline PetscErrorCode RatelGetGlobalDMPolytopeType(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt height, DMPolytopeType *cell_type)

Get global DMPlex topology type.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex holding mesh

• domain_label – [in] DMLabel for DMPlex domain

• label_value – [in] Stratum value

• height – [in] Height of DMPlex topology

• cell_type – [out] DMPlex topology type

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMPlexCeedBasisCreate(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt height, PetscInt dm_field, CeedBasis *basis)

Create CeedBasis from DMPlex.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex holding mesh

• domain_label – [in] DMLabel for DMPlex domain

• label_value – [in] Stratum value

• height – [in] Height of DMPlex topology

• dm_field – [in] Index of DMPlex field

• basis – [out] CeedBasis for DMPlex

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMPlexCeedBasisCoordinateCreate(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt height, CeedBasis *basis)

Create CeedBasis for coordinate space of DMPlex.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex holding mesh

• domain_label – [in] DMLabel for DMPlex domain

• label_value – [in] Stratum value

• height – [in] Height of DMPlex topology

• basis – [out] CeedBasis for coordinate space of DMPlex

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMPlexCeedBasisCellToFaceCreate(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt face, PetscInt dm_field, CeedBasis *basis)

Create CeedBasis for cell to face quadrature space evaluation from DMPlex.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex holding mesh

• domain_label – [in] DMLabel for DMPlex domain

• label_value – [in] Stratum value

• face – [in] Index of face

• dm_field – [in] Index of DMPlex field

• basis – [out] CeedBasis for cell to face evaluation

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMPlexCeedBasisCellToFaceCoordinateCreate(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt face, CeedBasis *basis)

Create CeedBasis for cell to face quadrature space evaluation on coordinate space of DMPlex.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex holding mesh

• domain_label – [in] DMLabel for DMPlex domain

• label_value – [in] Stratum value

• face – [in] Index of face

• basis – [out] CeedBasis for cell to face evaluation on coordinate space of DMPlex

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMPlexCeedElemRestrictionCreate(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt height, PetscInt dm_field, CeedElemRestriction *restriction)

Create CeedElemRestriction from DMPlex.

Not collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex holding mesh

• domain_label – [in] DMLabel for DMPlex domain

• label_value – [in] Stratum value

• height – [in] Height of DMPlex topology

• dm_field – [in] Index of DMPlex field

• restriction – [out] CeedElemRestriction for DMPlex

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMPlexCeedElemRestrictionCoordinateCreate(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt height, CeedElemRestriction *restriction)

Create CeedElemRestriction from DMPlex domain for mesh coordinates.

Not collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex holding mesh

• domain_label – [in] Label for DMPlex domain

• label_value – [in] Stratum value

• height – [in] Height of DMPlex topology

• restriction – [out] CeedElemRestriction for mesh

Returns:

An error code: 0 - success, otherwise - failure

static PetscErrorCode RatelDMPlexCeedElemRestrictionStridedCreate(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt height, PetscInt dm_field, PetscInt q_data_size, PetscBool is_collocated, CeedElemRestriction *restriction)

Create CeedElemRestriction from DMPlex domain for auxiliary QFunction data.

Not collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex holding mesh

• domain_label – [in] Label for DMPlex domain

• label_value – [in] Stratum value

• height – [in] Height of DMPlex topology

• dm_field – [in] Index of DMPlex field

• q_data_size – [in] Number of components for QFunction data

• is_collocated – [in] Boolean flag indicating if the data is collocated on the nodes (PETSC_TRUE) on on quadrature points (PETSC_FALSE)

• restriction – [out] Strided CeedElemRestriction for QFunction data

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMPlexCeedElemRestrictionQDataCreate(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt height, PetscInt q_data_size, CeedElemRestriction *restriction)

Create CeedElemRestriction from DMPlex domain for auxiliary QFunction data.

Not collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex holding mesh

• domain_label – [in] Label for DMPlex domain

• label_value – [in] Stratum value

• height – [in] Height of DMPlex topology

• q_data_size – [in] Number of components for QFunction data

• restriction – [out] Strided CeedElemRestriction for QFunction data

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMPlexCeedElemRestrictionCollocatedCreate(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt height, PetscInt dm_field, PetscInt q_data_size, CeedElemRestriction *restriction)

Create CeedElemRestriction from DMPlex domain for nodally collocated auxiliary QFunction data.

Not collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex holding mesh

• domain_label – [in] Label for DMPlex domain

• label_value – [in] Stratum value

• height – [in] Height of DMPlex topology

• dm_field – [in] Index of DMPlex field

• q_data_size – [in] Number of components for QFunction data

• restriction – [out] Strided CeedElemRestriction for QFunction data

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMGetFieldISLocal(Ratel ratel, DM dm, PetscInt field, PetscInt *block_size, IS *is)

Create index set for each local field DoF in DM, including constrained and ghost nodes.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DM holding mesh

• field – [in] Field index

• block_size – [out] Number of components of requested field

• is – [out] Output IS with one index for each DoF of field

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMHasFace(Ratel ratel, DM dm, PetscInt face_id, PetscBool *has_face)

Check whether a given face id exists in the “Face Sets” DMLabel.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DM holding mesh

• face_id – [in] Stratum value of face in “Face Sets” DM:abel

• has_face – [out] PETSC_TRUE if face exists in “Face Sets” DMLabel, otherwise PETSC_FALSE

Returns:

An error code: 0 - success, otherwise - failure

static PetscErrorCode RatelForcingBodyParamsFromOptions(Ratel ratel, const char *cl_prefix, const char *cl_message, RatelForcingBodyParams *params_forcing)

PetscErrorCode RatelMaterialForcingBodyDataFromOptions(RatelMaterial material, RatelForcingBodyParams *params_forcing)

Read user forcing vectors from options database.

Collective across MPI processes.

Parameters:

• material – [in] RatelMaterial context

• params_forcing – [out] Data structure to store traction data

PetscErrorCode RatelMaterialSetupForcingSuboperator(RatelMaterial material, CeedOperator op_residual)

Add forcing term to residual u term operator.

Collective across MPI processes.

Parameters:

• material – [in] RatelMaterial context

• op_residual – [out] Composite residual u term CeedOperator to add RatelMaterial sub-operator

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelMaterialSetupForcingEnergySuboperator(RatelMaterial material, DM dm_energy, CeedOperator op_external_energy)

Add energy of forcing term to external energy operator.

Collective across MPI processes.

Parameters:

• material – [in] RatelMaterial context

• dm_energy – [in] DM for external energy computation

• op_external_energy – [out] Composite external energy CeedOperator to add energy of body force sub-operator

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelCreateSubmatrix(Mat A, IS is_row, IS is_col, MatReuse reuse_submatrix, Mat *A_sub)

Create or return a submatrix for a given MATCEED corresponding to a Jacobian matrix from Ratel.

By default, this submatrix has a MatType of MATCEED. The MatType of the sub-matrix can be set with the command line option -ceed_field_[index]_mat_type or -ceed_[field name]_mat_type.

Collective across MPI processes.

Parameters:

• A – [in] Jacobian MATCEED created by Ratel

• is_row – [in] Parallel IS for rows of the submatrix on this process

• is_col – [in] Parallel IS for all columns of the submatrix

• reuse_submatrix – [in] Either MAT_INTITIAL_MATRIX or MAT_REUSE_MATRIX

• A_sub – [out] Submatrix MATCEED

Returns:

An error code: 0 - success, otherwise - failure

static PetscErrorCode RatelSNESSetJacobianMats(Ratel ratel, SNES snes)

Set Ratel Mats for a SNES.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• snes – [inout] Non-linear solver

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSNESFormResidual(SNES snes, Vec X, Vec Y, void *ctx)

Compute the non-linear residual for a SNES.

Collective across MPI processes.

Parameters:

• snes – [in] Non-linear solver

• X – [in] Current state vector

• Y – [out] Residual vector

• ctx – [in] Ratel context

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSNESFormJacobian(SNES snes, Vec X, Mat J, Mat J_pre, void *ctx)

Assemble SNES Jacobian.

Collective across MPI processes.

Parameters:

• snes – [in] Non-linear solver

• X – [in] Current non-linear residual

• J – [out] Jacobian operator

• J_pre – [out] Jacobian operator for preconditioning

• ctx – [inout] Ratel context

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelTSFormIResidual(TS ts, PetscReal time, Vec X, Vec X_t, Vec Y, void *ctx)

Compute the non-linear residual for a TS.

Collective across MPI processes.

Parameters:

• ts – [in] Time-stepper

• time – [in] Current time

• X – [in] Current state vector

• X_t – [in] Time derivative of current state vector

• Y – [out] Function vector

• ctx – [inout] Ratel context

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelTSFormIJacobian(TS ts, PetscReal time, Vec X, Vec X_t, PetscReal v, Mat J, Mat J_pre, void *ctx)

Assemble TS Jacobian.

Collective across MPI processes.

Parameters:

• ts – [in] Time-stepper

• time – [in] Current time

• X – [in] Current non-linear residual

• X_t – [in] Time derivative of current non-linear residual

• v – [in] Shift

• J – [out] Jacobian operator

• J_pre – [out] Jacobian operator for preconditioning

• ctx – [inout] Ratel context

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelTSFormI2Residual(TS ts, PetscReal time, Vec X, Vec X_t, Vec X_tt, Vec Y, void *ctx)

Compute the non-linear residual for a TS.

Collective across MPI processes.

Parameters:

• ts – [in] Time-stepper

• time – [in] Current time

• X – [in] Current state vector

• X_t – [in] Time derivative of current state vector

• X_tt – [in] Second time derivative of current state vector

• Y – [out] Function vector

• ctx – [inout] User context struct containing CeedOperator data

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelTSFormI2Jacobian(TS ts, PetscReal time, Vec X, Vec X_t, Vec X_tt, PetscReal v, PetscReal a, Mat J, Mat J_pre, void *ctx)

Assemble TS Jacobian.

Collective across MPI processes.

Parameters:

• ts – [in] Time-stepper

• time – [in] Current time

• X – [in] Current non-linear residual

• X_t – [in] Time derivative of current non-linear residual

• X_tt – [in] Second time derivative of current non-linear residual

• v – [in] Shift for X_t

• a – [in] Shift for X_tt

• J – [out] Jacobian operator

• J_pre – [out] Jacobian operator for preconditioning

• ctx – [inout] Ratel context

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelTSPreStep(TS ts)

Migrate from DMSwarm to DMPlex before TS step, if needed.

Collective across MPI processes.

Parameters:

• ts – [in] Time-stepper

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelTSPostEvaluate(TS ts)

TS post-evaluate routine to accept state values and migrate MPM points, if needed.

Collective across MPI processes.

Parameters:

• ts – [in] Time-stepper

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSNESFormObjective(SNES snes, Vec X, PetscReal *merit, void *ctx)

Compute SNES objective.

Collective across MPI processes.

Parameters:

• snes – [in] Non-linear solver

• X – [in] Current state vector

• merit – [out] SNES objective

• ctx – [in] Ratel context

Returns:

An error code: 0 - success, otherwise - failure

static PetscErrorCode RatelCeedOperatorClone_Single(Ratel ratel, CeedOperator op, CeedOperator *op_clone)

Clone a single, non-composite CeedOperator.

Not collective across MPI processes.

Parameters:

• ratel – [in] Ratel context object

• op – [in] CeedOperator to clone

• op_clone – [out] Cloned CeedOperator

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelCeedOperatorClone(Ratel ratel, CeedOperator op, CeedOperator *op_clone)

Create a new CeedOperator that references the same fields and inputs.

This new CeedOperator will reference the same input and output vectors, so using op_clone will change any passive output CeedVectors for op.

Not collective across MPI processes.

Parameters:

• ratel – [in] Ratel context object

• op – [in] CeedOperator to clone

• op_clone – [out] Cloned CeedOperator

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelIncrementalize(Ratel ratel, PetscBool scale_dt, CeedInt num_comp, CeedInt *num_times, CeedScalar *vectors, CeedScalar *times, RatelBCInterpolationType *interp_type)

Convert flexible loading values to incremental versions by differentiating vectors.

Not collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• scale_dt – [in] Should the time derivative be scaled by the time step?

• num_comp – [in] Number of components in vectors

• num_times – [inout] Length of times and vectors arrays, maybe edited in-place

• vectors – [inout] Array of size num_times x num_comp containing the values of the flexible loading vectors, edited in-place

• times – [inout] Array of size num_times containing the times at which the flexible loading transitions occur, edited in-place

• interp_type – [inout] Type of interpolation to use for the flexible loading vectors, edited in-place

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSetupStrainEnergyEvaluator(Ratel ratel, CeedEvaluator *evaluator_strain_energy)

Setup CeedEvaluator computing strain energy.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• evaluator_strain_energy – [out] Computed strain energy evaluator

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSetupExternalEnergyEvaluator(Ratel ratel, CeedEvaluator *evaluator_external_energy)

Setup CeedEvaluator computing external energy.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• evaluator_external_energy – [out] Computed external energy evaluator

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSetupDiagnosticEvaluators(Ratel ratel, KSP *ksp_diagnostic_projection, CeedEvaluator *evaluator_projected_diagnostic, CeedEvaluator *evaluator_dual_diagnostic, CeedEvaluator *evaluator_dual_nodal_scale)

Setup CeedEvaluator computing diagnostic values.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• ksp_diagnostic_projection – [out] KSP for solving projection

• evaluator_projected_diagnostic – [out] Diagnostic quantity projected values evaluator

• evaluator_dual_diagnostic – [out] Diagnostic quantity dual space values evaluator

• evaluator_dual_nodal_scale – [out] Diagnostic quantity dual space nodal scaling evaluator

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSetupMMSErrorEvaluator(Ratel ratel, CeedEvaluator *evaluator_mms_error)

Setup context computing MMS error.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• evaluator_mms_error – [out] MMS error computation evaluator

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSetupSurfaceCentroidEvaluators(Ratel ratel, CeedEvaluator evaluators_surface_centroid[])

Setup CeedEvaluator computing surface centroids.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• evaluators_surface_centroid – [out] Surface centroid computation evaluator

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSetupFaceRestrictionNodalOperator(Ratel ratel, DM dm, PetscInt face_label_value, const char face_name[], CeedOperator *op_face_restriction)

Setup CeedOperator restricting surface forces to specific face.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• dm – [in] DMPlex containing face

• face_label_value – [in] Face label value on dm

• face_name – [in] Name of face

• op_face_restriction – [out] CeedOperator restricting values to face

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSetupSurfaceForceEvaluator(Ratel ratel, CeedEvaluator *evaluator_surface_force, CeedEvaluator evaluators_surface_force_face[])

Setup CeedEvaluator for computing surface forces.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• evaluator_surface_force – [out] Surface force computation evaluator

• evaluators_surface_force_face – [out] Surface force face restriction evaluator

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSetupSurfaceForceCellToFaceEvaluators(Ratel ratel, CeedEvaluator evaluators_surface_force_cell_to_face[])

Setup CeedEvaluator computing surface forces with cell-to-face bases.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• evaluators_surface_force_cell_to_face – [out] Surface force computation evaluators

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelComputeDiagnosticQuantities_Internal(Ratel ratel, Vec U, PetscReal time, Vec D)

Internal code for computing diagnostic quantities.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• U – [in] Computed solution vector

• time – [in] Final time value, or 1.0 for SNES solution

• D – [out] Computed diagnostic quantities vector

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelComputeSurfaceForcesCellToFace_Internal(Ratel ratel, Vec U, PetscReal time, PetscScalar *surface_forces)

Internal routine for computing diagnostic quantities on mesh faces.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• U – [in] Computed solution vector

• time – [in] Current time value, or 1.0 for SNESsolution

• surface_forces – [out] Computed face forces

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelComputeSurfaceForces_Internal(Ratel ratel, Vec U, PetscReal time, PetscScalar *surface_forces)

Internal routine for computing surface forces on mesh faces.

Collective across MPI processes.

Parameters:

• ratel – [in] Ratel context

• U – [in] Computed solution vector

• time – [in] Current time value, or 1.0 for SNES solution

• surface_forces – [out] Computed face forces

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelComputeSurfaceCentroids_Internal(Ratel ratel, Vec U, PetscReal time, PetscScalar *surface_centroids)

Internal routine for computing centroids of mesh faces.

Parameters:

• ratel – [in] Ratel context

• U – [in] Computed solution vector

• time – [in] Current time value, or 1.0 for SNES solution

• surface_centroids – [out] Computed centroids

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDebugImpl256(MPI_Comm comm, PetscBool is_debug, PetscBool all_ranks, PetscInt comm_rank, const unsigned char color, const char *format, ...)

Print Ratel debugging information in color.

Forward RatelDebug* declarations & macros.

Parameters:

• comm – [in] MPI communicator for debugging

• is_debug – [in] Boolean flag for debugging

• all_ranks – [in] Boolean flag to force printing on all ranks

• comm_rank – [in] MPI communicator rank

• color – [in] Color to print

• format – [in] Printing format

PetscCallCeed(ceed_, ...)

do {                                                              \

int ierr_q_;                                                    \

PetscStackUpdateLine;                                           \

ierr_q_ = __VA_ARGS__;                                          \

if (PetscUnlikely(ierr_q_ != CEED_ERROR_SUCCESS)) {             \

const char *error_message;                                    \

CeedGetErrorMessage(ceed_, &error_message);                   \

SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "%s", error_message); \

}                                                               \

} while (0)

Calls a libCEED function and then checks the resulting error code.

If the error code is non-zero, then a PETSc error is set with the libCEED error message.

RatelCallCeed(ratel, ...) PetscCallCeed((ratel->ceed), __VA_ARGS__)

Calls a libCEED function and then checks the resulting error code.

If the error code is non-zero, then a PETSc error is set with the libCEED error message.

RATEL_MAX_FORCE_INTERP_POINTS 16

Maximum interpolation points for user-specified forcing function.

RATEL_PI_DOUBLE 3.14159265358979323846

RATEL_EPSILON_DOUBLE 1e-16

RatelMax(a, b) ((a) > (b) ? (a) : (b))

Generic max function.

Parameters:

• a – [in] First value

• b – [in] Second value

Returns:

Maximum of both values

RatelMin(a, b) ((a) < (b) ? (a) : (b))

Generic min function.

Parameters:

• a – [in] First value

• b – [in] Second value

Returns:

Minimum of both values

RatelValueInInterval(value, min, max) ((value) >= (min) && (value) <= (max))

Check if a value is in an interval.

Parameters:

• value – [in] Value to check

• min – [in] Minimum value of the interval

• max – [in] Maximum value of the interval

FLOPS_dXdxwdetJ (FLOPS_MatMatMult + 9)

FLOPS_ScaledMass 9

FLOPS_Log1pSeries 30

FLOPS_Expm1Series 15

FLOPS_AtanSeries (15 * 2 + 2)

FLOPS_Dot3 5

FLOPS_Norm3 6

FLOPS_ScalarVecMult 3

FLOPS_VecVecSubtract 3

FLOPS_VecVecCross 9

FLOPS_VecOuterMult 12

FLOPS_VecVecOuterMult 18

FLOPS_ScalarMatMultSymmetric 6

FLOPS_MatVecMult 18

FLOPS_MatTransposeVecMult 18

FLOPS_MatTrace 2

FLOPS_MatMatAdd 9

FLOPS_MatCopy 0

FLOPS_MatMatAddSymmetric 6

FLOPS_MatMatMatAddSymmetric 12

FLOPS_MatMatMatAdd (9 * 2)

FLOPS_MatDeviatoricSymmetric 6

FLOPS_MatMatMult 54

FLOPS_MatMatTransposeMult 54

FLOPS_MatMatMultPlusMatMatMult 108

FLOPS_MatMatContractSymmetric 15

FLOPS_MatDetA 14

FLOPS_MatDetAM1 29

FLOPS_MatMatMultSymmetric 36

FLOPS_MatNorm 19

FLOPS_MatMatMatMultSymmetric (FLOPS_MatMatMult + 18)

FLOPS_MatInverse 36

FLOPS_MatInverseSymmetric 24

FLOPS_CInverse (14 + FLOPS_MatInverseSymmetric)

FLOPS_CInverse_fwd (FLOPS_MatMatMatMultSymmetric + 6)

FLOPS_LinearStrain 6

FLOPS_LinearStrain_fwd 6

FLOPS_GreenEulerStrain (FLOPS_LinearStrain + 54)

FLOPS_GreenEulerStrain_fwd 90

FLOPS_GreenLagrangeStrain (FLOPS_LinearStrain + 54)

FLOPS_GreenLagrangeStrain_fwd 90

FLOPS_LinearStress 12

FLOPS_LinearStress_fwd 12

FLOPS_MixedLinearStress 15

FLOPS_MixedLinearStress_fwd 15

FLOPS_PoroElasticityLinearStress 18

FLOPS_PoroElasticityLinearStress_fwd 18

FLOPS_OrthogonalComplement (7 + FLOPS_VecVecCross)

FLOPS_ComputeEigenvector0 (4 + 3 * FLOPS_VecVecCross + 3 * FLOPS_Dot3 + FLOPS_ScalarVecMult)

FLOPS_ComputeEigenvector1   (FLOPS_OrthogonalComplement + 2 * FLOPS_MatVecMult + 3 * FLOPS_Dot3 + 2 * FLOPS_ScalarVecMult + FLOPS_VecVecSubtract + 8)

FLOPS_MatComputeEigenValueSymmetric   (FLOPS_MatTrace + 20 + 2 * FLOPS_MatDeviatoricSymmetric + FLOPS_MatMatMultSymmetric + 2 * FLOPS_MatMatAddSymmetric + 12 + 2 * FLOPS_MatNorm + 20 + \    FLOPS_AtanSeries)

FLOPS_MatComputeEigenVectorSymmetric (FLOPS_ComputeEigenvector0 + FLOPS_ComputeEigenvector1 + FLOPS_VecVecCross)

FLOPS_MatComputeEigensystemSymmetric (FLOPS_MatComputeEigenValueSymmetric + FLOPS_MatComputeEigenVectorSymmetric)

FLOPS_EigenValue_fwd (3 * FLOPS_MatVecMult + 3 * FLOPS_Dot3)

FLOPS_EigenVector_fwd (2 * FLOPS_MatVecMult + 2 * FLOPS_Dot3 + 20)

FLOPS_PrincipalStretch (9)

FLOPS_PrincipalStretch_fwd (FLOPS_EigenValue_fwd + 3)

FLOPS_EigenVectorOuterMult (3 * FLOPS_VecOuterMult)

FLOPS_EigenVectorOuterMult_fwd (3 * FLOPS_EigenVector_fwd + 6 * FLOPS_VecVecOuterMult + 3 * FLOPS_MatMatAdd)

FLOPS_VecSignum (3 + 1 + FLOPS_Norm3)

FLOPS_VecSignum_fwd (3 * 3 + 1 + FLOPS_VecSignum + FLOPS_Dot3)

RATEL_DIAGNOSTIC_PROJECTION_SCALAR_FIELD_NAME "scalar_field"

RATEL_DIAGNOSTIC_PROJECTION_SCALAR_LABEL_NAME "scalar_field_label"

RATEL_DIAGNOSTIC_PROJECTION_SCALAR_LABEL_VALUE 42

RATEL_DIAGNOSTIC_PROJECTION_SCALAR_LABEL_VALUE_STRING "42"

RATEL_DIAGNOSTIC_PROJECTION_SCALAR_OPTION_NAME "diagnostic_projection_scalar_mass_matrix"

RATEL_DIAGNOSTIC_PROJECTION_FULL_DIAGONAL_NAME "diagnostic_projection_full_diagonal"

struct RatelForcingBodyParams

#include <forcing.h>

User-specified forcing context.

Public Members

CeedScalar rho

Density.

CeedScalar acceleration[3 * RATEL_MAX_FORCE_INTERP_POINTS]

Acceleration vector.

CeedScalar times[RATEL_MAX_FORCE_INTERP_POINTS]

Transition times.

CeedInt num_times

Number of transition times.

RatelBCInterpolationType interpolation_type

Type of interpolation between points.

CeedScalar time

Current solver time.

struct RatelScaledMassParams

#include <scaled-mass.h>

Scaled mass context.

Public Members

CeedScalar rho

Scaling parameter.

CeedInt num_fields

Number of fields.

CeedInt field_sizes[RATEL_MAX_FIELDS]

Number of components in each field.

struct RatelViewer

#include <ratel-impl.h>

Ratel viewer information for monitor routines.

struct RatelViewers

#include <ratel-impl.h>

Ratel context for wrapping a collection of RatelViewer objects.

struct RatelCheckpointCtx

#include <ratel-impl.h>

Ratel context for checkpoints.

struct RatelCheckpointData

#include <ratel-impl.h>

Ratel data stored in checkpoint file.

struct RatelMPMOptions

#include <ratel-impl.h>

Ratel context for MPM options.

struct RatelPMGContext

#include <ratel-pmg-impl.h>

Ratel pMG context, does rely upon Ratel DM setup.