GEOSX
|
Implements kernels for solving quasi-static equilibrium. More...
#include <SolidMechanicsSmallStrainQuasiStaticKernel.hpp>
Classes | |
struct | NoOpFunctors |
Internal struct to provide no-op defaults used in the inclusion of lambda functions into kernel component functions. More... | |
class | StackVariables |
Kernel variables allocated on the stack. More... | |
Public Types | |
using | Base = finiteElement::ImplicitKernelBase< SUBREGION_TYPE, CONSTITUTIVE_TYPE, FE_TYPE, 3, 3 > |
Alias for the base class;. | |
Public Types inherited from geosx::finiteElement::ImplicitKernelBase< SUBREGION_TYPE, CONSTITUTIVE_TYPE, FE_TYPE, 3, 3 > | |
using | Base = KernelBase< SUBREGION_TYPE, CONSTITUTIVE_TYPE, FE_TYPE, NUM_DOF_PER_TEST_SP, NUM_DOF_PER_TRIAL_SP > |
Alias for the base class. (i.e. geosx::finiteElement::KernelBase) | |
Public Member Functions | |
QuasiStatic (NodeManager const &nodeManager, EdgeManager const &edgeManager, FaceManager const &faceManager, localIndex const targetRegionIndex, SUBREGION_TYPE const &elementSubRegion, FE_TYPE const &finiteElementSpace, CONSTITUTIVE_TYPE &inputConstitutiveType, arrayView1d< globalIndex const > const &inputDofNumber, globalIndex const rankOffset, CRSMatrixView< real64, globalIndex const > const &inputMatrix, arrayView1d< real64 > const &inputRhs, real64 const (&inputGravityVector)[3]) | |
Constructor. More... | |
GEOSX_HOST_DEVICE GEOSX_FORCE_INLINE void | setup (localIndex const k, StackVariables &stack) const |
Copy global values from primary field to a local stack array. More... | |
template<typename STRESS_MODIFIER = NoOpFunctors> | |
GEOSX_HOST_DEVICE GEOSX_FORCE_INLINE void | quadraturePointKernel (localIndex const k, localIndex const q, StackVariables &stack, STRESS_MODIFIER &&stressModifier=NoOpFunctors{}) const |
Performs a state update at a quadrature point. More... | |
GEOSX_HOST_DEVICE GEOSX_FORCE_INLINE real64 | complete (localIndex const k, StackVariables &stack) const |
Performs the complete phase for the kernel. More... | |
Public Member Functions inherited from geosx::finiteElement::ImplicitKernelBase< SUBREGION_TYPE, CONSTITUTIVE_TYPE, FE_TYPE, 3, 3 > | |
ImplicitKernelBase (NodeManager const &nodeManager, EdgeManager const &edgeManager, FaceManager const &faceManager, localIndex const targetRegionIndex, SUBREGION_TYPE const &elementSubRegion, FE_TYPE const &finiteElementSpace, CONSTITUTIVE_TYPE &inputConstitutiveType, arrayView1d< globalIndex const > const &inputDofNumber, globalIndex const rankOffset, CRSMatrixView< real64, globalIndex const > const &inputMatrix, arrayView1d< real64 > const &inputRhs) | |
Constructor. More... | |
GEOSX_HOST_DEVICE GEOSX_FORCE_INLINE void | setup (localIndex const k, StackVariables &stack) const |
Performs the setup phase for the kernel. More... | |
Public Member Functions inherited from geosx::finiteElement::KernelBase< SUBREGION_TYPE, CONSTITUTIVE_TYPE, FE_TYPE, NUM_DOF_PER_TEST_SP, NUM_DOF_PER_TRIAL_SP > | |
KernelBase (SUBREGION_TYPE const &elementSubRegion, FE_TYPE const &finiteElementSpace, CONSTITUTIVE_TYPE &inputConstitutiveType) | |
Constructor. More... | |
GEOSX_HOST_DEVICE GEOSX_FORCE_INLINE void | setup (localIndex const k, StackVariables &stack) const |
Performs the setup phase for the kernel. More... | |
GEOSX_HOST_DEVICE GEOSX_FORCE_INLINE void | quadraturePointKernel (localIndex const k, localIndex const q, StackVariables &stack) const |
Performs a state update at a quadrature point. More... | |
GEOSX_HOST_DEVICE GEOSX_FORCE_INLINE real64 | complete (localIndex const k, StackVariables &stack) const |
Performs the complete phase for the kernel. More... | |
Static Public Attributes | |
static constexpr int | numNodesPerElem = Base::numTestSupportPointsPerElem |
Static Public Attributes inherited from geosx::finiteElement::KernelBase< SUBREGION_TYPE, CONSTITUTIVE_TYPE, FE_TYPE, NUM_DOF_PER_TEST_SP, NUM_DOF_PER_TRIAL_SP > | |
static constexpr int | numTestSupportPointsPerElem = FE_TYPE::numNodes |
static constexpr int | numTrialSupportPointsPerElem = FE_TYPE::numNodes |
static constexpr int | numDofPerTestSupportPoint = NUM_DOF_PER_TEST_SP |
static constexpr int | numDofPerTrialSupportPoint = NUM_DOF_PER_TRIAL_SP |
static constexpr int | numQuadraturePointsPerElem = FE_TYPE::numQuadraturePoints |
Compile time value for the number of quadrature points per element. | |
Protected Attributes | |
arrayView2d< real64 const, nodes::REFERENCE_POSITION_USD > const | m_X |
The array containing the nodal position array. | |
arrayView2d< real64 const, nodes::TOTAL_DISPLACEMENT_USD > const | m_disp |
The rank-global displacement array. | |
arrayView2d< real64 const, nodes::INCR_DISPLACEMENT_USD > const | m_uhat |
The rank-global incremental displacement array. | |
real64 const | m_gravityVector [3] |
The gravity vector. | |
arrayView2d< real64 const > const | m_density |
The rank global density. | |
Protected Attributes inherited from geosx::finiteElement::ImplicitKernelBase< SUBREGION_TYPE, CONSTITUTIVE_TYPE, FE_TYPE, 3, 3 > | |
arrayView1d< globalIndex const > const | m_dofNumber |
The global degree of freedom number. | |
globalIndex const | m_dofRankOffset |
The global rank offset. | |
CRSMatrixView< real64, globalIndex const > const | m_matrix |
The global Jacobian matrix. | |
arrayView1d< real64 > const | m_rhs |
The global residaul vector. | |
Protected Attributes inherited from geosx::finiteElement::KernelBase< SUBREGION_TYPE, CONSTITUTIVE_TYPE, FE_TYPE, NUM_DOF_PER_TEST_SP, NUM_DOF_PER_TRIAL_SP > | |
traits::ViewTypeConst< typename SUBREGION_TYPE::NodeMapType::base_type > const | m_elemsToNodes |
The element to nodes map. | |
arrayView1d< integer const > const | m_elemGhostRank |
The element ghost rank array. | |
CONSTITUTIVE_TYPE::KernelWrapper const | m_constitutiveUpdate |
FE_TYPE const & | m_finiteElementSpace |
Additional Inherited Members | |
Static Public Member Functions inherited from geosx::finiteElement::KernelBase< SUBREGION_TYPE, CONSTITUTIVE_TYPE, FE_TYPE, NUM_DOF_PER_TEST_SP, NUM_DOF_PER_TRIAL_SP > | |
template<typename POLICY , typename KERNEL_TYPE > | |
static real64 | kernelLaunch (localIndex const numElems, KERNEL_TYPE const &kernelComponent) |
Kernel Launcher. More... | |
Implements kernels for solving quasi-static equilibrium.
Define the base interface for implicit finite element kernels. Define the base interface for finite element kernels.
SUBREGION_TYPE | The type of subregion that the kernel will act on. |
CONSTITUTIVE_TYPE | The type of constitutive model present in the subregion. |
NUM_TEST_SUPPORT_POINTS_PER_ELEM | The number of test space support points per element. |
NUM_TRIAL_SUPPORT_POINTS_PER_ELEM | The number of trial space support points per element. |
NUM_DOF_PER_TEST_SP | The number of DOF per test support point. |
NUM_DOF_PER_TRIAL_SP | The number of DOF per trial support point. |
KernelBase defines an interface for implementing finite element kernels that will be callable by the family of kernel launch functions. Specific physics kernels may or may not derive from KernelBase, but must follow the same interface in order to be callable from the generic launching functions.
The template parameters of KernelBase should be duplicated as part of the interface, EXCEPT for NUM_DOF_PER_TEST_SP
and NUM_DOF_PER_TRIAL_SP
. These values should be set internally by the physics solver since each physics discretization will have a constant intrinsic value for these quantities. For example, when solving or the heat equation with scalar temperature as the primary variable at the support point, these will have a value of 1. In contrast, when solving a solid mechanics problem, with vector displacement as the primary variable at the support point, these will have a value of 3. Note that the interface provided by geosx::finiteElement::RegionBasedKernelApplication will construct a kernel assuming only the first 4 template arguments.
Provides a common base for kernels that require the assembly of a system of equations. The types required to assemble the system, such as DOF information, the Matrix and Vector object, etc., are declared and set here.
NUM_NODES_PER_ELEM | The number of nodes per element for the SUBREGION_TYPE . |
UNUSED | An unused parameter since we are assuming that the test and trial space have the same number of support points. |
Implements the KernelBase interface functions required for solving the quasi-static equilibrium equations using one of the "finite element kernel application" functions such as geosx::finiteElement::RegionBasedKernelApplication.
In this implementation, the template parameter NUM_NODES_PER_ELEM
is used in place of both NUM_TEST_SUPPORT_POINTS_PER_ELEM
and NUM_TRIAL_SUPPORT_POINTS_PER_ELEM
, which are assumed to be equal. This results in the UNUSED
template parameter as only the NUM_NODES_PER_ELEM is passed to the ImplicitKernelBase template to form the base class.
Additionally, the number of degrees of freedom per support point for both the test and trial spaces are specified as 3
when specifying the base class.
Definition at line 56 of file SolidMechanicsSmallStrainQuasiStaticKernel.hpp.
|
inline |
Constructor.
Constructor.
nodeManager | Reference to the NodeManager object. |
edgeManager | Reference to the EdgeManager object. |
faceManager | Reference to the FaceManager object. |
targetRegionIndex | Index of the region the subregion belongs to. |
inputDofNumber | The dof number for the primary field. |
rankOffset | dof index offset of current rank |
inputMatrix | Reference to the Jacobian matrix. |
inputRhs | Reference to the RHS vector. Constructor. |
elementSubRegion | Reference to the SUBREGION_TYPE(class template parameter) object. |
inputConstitutiveType | The constitutive object. |
finiteElementSpace | Placeholder for the finite element space object, which currently doesn't do much. |
inputGravityVector | The gravity vector. |
Definition at line 90 of file SolidMechanicsSmallStrainQuasiStaticKernel.hpp.
|
inline |
Performs the complete phase for the kernel.
STACK_VARIABLE_TYPE | The type of StackVariable that holds the stack variables. This is most likely a defined in a type that derives from KernelBase. |
k | The element index. |
stack | The StackVariable object that hold the stack variables. |
The operations typically found in complete are the mapping of the local Jacobian and Residual into the global Jacobian and Residual.
Definition at line 281 of file SolidMechanicsSmallStrainQuasiStaticKernel.hpp.
|
inline |
Performs a state update at a quadrature point.
STACK_VARIABLE_TYPE | The type of StackVariable that holds the stack variables. This is most likely a defined in a type that derives from KernelBase. |
k | The element index. |
q | The quadrature point index. |
stack | The StackVariable object that hold the stack variables. |
The operations found here are the mapping from the support points to the quadrature point, calculation of gradients, etc. From this data the state of the constitutive model is updated if required by the physics package.
STRESS_MODIFIER | Type of optional functor to allow for the modification of stress prior to integration. |
stressModifier | An optional functor to allow for the modification of stress prior to integration. For solid mechanics kernels, the strain increment is calculated, and the constitutive update is called. In addition, the constitutive stiffness stack variable is filled by the constitutive model. |
Definition at line 239 of file SolidMechanicsSmallStrainQuasiStaticKernel.hpp.
|
inline |
Copy global values from primary field to a local stack array.
Performs the setup phase for the kernel.
STACK_VARIABLE_TYPE | The type of StackVariable that holds the stack variables. This is most likely a defined in a type that derives from KernelBase. |
k | The element index. |
stack | The StackVariable object that hold the stack variables. |
The operations typically found in setup are thing such as the collection of global data into local stack storage.
In this implementation, the element local Row and Column DOF stack arrays are filled for when we fill the global matrix and rhs.
For the QuasiStatic implementation, global values from the displacement, incremental displacement, and degree of freedom numbers are placed into element local stack storage.
Definition at line 172 of file SolidMechanicsSmallStrainQuasiStaticKernel.hpp.
|
static |
Number of nodes per element...which is equal to the numTestSupportPointPerElem and numTrialSupportPointPerElem by definition.
Definition at line 73 of file SolidMechanicsSmallStrainQuasiStaticKernel.hpp.