Boundary Conditions#
These functions add boundary terms to the residual evaluation.
-
typedef PetscErrorCode RatelBCFunction(PetscInt, PetscReal, const PetscReal*, PetscInt, PetscScalar*, void*)#
Function signature for BC function to set on DMPlex face.
-
int SetupClampBCs(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Setup coordinate and scaling data for Dirichlet (clamp) 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
-
int ApplyClampBCs(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute Dirichlet (clamp) boundary values.
- Parameters:
ctx – [in] QFunction context, holding
RatelBCClampData
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
-
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, returns
s = (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. If type isRATEL_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 ProjectOntoTangentBall(const CeedScalar n[3], CeedScalar radius, const CeedScalar v[3], CeedScalar v_out[3])#
-
int ApplyPlatenBCs(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
(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 platen, andgap = (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
RatelBCPressureData
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 coordinatesout – [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 ApplyPlatenBCs_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
(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
RatelBCPressureData
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 residualout – [out] Output array 0 - action of QFunction
- Returns:
An error code: 0 - success, otherwise - failure
-
int ApplyPressureBCs(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
RatelBCPressureData
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 ApplyPressureBCsJacobian(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
RatelBCPressureData
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 ApplySlipBCs(void *ctx, CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out)#
Compute Dirichlet (slip) boundary values.
- Parameters:
ctx – [in] QFunction context, holding
RatelBCSlipData
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 ApplyTractionBCs(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
RatelBCTractionData
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
-
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, RatelBCSlipData *slip_data)#
Read slip boundary conditions from options database.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contexti – [in] Field index
j – [in] Boundary index
slip_data – [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 RatelBoundaryClampDataFromOptions(Ratel ratel, PetscInt i, PetscInt j, RatelBCClampData *clamp_data)#
Read Dirichlet boundary conditions from options database.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contexti – [in] Field index
j – [in] Boundary index
clamp_data – [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, RatelBCTractionData *traction_data)#
Read traction boundary conditions from options database.
Collective across MPI processes.
- Parameters:
ratel – [in]
Ratel
contexti – [in] Boundary index
traction_data – [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, 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 contextface_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 RatelMaterialCeedAddBoundariesPlaten(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 RatelMaterialSetupPlatenJacobianMultigridLevel(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 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
-
RATEL_MAX_BC_INTERP_POINTS#
Maximum interpolation nodes for boundary conditions.
-
FLOPS_Platen_without_df1#
-
FLOPS_Jacobian_PressureBC#
-
struct RatelBCClampData#
- #include <boundary-data.h>
Clamp boundary condition context.
Public Members
-
CeedScalar translation[3 * RATEL_MAX_BC_INTERP_POINTS]#
Translation vectors.
-
CeedScalar rotation_axis[3 * 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.
-
CeedScalar time#
Current solver time.
-
CeedScalar translation[3 * RATEL_MAX_BC_INTERP_POINTS]#
-
struct RatelBCSlipData#
- #include <boundary-data.h>
Slip boundary condition context.
Public Members
-
CeedScalar translation[3 * 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[3]#
Components to constrain.
-
CeedInt num_components#
Number of constrained components.
-
CeedScalar time#
Current solver time.
-
CeedScalar translation[3 * RATEL_MAX_BC_INTERP_POINTS]#
-
struct RatelBCTractionData#
- #include <boundary-data.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 RatelBCPressureData#
- #include <pressure-boundary.h>
Pressure boundary condition context.