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_Residual#
PetscLogEvent RATEL_Jacobian#
PetscLogEvent RATEL_CeedOperatorApply#
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 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 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

    • 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:
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:
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:
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)

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 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 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 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 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 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 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

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

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

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 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

  • 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 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

  • 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 RatelOperatorApplyContextCreate(Ratel ratel, const char *name, DM dm_x, DM dm_y, Vec X_loc, Vec X_dot_loc, Vec Y_loc, CeedVector x_loc, CeedVector x_dot_loc, CeedVector y_loc, CeedOperator op, PetscLogEvent log_event, RatelOperatorApplyContext *ctx)#

Setup context data for operator application.

Collective across MPI processes.

Parameters:
  • ratel[in] Ratel context

  • name[in] Context name

  • dm_x[in] Input DM

  • dm_y[in] Output DM

  • X_loc[in] Input PETSc local vector, or NULL

  • X_dot_loc[in] Input PETSc local time derivative vector, or NULL

  • Y_loc[in] Output PETSc local vector, or NULL

  • x_loc[in] Input libCEED local vector, or NULL

  • x_dot_loc[in] Input libCEED local time derivative vector, or NULL

  • y_loc[in] Output libCEED local vector, or NULL

  • op[in] CeedOperator for evaluation

  • log_event[in] PetscLogEvent for performance metrics

  • ctx[out] Context data for operator evaluation

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelOperatorApplyContextReference(RatelOperatorApplyContext ctx)#

Increment reference counter for operator context.

Not collective across MPI processes.

Parameters:
  • ctx[inout] Context data for operator evaluation

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelOperatorApplyContextReferenceCopy(RatelOperatorApplyContext ctx, RatelOperatorApplyContext *ctx_copy)#

Copy reference for operator context.

Note: If ctx_copy is non-null, it is assumed to be a valid pointer to a RatelOperatorApplyContext.

Not collective across MPI processes.

Parameters:
  • ctx[in] Context data for operator evaluation

  • ctx_copy[out] Copy of pointer to context data for operator evaluation

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelOperatorApplyContextDestroy(RatelOperatorApplyContext ctx)#

Destroy context data for operator application.

Collective across MPI processes.

Parameters:
  • ctx[inout] Context data for operator evaluation

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 inline PetscErrorCode RatelApplyAddCeedOperator_Core(RatelOperatorApplyContext ctx, Vec X, Vec Y, Vec X_loc, Vec X_dot_loc, Vec Y_loc)#

Apply the local action of a CeedOperator and store result in PETSc vector.

Collective across MPI processes.

Parameters:
  • ctx[in] User context struct containing CeedOperator data

  • X[in] Input PETSc global vector, for logging only

  • Y[in] Output PETSc global vector, for logging only

  • X_loc[in] Input PETSc local vector

  • X_dot_loc[in] Input PETSc local time derivative vector, or NULL if and only if X_dot is NULL

  • Y_loc[out] Output PETSc local vector

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelApplyCeedOperatorLocalToLocal(RatelOperatorApplyContext ctx, Vec X_loc, Vec Y_loc)#

Apply the local action of a CeedOperator and store result in PETSc vector.

Collective across MPI processes.

Parameters:
  • ctx[in] User context struct containing CeedOperator data

  • X_loc[in] Input PETSc local vector

  • Y_loc[out] Output PETSc local vector

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelApplyCeedOperatorGlobalToLocal(RatelOperatorApplyContext ctx, Vec X, Vec Y_loc)#

Apply the local action of a CeedOperator and store result in PETSc vector.

Collective across MPI processes.

Parameters:
  • ctx[in] User context struct containing CeedOperator data

  • X[in] Input PETSc global vector

  • Y_loc[out] Output PETSc local vector

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelApplyCeedOperatorGlobalToGlobal(RatelOperatorApplyContext ctx, Vec X, Vec Y)#

Apply the local action of a CeedOperator and store result in PETSc vector.

Collective across MPI processes.

Parameters:
  • ctx[in] User context struct containing CeedOperator data

  • X[in] Input PETSc global vector

  • Y[out] Output PETSc global vector

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelApplyCeedOperatorVelocityGlobalToGlobal(RatelOperatorApplyContext ctx, Vec X, Vec X_dot, Vec Y)#

Apply the local action of a velocity-dependent CeedOperator and store result in PETSc vector.

Collective across MPI processes.

Parameters:
  • ctx[in] User context struct containing CeedOperator data

  • X[in] Input PETSc global vector

  • X_dot[in] Input PETSc global time derivative vector

  • Y[out] Output PETSc global vector

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelApplyAddCeedOperatorGlobalToGlobal(RatelOperatorApplyContext ctx, Vec X, Vec Y)#

Apply the local action of a CeedOperator and store result in PETSc vector.

Collective across MPI processes.

Parameters:
  • ctx[in] User context struct containing CeedOperator data

  • X[in] Input PETSc global vector

  • Y[out] Output PETSc global vector

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelGetCurrentTime(Ratel ratel, PetscReal *current_time)#

Get the current time used in all CeedOperator.

Not collective across MPI processes.

Parameters:
  • ratel[in] Ratel context to get time

  • current_time[out] Current time

Returns:

An error code: 0 - success, otherwise - failure

static PetscErrorCode RatelUpdateTimeContexts(Ratel ratel, PetscReal time)#

Update the current time used in all CeedOperator.

Not collective across MPI processes.

Parameters:
  • ratel[inout] Ratel context to update boundary conditions

  • time[in] Current time

Returns:

An error code: 0 - success, otherwise - failure

static PetscErrorCode RatelGetCurrentTimeStep(Ratel ratel, PetscReal *dt)#
PetscErrorCode RatelUpdateTimeStepContexts(Ratel ratel, PetscReal dt)#

Update the current time step size used in all CeedOperator.

Not collective across MPI processes.

Parameters:
  • ratel[inout] Ratel context to update boundary conditions

  • dt[in] Current time step

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelDMPlexInsertBoundaryValues(DM dm, PetscBool insert_essential, Vec U_loc, PetscReal time, Vec face_geometry_FVM, Vec cell_geometry_FVM, Vec grad_FVM)#

Compute the Dirichlet boundary values via CeedOperator.

Collective across MPI processes.

Parameters:
  • dm[in] DM to insert boundary values on

  • insert_essential[in] Boolean flag for Dirichlet or Neumann boundary conditions

  • U_loc[out] Local vector to insert the boundary values into

  • time[in] Current time

  • face_geometry_FVM[in] Face geometry data for finite volume discretizations

  • cell_geometry_FVM[in] Cell geometry data for finite volume discretizations

  • grad_FVM[in] Gradient reconstruction data for finite volume discretizations

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelUpdateTimeAndBoundaryValues(Ratel ratel, PetscReal time)#

Update time contexts and boundary values.

Collective across MPI processes.

Parameters:
  • ratel[inout] Ratel context to update boundary conditions

  • time[in] Current time

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelApplyOperator(Mat A, Vec X, Vec Y)#

Compute the action of the matrix-free CeedOperator.

Collective across MPI processes.

Parameters:
  • A[in] Operator MatShell

  • X[in] Input PETSc vector

  • Y[out] Output PETSc vector

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_residual_u)#

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_residual_u[in] User context struct containing CeedOperator data

Returns:

An error code: 0 - success, otherwise - failure

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

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_jacobian[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_residual_ut)#

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_residual_ut[inout] User context struct containing CeedOperator data

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_jacobian)#

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_jacobian[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_residual_utt)#

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_residual_utt[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_jacobian)#

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_jacobian[inout] Ratel context

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelTSPreStep(TS ts)#
PetscErrorCode RatelTSPostEvaluate(TS ts)#
PetscErrorCode RatelTSPostStage(TS ts, PetscReal stage_time, PetscInt stage_index, Vec *X)#

TS post-stage routine to accept state values.

Collective across MPI processes.

Parameters:
  • ts[in] Time-stepper

  • stage_time[in] Current time

  • stage_index[in] Index of current stage solution vector

  • X[in] Array of accepted stage solution vectors

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 RatelSetupStrainEnergyOperator(Ratel ratel, RatelOperatorApplyContext *ctx_strain_energy)#

Setup context computing strain energy.

Collective across MPI processes.

Parameters:
  • ratel[in] Ratel context

  • ctx_strain_energy[out] Computed strain energy context

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSetupDiagnosticOperators(Ratel ratel, KSP *ksp_diagnostic_projection, RatelOperatorApplyContext *ctx_projected_diagnostic, RatelOperatorApplyContext *ctx_dual_diagnostic)#

Setup contexts computing diagnostic values.

Collective across MPI processes.

Parameters:
  • ratel[in] Ratel context

  • ksp_diagnostic_projection[out] KSP for solving projection

  • ctx_projected_diagnostic[out] Diagnostic quantity projected values context

  • ctx_dual_diagnostic[out] Diagnostic quantity dual space values context

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSetupSurfaceForceCellToFaceOperators(Ratel ratel, RatelOperatorApplyContext *ctx_surface_force_cell_to_face, CeedOperator ops_surface_force[])#

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

Collective across MPI processes.

Parameters:
  • ratel[in] Ratel context

  • ctx_surface_force_cell_to_face[out] Surface force computation context

  • ops_surface_force[out] Surface force CeedOperator

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)#
PetscErrorCode RatelSetupSurfaceForceOperator(Ratel ratel, RatelOperatorApplyContext *ctx_surface_force)#

Setup context for computing surface forces.

Collective across MPI processes.

Parameters:
  • ratel[in] Ratel context

  • ctx_surface_force[out] Surface force computation context

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSetupSurfaceCentroidOperators(Ratel ratel, RatelOperatorApplyContext *ctx_surface_centroid, CeedOperator ops_surface_centroid[])#

Setup contexts computing surface centroids.

Collective across MPI processes.

Parameters:
  • ratel[in] Ratel context

  • ctx_surface_centroid[out] Surface centroid computation context

  • ops_surface_centroid[out] Surface centroid CeedOperator

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelSetupMMSErrorOperator(Ratel ratel, RatelOperatorApplyContext *ctx_mms_error)#

Setup context computing MMS error.

Collective across MPI processes.

Parameters:
  • ratel[in] Ratel context

  • ctx_mms_error[out] MMS error computation context

Returns:

An error code: 0 - success, otherwise - failure

PetscErrorCode RatelComputeDiagnosticQuantities_Internal(Ratel ratel, Vec U, PetscScalar 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, PetscScalar 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, PetscScalar 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, PetscScalar 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_ = __VA_ARGS__;                                      \

if (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 18#
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_MatTrace 2#
FLOPS_MatMatAdd 9#
FLOPS_MatMatAddSymmetric 6#
FLOPS_MatMatMatAddSymmetric 12#
FLOPS_MatDeviatoricSymmetric 6#
FLOPS_MatMatMult 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_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 + 21)#
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)#
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 RatelOperatorApplyContext#
#include <ratel-impl.h>

Ratel operator context for applying composite CeedOperator on a DM.

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.