Boundary Conditions#
These functions add boundary terms to the residual evaluation.
-
enum RatelBoundaryType#
Specify boundary condition.
Values:
-
enumerator RATEL_BOUNDARY_CLAMP#
Dirichlet clamped boundary.
-
enumerator RATEL_BOUNDARY_MMS#
Dirichlet MMS boundary.
-
enumerator RATEL_BOUNDARY_TRACTION#
Neumann traction boundary,.
-
enumerator RATEL_BOUNDARY_PLATEN#
Platen boundary integral.
-
enumerator RATEL_BOUNDARY_CLAMP#
-
enum RatelBCInterpolationType#
Specify interpolation type for flexible boundary condition.
Values:
-
enumerator RATEL_BC_INTERP_NONE#
No interpolation, piecewise constant.
-
enumerator RATEL_BC_INTERP_LINEAR#
Linear interpolation.
-
enumerator RATEL_BC_INTERP_NONE#
-
enum RatelPlatenType#
Specify platen type.
Values:
-
enumerator RATEL_PLATEN_NITSCHE#
-
enumerator RATEL_PLATEN_PENALTY#
-
enumerator RATEL_PLATEN_NITSCHE#
-
typedef PetscErrorCode RatelBCFunction(PetscInt, PetscReal, const PetscReal*, PetscInt, PetscScalar*, void*)#
Function signature for BC function to set on DMPlex face.
-
int MMSBCs_CEED_ScalarBPs(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute Dirichlet MMS boundary values for CEED scalar BPs.
- Parameters:
ctx – [in] QFunction context, holding
RatelMMSCEEDBPsParams
Q – [in] Number of quadrature points
in – [in] Input arrays
0 - stored coordinates
1 - stored multiplicity scaling factor
out – [out] Output array
0 - Dirichlet values
- Returns:
An error code: 0 - success, otherwise - failure
-
int MMSBCs_CEED_VectorBPs(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute Dirichlet MMS boundary values for CEED vector BPs.
- Parameters:
ctx – [in] QFunction context, holding
RatelMMSCEEDBPsParams
Q – [in] Number of quadrature points
in – [in] Input arrays
0 - stored coordinates
1 - stored multiplicity scaling factor
out – [out] Output array
0 - Dirichlet values
- Returns:
An error code: 0 - success, otherwise - failure
-
int ClampBCs(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute Dirichlet (clamp) boundary values.
- Parameters:
ctx – [in] QFunction context, holding
RatelBCClampParams
Q – [in] Number of quadrature points
in – [in] Input arrays
0 - stored coordinates
1 - stored multiplicity scaling factor
out – [out] Output array
0 - Dirichlet values
- Returns:
An error code: 0 - success, otherwise - failure
-
int SetupDirichletBCs(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Setup coordinate and scaling data for Dirichlet boundary.
- Parameters:
ctx – [in] QFunction context, unused
Q – [in] Number of quadrature points
in – [in] Input arrays
0 - nodal coordinates
1 - nodal multiplicity
out – [out] Output array
0 - stored coordinates
1 - stored multiplicity scaling factor
- Returns:
An error code: 0 - success, otherwise - failure
-
CeedInt RatelBCInterpGetPreviousKnotIndex(CeedScalar t, CeedInt num_knots, const CeedScalar *knots)#
Compute the index
j \in [0, num_knots]
such thatknots[j] <= t < knots[j+1]
.If
t < knots[0]
, returns -1.- Parameters:
t – [in] Current time
num_knots – [in] Number of knots
knots – [in] Values of knots
- Returns:
A scalar value: the index of the previous knot
-
CeedScalar RatelBCInterpScaleTime(CeedScalar t, CeedInt prev_knot_index, CeedInt num_knots, const CeedScalar *knots)#
Compute rescaled time between knots, \(s \in [0,1)\).
If
0 <= prev_knot_index < num_knots-1
, returnss = (t - knots[prev_knot_index]) / (knots[prev_knot_index + 1])
. Ift
is before the first knot, i.e. prev_knot_index == -1, implicitly inserts a knot at time 0 and returnss = t / knots[0]
. Ift
is past the final knot, returnss = 1
. Ift > knots[num_knots - 1]
, returnspoints[num_knots - 1]
.- Parameters:
t – [in] Current time
prev_knot_index – [in] Index of last knot
num_knots – [in] Number of knots
knots – [in] Values of knots
- Returns:
A scalar value: the rescaled time between knots
-
int RatelBCInterpolate(RatelBCInterpolationType type, CeedScalar t, CeedInt prev_knot_index, CeedInt num_knots, const CeedScalar *knots, CeedInt point_dim, const CeedScalar *points, CeedScalar *out_point)#
Compute the interpolated value between the given points and knots.
For
t > knots[num_knots-1]
, returnspoints[num_knots-1]
.- Parameters:
type – [in] If type is
RATEL_BC_INTERP_NONE
- returns the point corresponding to the first knot after time t.RATEL_BC_INTERP_LINEAR
- returns the result of linearly interpolating between adjacent points.
t – [in] Current time
prev_knot_index – [in] Index of last knot
num_knots – [in] Number of knots
knots – [in] Values of knots
point_dim – [in] Dimension of coordinates
points – [in] Coordinates of knots
out_point – [out] Interpolated value
- Returns:
An error code: 0 - success, otherwise - failure
-
int MMSBCs_Linear(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute Dirichlet MMS boundary values 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 - stored coordinates
1 - stored multiplicity scaling factor
out – [out] Output array
0 - Dirichlet values
- Returns:
An error code: 0 - success, otherwise - failure
-
int MMSBCs_MixedLinear(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute Dirichlet MMS boundary values 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 - stored coordinates
1 - stored multiplicity scaling factor
out – [out] Output array
0 - Dirichlet values
- Returns:
An error code: 0 - success, otherwise - failure
-
int PlatenBCs(void *ctx, CeedInt Q, RatelComputef1 compute_f1, bool has_state_values, bool has_stored_values, CeedInt num_active_field_eval_modes, const CeedScalar *const *in, CeedScalar *const *out)#
Compute the surface integral of the platen contact boundary condition.
Updated to support initial gap using https://hal.archives-ouvertes.fr/hal-00722114/document.
In the frictionless case, compute
where(f1_gamma < 0) ? (f1_gamma * n_platen * (w detNb)) : 0
f1_gamma = f1_n + gamma*gap
f1_n = (f1·n_platen)·n_platen
is the contact pressure between the surface and platengap = (X + u - C)·n_platen
is the gap between the surface and the platen in the current configuration
For linear elasticity
f1 = sigma
, wheresigma
is Cauchy stress tensor. For hyperelastic in initial configurationf1 = P
, whereP
is First Piola Kirchhoff stress. For hyperelastic in current configurationf1 = tau
, where tau is Kirchhoff stress.For the frictional case, see the theory section of the documentation
- Parameters:
ctx – [in] QFunction platen, holding
RatelBCPlatenParams
Q – [in] Number of quadrature points
compute_f1 – [in] Function to compute f1
has_state_values – [in] Boolean flag indicating model state values in residual evaluation
has_stored_values – [in] Boolean flag indicating model stores values in residual evaluation
num_active_field_eval_modes – [in] Number of active field evaluation modes
in – [in] Input arrays
0 - qdata
input_data_offset
-u
input_data_offset + 1
-du/dt
(oru
in static case)input_data_offset + 2
- quadrature point coordinates
out – [out] Output array
output_data_offset
- action of QFunction onu
output_data_offset +
- storedf1_n + gamma*gap
andf1_T - gamma*(du/dt)_T
- Returns:
An error code: 0 - success, otherwise - failure
-
int PlatenBCs_Jacobian(void *ctx, CeedInt Q, RatelComputef1_fwd compute_df1, bool has_state_values, bool has_stored_values, CeedInt num_active_field_eval_modes, const CeedScalar *const *in, CeedScalar *const *out)#
Compute the incremental traction due to platen contact boundary condition.
In the frictionless case compute
where(f1_gamma < 0) ? ((df1_N + gamma * du_N) * n_platen * (w detNb)) : 0
df1_N = (df1·n_platen)·n_platen
is the incremental surface force away from platen.du_N = du·n_platen
is the incremental displacement away from platen.
For linear elasticity
df1 = dsigma
. For hyperelastic in initial configurationdf1 = dP
. For hyperelastic in current configurationdf1 = dtau - tau * grad_du^T
.For the frictional case, see the theory section of the documentation
- Parameters:
ctx – [in] QFunction platen, holding
RatelBCPlatenParams
Q – [in] Number of quadrature points
compute_df1 – [in] Function to compute f1
has_state_values – [in] Boolean flag indicating model state values in residual evaluation
has_stored_values – [in] Boolean flag indicating model stores values in residual evaluation
num_active_field_eval_modes – [in] Number of active fields
in – [in] Input arrays
0 - qdata
input_data_offset
-du
, incremental change to uinput_data_offset + 1
- stored values from residual
out – [out] Output array
0 - action of QFunction
- Returns:
An error code: 0 - success, otherwise - failure
-
int SetupPenaltyPlatens(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute geometric factors for integration, gradient transformations, and coordinate transformations on element faces.
Reference (parent) 2D coordinates are given by
X
and physical (current) 3D coordinates are given byx
. The change of coordinate matrix is given bydxdX_{i,j} = dx_i/dX_j (indicial notation) [3 * 2]
.(N_1, N_2, N_3)
is given by the cross product of the columns ofdxdX_{i,j}
.detNb
is the magnitude of(N_1, N_2, N_3)
.- Parameters:
ctx – [in] QFunction context, unused
Q – [in] Number of quadrature points
in – [in] Input arrays
0 - Jacobian of face coordinates
1 - quadrature weights
out – [out] Output array
0 - qdata,
w
,detNb
, andN
- Returns:
An error code: 0 - success, otherwise - failure
-
int PlatenPenaltyBCs(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute the surface integral for penalty method based contact on the constrained faces.
In the frictionless case, compute
where(gap < 0) ? (w detNb) * gap / epsilon * N : 0
epsilon = sqrt(detNb) / platen->gamma
is the penalty parameter- Parameters:
ctx – [in] QFunction context, holding
RatelBCPlatenParams
Q – [in] Number of quadrature points
in – [in] Input arrays
0 - quadrature data
1 -
u
2 -
du/dt
(oru
in static case)3 - quadrature point coordinates
out – [out] Output array
0 -
v
1 - stored values from residual
- Returns:
An error code: 0 - success, otherwise - failure
-
int PlatenPenaltyBCsJacobian(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute the Jacobian of the surface integral for penalty method based contact on the constrained faces.
- Parameters:
ctx – [in] QFunction context, holding
RatelBCPlatenParams
Q – [in] Number of quadrature points
in – [in] Input arrays
0 - quadrature data
1 - incremental change in
u
2 - stored values from residual
out – [out] Output array
0 - action of QFunction on
du
- Returns:
An error code: 0 - success, otherwise - failure
-
int PressureBCs(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute the surface integral for pressure on the constrained faces.
- Parameters:
ctx – [in] QFunction context, holding
RatelBCPressureParams
Q – [in] Number of quadrature points
in – [in] Input arrays
0 - gradient of coordinates
1 - quadrature weights
2 - gradient of u with respect to reference coordinates
out – [out] Output array
0 -
\int v^T . (pn) ds
- Returns:
An error code: 0 - success, otherwise - failure
-
int PressureBCsJacobian(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute the Jacobian of the surface integral for pressure on the constrained faces.
- Parameters:
ctx – [in] QFunction context, holding
RatelBCPressureParams
Q – [in] Number of quadrature points
in – [in] Input arrays
0 - gradient of coordinates
1 - quadrature weights
2 - gradient of incremental change in
u
out – [out] Output array
0 - action of QFunction on
dug
- Returns:
An error code: 0 - success, otherwise - failure
-
int SetupSlipBCs(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Setup scaling data for Dirichlet (clamp) boundary.
- Parameters:
ctx – [in] QFunction context, unused
Q – [in] Number of quadrature points
in – [in] Input arrays
0 - nodal multiplicity
out – [out] Output array
0 - stored multiplicity scaling factor
- Returns:
An error code: 0 - success, otherwise - failure
-
int SlipBCs(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute Dirichlet (slip) boundary values.
- Parameters:
ctx – [in] QFunction context, holding
RatelBCSlipParams
Q – [in] Number of quadrature points
in – [in] Input arrays
0 - stored multiplicity scaling factor
out – [out] Output array
0 - Dirichlet values
- Returns:
An error code: 0 - success, otherwise - failure
-
int TractionBCs(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute surface integral of the traction vector on the constrained faces.
- Parameters:
ctx – [in] QFunction context, holding
RatelBCTractionParams
Q – [in] Number of quadrature points
in – [in] Input arrays
0 - surface qdata
out – [out] Output array
0 -
v = t * (w detNb)
- Returns:
An error code: 0 - success, otherwise - failure
-
int MMSError_MixedLinear(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute error from true solution for mixed linear elasticity MMS.
- Parameters:
ctx – [in] QFunction context, holding
RatelLinearElasticityParams
Q – [in] Number of quadrature points
in – [in] Input arrays
0 - qdata
1 - quadrature point coordinates
2 - computed solution, displacement field
3 - computed solution, pressure field
out – [out] Output array
0 - error from true solution, displacement field
0 - error from true solution, pressure field
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelCreateBCLabel(DM dm, const char name[])#
Create boundary label.
Not collective across MPI processes.
- Parameters:
dm – [inout]
DM
to add boundary labelname – [in] Name for new boundary label
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelDMAddBoundariesDirichlet(Ratel ratel, DM dm)#
Add Dirichlet boundary conditions to
DM
.Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextdm – [out]
DM
to update with Dirichlet boundaries
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelDMAddBoundariesSlip(Ratel ratel, DM dm)#
Add Dirichlet slip boundaries to
DM
.Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextdm – [out]
DM
to update with dirichlet slip boundaries
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelBoundarySlipDataFromOptions(Ratel ratel, PetscInt i, PetscInt j, RatelBCSlipParams *params_slip)#
Read slip boundary conditions from options database.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contexti – [in] Field index
j – [in] Boundary index
params_slip – [out] Data structure to store slip data
-
PetscErrorCode RatelCeedAddBoundariesDirichletSlip(Ratel ratel, DM dm, CeedOperator op_dirichlet)#
Add Dirichlet slip boundary conditions to Dirichlet boundary condition operator.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextdm – [in]
DM
to use for Dirichlet boundariesop_dirichlet – [out] Composite Dirichlet
CeedOperator
to add sub-operators to
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelCeedAddBoundariesDirichletMMS(Ratel ratel, DM dm, CeedOperator op_dirichlet)#
Add Dirichlet MMS boundary conditions to Dirichlet boundary condition operator.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextdm – [in]
DM
to use for MMS boundariesop_dirichlet – [out] Composite Dirichlet
CeedOperator
to add sub-operators to
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelBoundaryClampDataFromOptions(Ratel ratel, PetscInt i, PetscInt j, RatelBCClampParams *params_clamp)#
Read Dirichlet boundary conditions from options database.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contexti – [in] Field index
j – [in] Boundary index
params_clamp – [out] Data structure to store clamp data
-
PetscErrorCode RatelCeedAddBoundariesDirichletClamp(Ratel ratel, DM dm, CeedOperator op_dirichlet)#
Add Dirichlet clamp boundary conditions to Dirichlet boundary condition operator.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextdm – [in]
DM
to use for Dirichlet boundariesop_dirichlet – [out] Composite Dirichlet
CeedOperator
to add sub-operators to
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelBoundaryTractionDataFromOptions(Ratel ratel, PetscInt i, RatelBCTractionParams *params_traction)#
Read traction boundary conditions from options database.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contexti – [in] Boundary index
params_traction – [out] Data structure to store traction data
-
PetscErrorCode RatelCeedAddBoundariesNeumann(Ratel ratel, DM dm, CeedOperator op_residual)#
Add Neumann boundary conditions to residual operator.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextdm – [in]
DM
to use for Neumann boundariesop_residual – [out] Composite residual u term
CeedOperator
to add Neumann sub-operators
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelSetupSurfaceForceCentroids(Ratel ratel, DM dm)#
Setup surface force face centroid computation.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextdm – [in]
DM
with surface force faces
- Returns:
An error code: 0 - success, otherwise - failure
-
static PetscErrorCode RatelMaterialPlatenContextRegisterFields(RatelMaterial material, CeedQFunctionContext ctx)#
Register material specific fields for platen boundary conditions
CeedQFunctionContext
Not collective across MPI processes.
- Parameters:
material – [in]
RatelMaterial
to setup platen contextctx – [inout]
CeedQFunctionContext
for platen
- Returns:
An error code: 0 - success, otherwise - failure
-
static PetscErrorCode RatelCeedPlatenContextCreate(RatelMaterial material, CeedScalar normal[3], CeedScalar center[3], CeedScalar distance[RATEL_MAX_BC_INTERP_POINTS], CeedScalar times[RATEL_MAX_BC_INTERP_POINTS], PetscInt num_times, RatelBCInterpolationType interpolation, CeedScalar gamma, CeedScalar friction, CeedInt face_id, CeedInt face_orientation, CeedQFunctionContext *ctx)#
Create libCEED
CeedQFunctionContext
for platen boundary conditions.Not collective across MPI processes.
- Parameters:
material – [in]
RatelMaterial
to setup platen contextnormal – [in] Normal vector to platen surface
center – [in] Center of platen surface
distance – [in] Displacement of platen along normal
times – [in] Transition times between distances
num_times – [in] Number of transition times
interpolation – [in] Interpolation type for platen BC
gamma – [in] Nitsche’s method penalty parameter
friction – [in] Coefficient of friction for frictional contact
face_id – [in]
DMPlex
face domain idface_orientation – [in]
DMPlex
face domain stratum valuectx – [out]
CeedQFunctionContext
for platen
- Returns:
An error code: 0 - success, otherwise - failure
-
static PetscErrorCode RatelMaterialPlatenContextCreate(RatelMaterial material, const char *platen_name, PetscInt face_id, PetscInt face_orientation, CeedQFunctionContext *ctx)#
Process command line options for platen boundary conditions.
Collective across MPI processes.
- Parameters:
material – [in]
RatelMaterial
to setup platen contextplaten_name – [in]
const char *
name of platenface_id – [in]
DMPlex
face domain idface_orientation – [in]
DMPlex
face domain stratum valuectx – [out]
CeedQFunctionContext
for platen
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelMaterialCeedAddBoundariesPlatenNitsche(RatelMaterial material, CeedOperator op_residual, CeedVector u_dot_loc, CeedOperator op_jacobian)#
Add platen contact boundary conditions to the residual u term operator.
Collective across MPI processes.
- Parameters:
material – [in]
RatelMaterial
contextop_residual – [out] Composite residual u term
CeedOperator
to add RatelMaterial sub-operatoru_dot_loc – [out]
CeedVector
for passive input velocity fieldop_jacobian – [out] Composite Jacobian
CeedOperator
to addRatelMaterial
sub-operator
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelMaterialSetupPlatenNitscheJacobianMultigridLevel(RatelMaterial material, DM dm_level, CeedVector m_loc, CeedOperator sub_op_jacobian_fine, CeedOperator op_jacobian_coarse, CeedOperator op_prolong, CeedOperator op_restrict)#
Setup multigrid level suboperator for
RatelMaterial
platen boundary conditions Jacobian.Collective across MPI processes.
- Parameters:
material – [in]
RatelMaterial
contextdm_level – [in]
DMPlex
for multigrid level to setupm_loc – [in]
CeedVector
holding multiplicitysub_op_jacobian_fine – [in] Fine grid Jacobian
CeedOperator
op_jacobian_coarse – [inout] Composite Jacobian
CeedOperator
to addRatelMaterial
suboperatorsop_prolong – [inout] Composite prolongation
CeedOperator
to addRatelMaterial
suboperators, unused but needed for compatable function signatureop_restrict – [inout] Composite restriction
CeedOperator
to addRatelMaterial
suboperators, unused but needed for compatable function signature
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelCeedAddBoundariesPlatenPenalty(Ratel ratel, DM dm, CeedOperator op_residual, CeedVector u_dot_loc, CeedOperator op_jacobian)#
-
PetscErrorCode RatelMaterialSetupPlatenPenaltyJacobianMultigridLevel(RatelMaterial material, DM dm_level, CeedVector m_loc, CeedOperator sub_op_jacobian_fine, CeedOperator op_jacobian_coarse, CeedOperator op_prolong, CeedOperator op_restrict)#
-
PetscErrorCode RatelCeedAddBoundariesPressure(Ratel ratel, DM dm, CeedOperator op_residual, CeedOperator op_jacobian)#
Add pressure boundary conditions to residual and Jacobian operators.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextdm – [in]
DM
to use for pressure boundariesop_residual – [out] Composite residual u term
CeedOperator
to add pressure sub-operatorsop_jacobian – [out] Composite Jacobian
CeedOperator
to add pressure sub-operators
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelSetupPressureJacobianMultigridLevel(Ratel ratel, DM dm_level, CeedVector m_loc, CeedOperator op_jacobian_fine, CeedOperator op_jacobian_coarse)#
Setup multigrid level suboperators for
RatelMaterial
pressure boundary conditions Jacobian.Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextdm_level – [in]
DMPlex
for multigrid level to setupm_loc – [in]
CeedVector
holding multiplicityop_jacobian_fine – [in] Fine grid Jacobian
CeedOperator
op_jacobian_coarse – [inout] Composite Jacobian
CeedOperator
to addRatelMaterial
suboperators
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelBoundingBoxParamsFromOptions(Ratel ratel, const char option_prefix[], const char name[], RatelBoundingBoxParams *bounding_box)#
Read bounding box from options.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextoption_prefix – [in] Prefix for options
name – [in] Name of face
bounding_box – [out] Output bounding box
- Returns:
An error code: 0 - success, otherwise - failure
-
PetscErrorCode RatelFaceLabelValueFromOptions(Ratel ratel, const char option_prefix[], const char name[], PetscInt *label_value)#
Read label value from options.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contextoption_prefix – [in] Prefix for options
name – [in] Name of face
label_value – [out] Label value read from options
- Returns:
An error code: 0 - success, otherwise - failure
-
RATEL_MAX_BC_INTERP_POINTS 256#
Maximum interpolation nodes for boundary conditions.
-
FLOPS_Platen_without_df1 (7 * FLOPS_Dot3 + 64)#
-
FLOPS_Jacobian_PlatenPenaltyBC (3 * FLOPS_Dot3 + FLOPS_Norm3 + 57)#
-
FLOPS_Jacobian_PressureBC 30#
-
struct RatelBCClampParams#
- #include <boundary.h>
Clamp boundary condition context.
Public Members
-
CeedScalar translation[4 * RATEL_MAX_BC_INTERP_POINTS]#
Translation vectors.
-
CeedScalar rotation_axis[4 * RATEL_MAX_BC_INTERP_POINTS]#
Rotation axes.
-
CeedScalar rotation_polynomial[2 * RATEL_MAX_BC_INTERP_POINTS]#
Rotation polynomial coefficients.
-
CeedScalar times[RATEL_MAX_BC_INTERP_POINTS]#
Transition times.
-
CeedInt num_times#
Number of transition times.
-
RatelBCInterpolationType interpolation_type#
Type of interpolation between points.
-
CeedInt num_comp_u#
Number of components in field.
-
CeedScalar time#
Current solver time.
-
CeedScalar translation[4 * RATEL_MAX_BC_INTERP_POINTS]#
-
struct RatelBCSlipParams#
- #include <boundary.h>
Slip boundary condition context.
Public Members
-
CeedScalar translation[4 * RATEL_MAX_BC_INTERP_POINTS]#
Translation vectors.
-
CeedScalar times[RATEL_MAX_BC_INTERP_POINTS]#
Transition times.
-
CeedInt num_times#
Number of transition times.
-
RatelBCInterpolationType interpolation_type#
Type of interpolation between points.
-
CeedInt components[4]#
Components to constrain.
-
CeedInt num_comp_slip#
Number of constrained components.
-
CeedInt num_comp_u#
Number of components in field.
-
CeedScalar time#
Current solver time.
-
CeedScalar translation[4 * RATEL_MAX_BC_INTERP_POINTS]#
-
struct RatelBCTractionParams#
- #include <boundary.h>
Traction boundary condition context.
Public Members
-
CeedScalar direction[3 * RATEL_MAX_BC_INTERP_POINTS]#
Traction vector.
-
CeedScalar times[RATEL_MAX_BC_INTERP_POINTS]#
Transition times.
-
CeedInt num_times#
Number of transition times.
-
RatelBCInterpolationType interpolation_type#
Type of interpolation between points.
-
CeedScalar time#
Current solver time.
-
CeedScalar direction[3 * RATEL_MAX_BC_INTERP_POINTS]#
-
struct RatelBoundingBoxParams#
- #include <bounding-box.h>
-
struct RatelBCPlatenParamsCommon#
- #include <platen.h>
Common Platen BC data.
Public Members
-
CeedScalar normal[3]#
Exterior normal for platen.
-
CeedScalar center[3]#
Center of the platen.
-
CeedScalar distance[RATEL_MAX_BC_INTERP_POINTS]#
Displacement of the platen along normal vector.
-
CeedScalar times[RATEL_MAX_BC_INTERP_POINTS]#
Transition times.
-
CeedInt num_times#
Number of transition times.
-
RatelBCInterpolationType interpolation_type#
Type of interpolation between points.
-
CeedScalar gamma#
Nitsche’s method parameter.
-
CeedScalar f#
Coefficient of friction.
-
CeedInt face_id#
DM face id.
-
CeedInt face_domain_value#
DM face orientation.
-
CeedScalar normal[3]#
-
struct RatelBCPlatenParams#
- #include <platen.h>
Platen BC data, including material data.
Public Members
-
RatelBCPlatenParamsCommon platen#
Common platen data.
-
CeedScalar time#
Current solver time.
-
CeedScalar material[RATEL_MAX_MATERIAL_SIZE]#
Material data.
-
RatelBCPlatenParamsCommon platen#
-
struct RatelBCPressureParams#
- #include <pressure.h>
Pressure boundary condition context.