Internal Functions#
These are additional functions that do not clearly fit into another category.
-
static inline PetscInt RatelIntMin(PetscInt a, PetscInt b)#
Return minimum of two integers.
- Parameters:
a – [in] The first integer to compare
b – [in] The second integer to compare
- Returns:
The minimum of the two integers
-
static inline PetscInt RatelIntMax(PetscInt a, PetscInt b)#
Return maximum of two integers.
- Parameters:
a – [in] The first integer to compare
b – [in] The second integer to compare
- Returns:
The maximum of the two integers
-
static inline CeedInt RatelGetQuadratureSize(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
-
static inline CeedMemType RatelMemTypeP2C(PetscMemType mem_type)#
Translate PetscMemType to CeedMemType.
- Parameters:
mem_type – [in] PetscMemType
- Returns:
Equivalent CeedMemType
-
PetscErrorCode RatelProcessCommandLineOptions(Ratel ratel)#
Process general command line options.
- Parameters:
ratel – [inout] Ratel context object
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelDMCreateFromOptions(Ratel ratel, VecType vec_type, DM *dm)#
Read mesh and set from options.
- 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 DM
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelDMFaceLabelCreate(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.
- Parameters:
ratel – [inout] Ratel context
dm – [inout] DM holding mesh
material – [in] RatelMaterial to setup operators for
dm_face – [in] Face number on PETSc 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_index
anddm_face
- 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.
- Parameters:
ratel – [in] Ratel context
dm – [in] DM holding mesh
domain_label – [in] Label for DMPlex domain
dm_field – [in] Index of DMPlex field
ds_field – [out] Index of DS field
- Returns:
An error code: 0 - success, otherwise - failure
-
static PetscErrorCode RatelBasisCreateFromTabulation(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt height, PetscInt face, PetscInt dm_field, PetscFE fe, PetscTabulation basis_tabulation, PetscQuadrature quadrature, CeedBasis *basis)#
Create libCEED Basis from PetscTabulation.
- Parameters:
ratel – [in] Ratel context
dm – [in] DM holding mesh
domain_label – [in] Label for DMPlex domain
label_value – [in] Stratum value
height – [in] Height of DMPlex topology
face – [in] Index of basis face, or 0
dm_field – [in] Index of DMPlex field
fe – [in] PETSc basis FE object
basis_tabulation – [in] PETSc basis tabulation
quadrature – [in] PETSc basis quadrature
basis – [out] libCEED basis
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelBasisCreateFromPlex(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt height, PetscInt dm_field, CeedBasis *basis)#
Create libCEED Basis from DMPlex.
- Parameters:
ratel – [in] Ratel context
dm – [in] DM 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
basis – [out] libCEED basis
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelOrientedBasisCreateCellToFaceFromPlex(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt face, PetscInt dm_field, CeedBasis *basis)#
Create libCEED Basis for cell to face from DMPlex.
- Parameters:
ratel – [in] Ratel context
dm – [in] DM holding mesh
domain_label – [in] Label for DMPlex domain
label_value – [in] Stratum value
face – [in] Index of face
dm_field – [in] Index of DMPlex field
basis – [out] libCEED basis
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelRestrictionCreateFromPlex(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt height, PetscInt dm_field, CeedElemRestriction *restriction)#
Create local libCEED ElemRestriction from DMPlex.
- Parameters:
ratel – [in] Ratel context
dm – [in] DM 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
restriction – [out] libCEED element restriction
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelGetRestrictionForDomain(Ratel ratel, DM dm, DMLabel domain_label, PetscInt label_value, PetscInt height, PetscInt dm_field, CeedInt Q, CeedInt q_data_size, CeedElemRestriction *restriction_u, CeedElemRestriction *restriction_x, CeedElemRestriction *restriction_qd)#
Create local libCEED restriction from DMPlex domain.
- Parameters:
ratel – [in] Ratel context
dm – [in] DM 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 – [in] Number of quadrature points
q_data_size – [in] Size of stored quadrature point data
restriction_u – [out] libCEED element restriction for solution
restriction_x – [out] libCEED element restriction for mesh
restriction_qd – [out] libCEED element restriction for quadrature point data
- 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 Y_loc, CeedVector x_loc, CeedVector y_loc, CeedOperator op, RatelOperatorApplyContext *ctx)#
Setup context data for operator application.
- 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
Y_loc – [in] Output PETSc local vector, or NULL
x_loc – [in] Input libCEED local vector, or NULL
y_loc – [in] Output libCEED local vector, or NULL
op – [in] libCEED operator
ctx – [out] Context data for operator evaluation
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelOperatorApplyContextDestroy(RatelOperatorApplyContext ctx)#
Destroy context data for operator application.
- Parameters:
ctx – [inout] Context data for operator evaluation
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelProlongRestrictContextCreate(Ratel ratel, const char *name, DM dm_c, DM dm_f, Vec X_loc_c, Vec X_loc_f, CeedVector x_loc_c, CeedVector x_loc_f, CeedVector y_loc_c, CeedVector y_loc_f, CeedOperator op_prolong, CeedOperator op_restrict, RatelProlongRestrictContext *ctx)#
Setup context data for prolongation/restriction.
- Parameters:
ratel – [in] Ratel context
name – [in] Context name
dm_c – [in] Coarse grid DM
dm_f – [in] Fine grid DM
X_loc_c – [in] Coarse grid input PETSc local vector
X_loc_f – [in] Fine grid input PETSc local vector
x_loc_c – [in] Coarse grid input libCEED local vector, or NULL
x_loc_f – [in] Fine grid input libCEED local vector, or NULL
y_loc_c – [in] Coarse grid output libCEED local vector, or NULL
y_loc_f – [in] Fine grid output libCEED local vector, or NULL
op_prolong – [in] libCEED prolongation operator
op_restrict – [in] libCEED restriction operator
ctx – [out] Context data for prolongation/restriction evaluation
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelProlongRestrictContextDestroy(RatelProlongRestrictContext ctx)#
Destroy context data for prolongation/restriction.
- Parameters:
ctx – [inout] Context data for prolongation/restriction evaluation
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelGetDiagonal(Mat A, Vec D)#
Compute the diagonal of an operator via libCEED.
- Parameters:
A – [in] Operator MatShell
D – [out] Vector holding operator diagonal
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelMatSetPreallocationCOO(Ratel ratel, CeedOperator op, CeedVector *coo_values, Mat A)#
Set COO preallocation for Mat associated with a CeedOperator.
- Parameters:
ratel – [in] Ratel context
op – [in] CeedOperator to assemble sparsity pattern
coo_values – [out] CeedVector to hold COO values
A – [inout] Sparse matrix to hold assembled CeedOperator
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelMatSetValuesCOO(Ratel ratel, CeedOperator op, CeedVector coo_values, CeedMemType mem_type, Mat A)#
Set COO values for Mat associated with a CeedOperator.
- Parameters:
ratel – [in] Ratel context
op – [in] CeedOperator to assemble
coo_values – [in] CeedVector to hold COO values
mem_type – [in] MemType to pass COO values
A – [inout] Sparse matrix to hold assembled CeedOperator
- Returns:
An error code: 0 - success, otherwise - failure
-
static inline PetscErrorCode RatelApplyAddCeedOperator_Core(Vec X_loc, Vec X, Vec Y_loc, Vec Y, RatelOperatorApplyContext ctx)#
Apply the local action of a libCEED operator and store result in PETSc vector i.e.
compute A X = Y
- Parameters:
X_loc – [in] Input PETSc local vector
X – [in] Input PETSc global vector
Y_loc – [out] Output PETSc local vector
Y – [out] Output PETSc global vector, or NULL to only write into local vector
ctx – [in] User context struct containing libCEED operator data
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelApplyCeedOperatorGlobalToLocal(Vec X, Vec Y_loc, RatelOperatorApplyContext ctx)#
Apply the local action of a libCEED operator and store result in PETSc vector i.e.
compute A X = Y
- Parameters:
X – [in] Input PETSc global vector
Y_loc – [out] Output PETSc local vector
ctx – [in] User context struct containing libCEED operator data
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelApplyCeedOperatorGlobalToGlobal(Vec X, Vec Y, RatelOperatorApplyContext ctx)#
Apply the local action of a libCEED operator and store result in PETSc vector i.e.
compute A X = Y
- Parameters:
X – [in] Input PETSc global vector
Y – [out] Output PETSc global vector
ctx – [in] User context struct containing libCEED operator data
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelApplyAddCeedOperatorGlobalToGlobal(Vec X, Vec Y, RatelOperatorApplyContext ctx)#
Apply the local action of a libCEED operator and store result in PETSc vector i.e.
compute A X = Y
- Parameters:
X – [in] Input PETSc global vector
Y – [out] Output PETSc global vector
ctx – [in] User context struct containing libCEED operator data
- Returns:
An error code: 0 - success, otherwise - failure
-
static PetscErrorCode RatelGetCurrentTime(Ratel ratel, PetscReal *current_time)#
Get the current time used in libCEED operators.
- 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 libCEED operators.
- Parameters:
ratel – [inout] Ratel context to update boundary conditions
time – [in] Current time
- 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 libCEED.
- 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.
- 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 operator libCEED.
- Parameters:
A – [in] Operator MatShell
X – [in] Input PETSc vector
Y – [out] Output PETSc vector
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelApplyJacobian(Mat A, Vec X, Vec Y)#
Compute the action of the Jacobian for a SNES operator via libCEED.
- Parameters:
A – [in] Jacobian operator MatShell
X – [in] Input PETSc vector
Y – [out] Output PETSc vector
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelApplyProlongation(Mat A, Vec C, Vec F)#
Compute action of a p-multigrid prolongation operator via libCEED.
- Parameters:
A – [in] Prolongation operator MatShell
C – [in] Input PETSc vector on coarse grid
F – [out] Output PETSc vector on fine grid
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelApplyRestriction(Mat A, Vec F, Vec C)#
Compute action of a p-multigrid restriction operator via libCEED.
- Parameters:
A – [in] Restriction operator MatShell
F – [in] Input PETSc vector on fine grid
C – [out] Output PETSc vector on coarse grid
- 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 operator via libCEED.
- Parameters:
snes – [in] Non-linear solver
X – [in] Current state vector
Y – [out] Residual vector
ctx_residual_u – [in] User context struct containing libCEED operator data
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelSNESFormJacobian(SNES snes, Vec X, Mat J, Mat J_pre, void *ctx_form_jacobian)#
Assemble SNES Jacobian for each level of multigrid hierarchy.
- Parameters:
snes – [in] Non-linear solver
X – [in] Current non-linear residual
J – [out] Fine grid Jacobian operator
J_pre – [out] Fine grid Jacobian operator for preconditioning
ctx_form_jacobian – [inout] Jacobian context, holding Jacobian operators for each level of multigrid hierarchy
- 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 operator via libCEED.
- 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 libCEED operator 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_form_jacobian)#
Assemble TS Jacobian for each level of multigrid hierarchy.
- 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] Fine grid Jacobian operator
J_pre – [out] Fine grid Jacobian operator for preconditioning
ctx_form_jacobian – [inout] Jacobian context, holding Jacobian operators for each level of multigrid hierarchy
- 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 operator via libCEED.
- 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 libCEED operator 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_form_jacobian)#
Assemble TS Jacobian for each level of multigrid hierarchy.
- 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] Fine grid Jacobian operator
J_pre – [out] Fine grid Jacobian operator for preconditioning
ctx_form_jacobian – [inout] Jacobian context, holding Jacobian operators for each level of multigrid hierarchy
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelTSPostStage(TS ts, PetscReal stage_time, PetscInt stage_index, Vec *X)#
TS post-stage routine to accept state values.
- 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
-
PetscErrorCode RatelIntArrayC2P(PetscInt num_entries, CeedInt **array_ceed, PetscInt **array_petsc)#
Translate array of CeedInt to PetscInt.
If the types differ,
array_ceed
is freed withfree()
andarray_petsc
is allocated withmalloc()
. Caller is responsible for freeingarray_petsc
withfree()
.- Parameters:
num_entries – [in] Number of array entries
array_ceed – [inout] Array of CeedInts
array_petsc – [out] Array of PetscInts
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelIntArrayP2C(PetscInt num_entries, PetscInt **array_petsc, CeedInt **array_ceed)#
Translate array of PetscInt to CeedInt.
If the types differ,
array_petsc
is freed withPetscFree()
andarray_ceed
is allocated withPetscMalloc1()
. Caller is responsible for freeingarray_ceed
withPetscFree()
.- Parameters:
num_entries – [in] Number of array entries
array_petsc – [inout] Array of PetscInts
array_ceed – [out] Array of CeedInts
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelVecP2C(Ratel ratel, Vec X_petsc, PetscMemType *mem_type, CeedVector x_ceed)#
Transfer array from PETSc Vec to CeedVector.
- Parameters:
ratel – [in] Ratel context object
X_petsc – [in] PETSc Vec
mem_type – [out] PETSc MemType
x_ceed – [out] CeedVector
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelVecC2P(Ratel ratel, CeedVector x_ceed, PetscMemType mem_type, Vec X_petsc)#
Transfer array from CeedVector to PETSc Vec.
- Parameters:
ratel – [in] Ratel context object
x_ceed – [in] CeedVector
mem_type – [in] PETSc MemType
X_petsc – [out] PETSc Vec
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelVecReadP2C(Ratel ratel, Vec X_petsc, PetscMemType *mem_type, CeedVector x_ceed)#
Transfer read only array from PETSc Vec to CeedVector.
- Parameters:
ratel – [in] Ratel context object
X_petsc – [in] PETSc Vec
mem_type – [out] PETSc MemType
x_ceed – [out] CeedVector
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelVecReadC2P(Ratel ratel, CeedVector x_ceed, PetscMemType mem_type, Vec X_petsc)#
Transfer read only array from CeedVector to PETSc Vec.
- Parameters:
ratel – [in] Ratel context object
x_ceed – [in] CeedVector
mem_type – [in] PETSc MemType
X_petsc – [out] PETSc Vec
- 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.
- 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 RatelComputeSurfaceForces_Internal(Ratel ratel, Vec U, PetscScalar time, PetscScalar *surface_forces)#
Internal routine for computing diagnostic quantities on mesh faces.
- 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 RatelComputeReactionForces_Internal(Ratel ratel, Vec U, PetscScalar time, PetscScalar *reaction_forces)#
Internal routine for computing reaction forces on mesh faces.
- Parameters:
ratel – [in] Ratel context
U – [in] Computed solution vector
time – [in] Current time value, or 1.0 for SNES solution
reaction_forces – [out] Computed face forces
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelTSMonitor(TS ts, PetscInt steps, PetscReal time, Vec U, void *ctx)#
RatelTSMonitor to allow multiple monitor outputs.
- Parameters:
ts – [in] TS object to monitor
steps – [in] Iteration number
time – [in] Current time
U – [in] Current state vector
ctx – [in] Monitor context holding Ratel context object
- Returns:
An error code: 0 - success, otherwise - failure
-
void 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.
- 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