Core Ratel Functions#
These functions are publicly exposed for users.
-
enum RatelSolverType#
Specify numerical method.
Values:
-
enumerator RATEL_SOLVER_STATIC#
Static solver, to be solved with a PETSc SNES.
-
enumerator RATEL_SOLVER_QUASISTATIC#
Quasistatic solver, to be solved with a PETSc TS.
-
enumerator RATEL_SOLVER_DYNAMIC#
Fully dynamic solver, to be solved with a PETSc TS.
-
enumerator RATEL_SOLVER_STATIC#
-
enum RatelForcingType#
Specify forcing term.
Values:
-
enumerator RATEL_FORCING_NONE#
No forcing term.
-
enumerator RATEL_FORCING_CONSTANT#
Constant forcing given by forcing vector.
-
enumerator RATEL_FORCING_MMS#
Forcing for linear elasticity manufactured solution.
-
enumerator RATEL_FORCING_NONE#
-
enum RatelInitialConditionType#
Specify initial condition.
Values:
-
enumerator RATEL_INITIAL_CONDITION_ZERO#
Initial condition zero vector.
-
enumerator RATEL_INITIAL_CONDITION_CONTINUE#
Initial condition given by binary file.
-
enumerator RATEL_INITIAL_CONDITION_ZERO#
-
enum RatelPMultigridCoarseningType#
Specify p-multigrid coarsening strategy.
Values:
-
enumerator RATEL_P_MULTIGRID_COARSENING_LOGARITHMIC#
P-multigrid, coarsen logarithmically, decreasing basis order to next lowest power of 2.
-
enumerator RATEL_P_MULTIGRID_COARSENING_UNIFORM#
P-multigrid, coarsen uniformly, decreasing basis order by 1.
-
enumerator RATEL_P_MULTIGRID_COARSENING_USER#
P-multigrid, user defined coarsening strategy.
-
enumerator RATEL_P_MULTIGRID_COARSENING_LOGARITHMIC#
-
typedef struct Ratel_private *Ratel#
Library context created by RatelInit()
-
PetscErrorCode RatelHasMMS(Ratel ratel, PetscBool *has_mms)#
Determine if
Ratel
context has MMS solution.Diagnostic quantities.
Not collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contexthas_mms – [out] Boolean flag indicating if
Ratel
context has MMS solution
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelComputeMMSL2Error(Ratel ratel, Vec U, PetscInt *num_fields, PetscScalar **l2_error)#
Compute the L2 error from the manufactured solution.
Note: User is responsible for freeing array of L2 error values.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextU – [in] Computed solution vector
num_fields – [out] Number of L2 error values
l2_error – [out] Computed L2 error for solution compared to MMS
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelGetExpectedStrainEnergy(Ratel ratel, PetscBool *has_expected, PetscScalar *expected_strain_energy)#
Retrieve expected strain energy.
Not collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contexthas_expected – [out]
PETSC_TRUE
if expected strain energy was setexpected_strain_energy – [out] Expected strain energy
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelGetExpectedMaxDisplacement(Ratel ratel, PetscInt *num_components, PetscBool has_expected[], PetscScalar expected_max_displacement[])#
Retrieve expected max displacements.
Not collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextnum_components – [inout] Size of input arrays / number of output components
has_expected – [out] Boolean array,
PETSC_TRUE
at index if expected value set for componentexpected_max_displacement – [out] Expected max displacement for component
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelGetExpectedFaceSurfaceForce(Ratel ratel, PetscInt face, PetscInt *num_components, PetscBool has_expected[], PetscScalar expected_surface_force[])#
Retrieve expected surface forces on face.
Not collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextface – [in] Face number
num_components – [inout] Size of input arrays / number of output components
has_expected – [out] Boolean array,
PETSC_TRUE
at index if expected value set for componentexpected_surface_force – [out] Expected surface forces
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelGetExpectedFaceCentroid(Ratel ratel, PetscInt face, PetscInt *num_components, PetscBool has_expected[], PetscScalar expected_centroid[])#
Retrieve expected centroid of face.
Not collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextface – [in] Face number
num_components – [inout] Size of input arrays / number of output components
has_expected – [out] Boolean array,
PETSC_TRUE
at index if expected value set for componentexpected_centroid – [out] Expected centroid values
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelComputeStrainEnergyError(Ratel ratel, PetscScalar strain_energy, PetscBool *has_expected, PetscScalar *strain_energy_error)#
Compute the error between the expected and final strain energy.
Not collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
Contextstrain_energy – [in] Computed final strain energy
has_expected – [out]
PETSC_TRUE
if expected strain energy providedstrain_energy_error – [out] Computed error in the final strain energy
-
PetscErrorCode RatelComputeMaxDisplacementError(Ratel ratel, PetscInt num_components, const PetscScalar max_displacement[], PetscBool *has_any_expected, PetscScalar max_displacement_error[])#
Compute the error between the expected and final maximum displacements.
Not collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
Contextnum_components – [in] Number of components of the input and output arrays
max_displacement – [in] Computed final maximum displacements
has_any_expected – [out]
PETSC_TRUE
if expected maximum displacement provided for any componentmax_displacement_error – [out] Computed error in the final maximum displacements
-
PetscErrorCode RatelComputeFaceForceErrors(Ratel ratel, PetscInt face, PetscInt num_components, const PetscScalar centroid[], const PetscScalar surface_force[], PetscBool *has_any_expected_centroid, PetscBool *has_any_expected_surface_force, PetscScalar centroid_error[], PetscScalar surface_force_error[])#
Compute the error between the expected and final face forces and centroids.
Not collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
Contextface – [in] Face number
num_components – [in] Number of components of the input and output arrays
centroid – [in] Computed final centroid for face
surface_force – [in] Computed final surface forces for face
has_any_expected_centroid – [out]
PETSC_TRUE
if expected centroid provided for facehas_any_expected_surface_force – [out]
PETSC_TRUE
if expected surface force provided for facecentroid_error – [out] Computed error in the final centroid for face
surface_force_error – [out] Computed error in the final surface force for face
-
PetscErrorCode RatelComputeStrainEnergy(Ratel ratel, Vec U, PetscScalar time, PetscScalar *strain_energy)#
Compute the final strain energy in the computed solution.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextU – [in] Computed solution vector
time – [in] Final time value, or
1.0
forSNES
solutionstrain_energy – [out] Computed strain energy
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelComputeMaxDisplacement(Ratel ratel, Vec U, PetscScalar time, PetscInt *num_components, PetscScalar max_displacement[])#
Compute maximum displacement of solution.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextU – [in] Computed solution vector
time – [in] Final time value, or
1.0
forSNES
solutionnum_components – [inout] Size of output buffer/number of components
max_displacement – [out] Computed max displacements
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelGetSurfaceForceFaces(Ratel ratel, PetscInt *num_faces, const PetscInt **face_numbers)#
List mesh face numbers for diagnostic force computation.
Not collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextnum_faces – [out] Number of faces where forces are computed
face_numbers – [out]
DMPlex
mesh face numbers
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelComputeSurfaceForcesCellToFace(Ratel ratel, Vec U, PetscScalar time, PetscInt *num_components, PetscScalar **surface_forces)#
Compute surface forces on mesh faces.
Note: The component
i
of the force for facef
is at indexsurface_forces[num_comp * f + i]
.Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextU – [in] Computed solution vector
time – [in] Final time value, or
1.0
forSNES
solutionnum_components – [out] Number of force components in output buffer
surface_forces – [out] Computed face forces - caller is responsible for freeing
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelComputeSurfaceForces(Ratel ratel, Vec U, PetscScalar time, PetscInt *num_components, PetscScalar **surface_forces)#
Compute surface forces on mesh faces.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextU – [in] Computed solution vector
time – [in] Final time value, or
1.0
forSNES
solutionnum_components – [out] Number of force components in output buffer
surface_forces – [out] Computed face surface forces - caller is responsible for freeing
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelComputeSurfaceCentroids(Ratel ratel, Vec U, PetscScalar time, PetscInt *num_components, PetscScalar **surface_centroids)#
Compute centroids of mesh faces.
Collective across MPI processes.
- Parameters:
ratel – [in] Ratel context
U – [in] Computed solution vector
time – [in] Final time value, or
1.0
forSNES
solutionnum_components – [out] Number of centroid components in output buffer
surface_centroids – [out] Computed surface centroids - caller is responsible for freeing
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelComputeDiagnosticQuantities(Ratel ratel, Vec U, PetscScalar time, Vec *D)#
Compute diagnostic quantities.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextU – [in] Computed solution vector
time – [in] Final time value, or
1.0
forSNES
solutionD – [out] Computed diagnostic quantities vector
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelTSMonitorStrainEnergy(TS ts, PetscInt steps, PetscReal time, Vec U, void *ctx)#
TSMonitor
function for strain energy.Collective across MPI processes.
- Parameters:
ts – [in]
TS
object to monitorsteps – [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
-
PetscErrorCode RatelTSMonitorSurfaceForceCellToFace(TS ts, PetscInt steps, PetscReal time, Vec U, void *ctx)#
TSMonitor
function for surface forces.Collective across MPI processes.
- Parameters:
ts – [in]
TS
object to monitorsteps – [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
-
PetscErrorCode RatelTSMonitorSurfaceForce(TS ts, PetscInt steps, PetscReal time, Vec U, void *ctx)#
TSMonitor
function for face surface forces.Collective across MPI processes.
- Parameters:
ts – [in]
TS
object to monitorsteps – [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
-
PetscErrorCode RatelTSMonitorDiagnosticQuantities(TS ts, PetscInt steps, PetscReal time, Vec U, void *ctx)#
TSMonitor
function for diagnostic quantities.Collective across MPI processes.
- Parameters:
ts – [in]
TS
object to monitorsteps – [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
-
PetscErrorCode RatelTSMonitorCheckpoint(TS ts, PetscInt steps, PetscReal time, Vec U, void *ctx)#
TSMonitor
function for binary checkpoints.Collective across MPI processes.
- Parameters:
ts – [in]
TS
object to monitorsteps – [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
-
PetscErrorCode RatelCheckpointCtxCreateFromOptions(Ratel ratel, const char opt[], const char text[], RatelCheckpointCtx *ctx, PetscBool *set)#
Create a
RatelCheckpointCtx
object fromPetscOptions
command line handling.Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextopt – [in]
TS
monitor flagtext – [in] Description of the flag
ctx – [out]
RatelCheckpointCtx
object to createset – [out]
PetscBool
indicating if the flag was set
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelCheckpointCtxCreate(Ratel ratel, PetscInt interval, const char file_name[], RatelCheckpointCtx *ctx)#
Create a
RatelCheckpointCtx
object.Not collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextinterval – [in] Checkpoint interval
file_name – [in] Base file name for checkpoint files
ctx – [out]
RatelCheckpointCtx
object to create
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelCheckpointCtxDestroy(RatelCheckpointCtx *ctx)#
Destroy
RatelCheckpointCtx
.Not collective across MPI processes.
- Parameters:
ctx – [inout]
RatelCheckpointCtx
object to destroy
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelViewerCreateFromOptions(Ratel ratel, const char opt[], const char text[], RatelViewer *monitor_viewer, PetscBool *set)#
Create a
RatelViewer
object fromPetscOptions
command line handling.Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextopt – [in]
TS
monitor flagtext – [in] Description of the flag
monitor_viewer – [inout]
RatelViewer
object to createset – [inout]
PetscBool
indicating if the flag was set
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelViewerCreate(Ratel ratel, PetscViewer viewer, PetscViewerFormat viewer_format, RatelViewer *monitor_viewer)#
Create a
RatelViewer
object.Not collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextviewer – [in]
PetscViewer
object to useviewer_format – [in]
PetscViewerFormat
object to usemonitor_viewer – [out]
RatelViewer
object to create
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelViewerDestroy(RatelViewer *monitor_viewer)#
Destroy viewer data for
RatelViewer
.Not collective across MPI processes.
- Parameters:
monitor_viewer – [inout]
RatelViewer
object to destroy
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelTSSetup(Ratel ratel, TS ts)#
Setup default
TS
options, theDM
, and options from command line.Solver management.
Note: Sets
SNES
defaults fromRatelSNESSetup()
.Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextts – [inout]
TS
object to setup
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelTSSetFromOptions(Ratel ratel, TS ts)#
Set additional
Ratel
specificTS
options.Collective across MPI processes.`-ts_monitor_strain_energy viewer:filename.extension` sets TSMonitor function for strain energy `-ts_monitor_surface_force viewer:filename.extension` sets TSMonitor function for surface forces `-ts_monitor_surface_force_cell_to_face viewer:filename.extension` sets TSMonitor function for cell-to-face surface forces `-ts_monitor_diagnostic_quantities viewer:filename.extension` sets TSMonitor function for diagnostic quantities `-ts_monitor_checkpoint viewer:filename.extension` sets TSMonitor function for solution checkpoint binaries
- Parameters:
ratel – [in]
Ratel
contextts – [inout]
TS
object to setup
- Returns:
An error code: 0 - success, otherwise - failure
-
static PetscErrorCode RatelKSPSetup(Ratel ratel, KSP ksp)#
Setup default
KSP
options.Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextksp – [inout]
KSP
object to setup
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelSNESSetup(Ratel ratel, SNES snes)#
Setup default
SNES
options, the DM, and options from command line.Note: Sets default
SNES
line search to critical point method. Note: Sets defaultKSP
toCG
with natural norm for most problems. Sets defaultKSP
toGMRES
with preconditioned norm for problems with multiple active fields or platen boundary conditions. Note: Sets defaultPC
toPCPMG
most problems. Sets defaultPC
toPCGAMG
for problems one material and linear elements. Sets defaultPC
toPCJACOBI
for problems with multiple materials.Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextsnes – [inout]
SNES
object to setup
- Returns:
An error code: 0 - success, otherwise - failure
-
static inline PetscErrorCode RatelSetupZeroInitialCondition(Ratel ratel, Vec U)#
Setup zero initial condition with random noise, if needed.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextU – [out] Global vector to set with zero initial condition
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelTSSetupInitialCondition(Ratel ratel, TS ts, Vec U)#
Setup inital conditions for
TS
based upon command line options.Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextts – [inout]
TS
to setupU – [out] Global vector to set with initial conditions
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelSNESSetupInitialCondition(Ratel ratel, SNES snes, Vec U)#
Setup inital conditions for SNES based upon command line options.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextsnes – [inout]
SNES
to setupU – [out] Global vector to set with initial conditions
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelTSCheckpointFinalSolutionFromOptions(Ratel ratel, TS ts, Vec U)#
Checkpoint final solution for
TS
.Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextts – [in]
TS
to checkpointU – [out] Global vector to checkpoint
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelSNESCheckpointFinalSolutionFromOptions(Ratel ratel, SNES snes, Vec U)#
Checkpoint final solution for
SNES
.Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextsnes – [in]
SNES
to checkpointU – [out] Global vector to checkpoint
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelGetVersion(int *major, int *minor, int *patch, PetscBool *release)#
Get
Ratel
library version info.Ratel
version numbers have the form major.minor.patch. Non-release versions may contain unstable interfaces..Not collective across MPI processes.
The caller may pass NULL for any arguments that are not needed.
See also
- Parameters:
major – [out] Major version of the library
minor – [out] Minor version of the library
patch – [out] Patch (subminor) version of the library
release – [out] True for releases; false for development branches.
- Returns:
An error code: 0 - success, otherwise - failure
-
static PetscErrorCode RatelRegisterLogEvents()#
Register core Ratel log events.
Not collective across MPI processes.
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelInit(MPI_Comm comm, Ratel *ratel)#
Setup
Ratel
context object.Library context management.
Note: This function call initializes the
libCEED
context.Collective across MPI processes.
- Parameters:
comm – [in] MPI communication object
ratel – [out]
Ratel
context object
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelView(Ratel ratel, PetscViewer viewer)#
View a
Ratel
context.Note: This function can be called with the command line option
-ratel_view
.Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
context to viewviewer – [in] Optional visualization context
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelDestroy(Ratel *ratel)#
Destroy a
Ratel
context.Collective across MPI processes.
- Parameters:
ratel – [inout]
Ratel
context object to destroy
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelDMCreate(Ratel ratel, RatelSolverType solver_type, DM *dm)#
Create
DM
withSNES
orTS
hooks fromRatel
context.User DM setup.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextsolver_type – [in] Numerical method to use
dm – [out]
DMPlex
object withSNES
orTS
hooks
- Returns:
An error code: 0 - success, otherwise - failure
-
RATEL_VERSION_MAJOR#
Ratel library version numbering.
-
RATEL_VERSION_GE(major, minor, patch)#
Compile-time check that the the current library version is at least as recent as the specified version.
This macro is typically used in
#if RATEL_VERSION_GE(0, 1, 0) code path that needs at least 0.1.0 #else fallback code for older versions #endif
A non-release version always compares as positive infinity.
See also
- Parameters:
major – Major version
minor – Minor version
patch – Patch (subminor) version