GEOS
SolidMechanicsLagrangianFEM.hpp
Go to the documentation of this file.
1 /*
2  * ------------------------------------------------------------------------------------------------------------
3  * SPDX-License-Identifier: LGPL-2.1-only
4  *
5  * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC
6  * Copyright (c) 2018-2024 TotalEnergies
7  * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University
8  * Copyright (c) 2023-2024 Chevron
9  * Copyright (c) 2019- GEOS/GEOSX Contributors
10  * All rights reserved
11  *
12  * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details.
13  * ------------------------------------------------------------------------------------------------------------
14  */
15 
20 #ifndef GEOS_PHYSICSSOLVERS_SOLIDMECHANICS_SOLIDMECHANICSLAGRANGIANFEM_HPP_
21 #define GEOS_PHYSICSSOLVERS_SOLIDMECHANICS_SOLIDMECHANICSLAGRANGIANFEM_HPP_
22 
24 #include "common/TimingMacros.hpp"
26 #include "kernels/StrainHelper.hpp"
27 #include "mesh/mpiCommunications/CommunicationTools.hpp"
28 #include "mesh/mpiCommunications/MPI_iCommData.hpp"
30 
32 
33 namespace geos
34 {
35 
42 {
43 public:
44 
46  static string coupledSolverAttributePrefix() { return "solid"; }
47 
54  {
55  QuasiStatic,
58  };
59 
65  SolidMechanicsLagrangianFEM( const string & name,
66  Group * const parent );
67 
68 
71 
72  SolidMechanicsLagrangianFEM & operator=( SolidMechanicsLagrangianFEM const & ) = delete;
74 
78  virtual ~SolidMechanicsLagrangianFEM() override;
79 
83  static string catalogName() { return "SolidMechanicsLagrangianFEM"; }
87  string getCatalogName() const override { return catalogName(); }
88 
89  virtual void initializePreSubGroups() override;
90 
91  virtual void registerDataOnMesh( Group & meshBodies ) override;
92 
99  virtual
100  real64 solverStep( real64 const & time_n,
101  real64 const & dt,
102  integer const cycleNumber,
103  DomainPartition & domain ) override;
104 
105  virtual
106  real64 explicitStep( real64 const & time_n,
107  real64 const & dt,
108  integer const cycleNumber,
109  DomainPartition & domain ) override;
110 
111  virtual void
112  implicitStepSetup( real64 const & time_n,
113  real64 const & dt,
114  DomainPartition & domain ) override;
115 
116  virtual void
117  setupDofs( DomainPartition const & domain,
118  DofManager & dofManager ) const override;
119 
120  virtual void
122  DofManager & dofManager,
123  CRSMatrix< real64, globalIndex > & localMatrix,
124  ParallelVector & rhs,
125  ParallelVector & solution,
126  bool const setSparsity = false ) override;
127 
128  virtual void
129  assembleSystem( real64 const time,
130  real64 const dt,
131  DomainPartition & domain,
132  DofManager const & dofManager,
133  CRSMatrixView< real64, globalIndex const > const & localMatrix,
134  arrayView1d< real64 > const & localRhs ) override;
135 
136  virtual void
137  applySystemSolution( DofManager const & dofManager,
138  arrayView1d< real64 const > const & localSolution,
139  real64 const scalingFactor,
140  real64 const dt,
141  DomainPartition & domain ) override;
142 
143  virtual void updateState( DomainPartition & domain ) override
144  {
145  // There should be nothing to update
146  GEOS_UNUSED_VAR( domain );
147  };
148 
149  virtual void applyBoundaryConditions( real64 const time,
150  real64 const dt,
151  DomainPartition & domain,
152  DofManager const & dofManager,
153  CRSMatrixView< real64, globalIndex const > const & localMatrix,
154  arrayView1d< real64 > const & localRhs ) override;
155 
156  virtual real64
157  calculateResidualNorm( real64 const & time_n,
158  real64 const & dt,
159  DomainPartition const & domain,
160  DofManager const & dofManager,
161  arrayView1d< real64 const > const & localRhs ) override;
162 
163  virtual void resetStateToBeginningOfStep( DomainPartition & domain ) override;
164 
165  virtual void implicitStepComplete( real64 const & time,
166  real64 const & dt,
167  DomainPartition & domain ) override;
168 
172  template< typename TYPE_LIST,
173  typename KERNEL_WRAPPER,
174  typename ... PARAMS >
175  real64 assemblyLaunch( MeshLevel & mesh,
176  DofManager const & dofManager,
177  string_array const & regionNames,
178  string const & materialNamesString,
179  CRSMatrixView< real64, globalIndex const > const & localMatrix,
180  arrayView1d< real64 > const & localRhs,
181  real64 const dt,
182  PARAMS && ... params );
183 
184  real64 explicitKernelDispatch( MeshLevel & mesh,
185  string_array const & targetRegions,
186  string const & finiteElementName,
187  real64 const dt,
188  std::string const & elementListName );
189 
200  DofManager const & dofManager,
201  DomainPartition & domain,
202  CRSMatrixView< real64, globalIndex const > const & localMatrix,
203  arrayView1d< real64 > const & localRhs );
204 
205  void applyTractionBC( real64 const time,
206  DofManager const & dofManager,
207  DomainPartition & domain,
208  arrayView1d< real64 > const & localRhs );
209 
210  void applyChomboPressure( DofManager const & dofManager,
211  DomainPartition & domain,
212  arrayView1d< real64 > const & localRhs );
213 
214 
215  void applyContactConstraint( DofManager const & dofManager,
216  DomainPartition & domain,
217  CRSMatrixView< real64, globalIndex const > const & localMatrix,
218  arrayView1d< real64 > const & localRhs );
219 
220  virtual real64
222  DofManager const & dofManager,
223  arrayView1d< real64 const > const & localSolution ) override;
224 
225  void enableFixedStressPoromechanicsUpdate();
226 
227  virtual void saveSequentialIterationState( DomainPartition & domain ) override;
228 
230  {
231  static constexpr char const * newmarkGammaString() { return "newmarkGamma"; }
232  static constexpr char const * newmarkBetaString() { return "newmarkBeta"; }
233  static constexpr char const * massDampingString() { return "massDamping"; }
234  static constexpr char const * stiffnessDampingString() { return "stiffnessDamping"; }
235  static constexpr char const * timeIntegrationOptionString() { return "timeIntegrationOption"; }
236  static constexpr char const * maxNumResolvesString() { return "maxNumResolves"; }
237  static constexpr char const * strainTheoryString() { return "strainTheory"; }
238  static constexpr char const * solidMaterialNamesString() { return "solidMaterialNames"; }
239  static constexpr char const * contactRelationNameString() { return "contactRelationName"; }
240  static constexpr char const * noContactRelationNameString() { return "NOCONTACT"; }
241  static constexpr char const * maxForceString() { return "maxForce"; }
242  static constexpr char const * elemsAttachedToSendOrReceiveNodesString() { return "elemsAttachedToSendOrReceiveNodes"; }
243  static constexpr char const * elemsNotAttachedToSendOrReceiveNodesString() { return "elemsNotAttachedToSendOrReceiveNodes"; }
244  static constexpr char const * surfaceGeneratorNameString() { return "surfaceGeneratorName"; }
245 
246  static constexpr char const * sendOrReceiveNodesString() { return "sendOrReceiveNodes";}
247  static constexpr char const * nonSendOrReceiveNodesString() { return "nonSendOrReceiveNodes";}
248  static constexpr char const * targetNodesString() { return "targetNodes";}
249  static constexpr char const * forceString() { return "Force";}
250 
251  static constexpr char const * contactPenaltyStiffnessString() { return "contactPenaltyStiffness"; }
252 
253  };
254 
255  SortedArray< localIndex > & getElemsAttachedToSendOrReceiveNodes( ElementSubRegionBase & subRegion )
256  {
257  return subRegion.getReference< SortedArray< localIndex > >( viewKeyStruct::elemsAttachedToSendOrReceiveNodesString() );
258  }
259 
260  SortedArray< localIndex > & getElemsNotAttachedToSendOrReceiveNodes( ElementSubRegionBase & subRegion )
261  {
262  return subRegion.getReference< SortedArray< localIndex > >( viewKeyStruct::elemsNotAttachedToSendOrReceiveNodesString() );
263  }
264 
265  real64 & getMaxForce() { return m_maxForce; }
266  real64 const & getMaxForce() const { return m_maxForce; }
267 
268  arrayView1d< ParallelVector > const & getRigidBodyModes() const
269  {
270  return m_rigidBodyModes;
271  }
272 
273  array1d< ParallelVector > & getRigidBodyModes()
274  {
275  return m_rigidBodyModes;
276  }
277 
278  /*
279  * @brief Utility function to set the stress initialization flag
280  * @param[in] performStressInitialization true if the solver has to initialize stress, false otherwise
281  */
282  void setStressInitialization( bool const performStressInitialization )
283  {
284  m_performStressInitialization = performStressInitialization;
285  }
286 
287 protected:
288  virtual void postInputInitialization() override;
289 
291 
292  virtual void setConstitutiveNamesCallSuper( ElementSubRegionBase & subRegion ) const override;
293 
294  real64 m_newmarkGamma;
295  real64 m_newmarkBeta;
296  real64 m_massDamping;
297  real64 m_stiffnessDamping;
298  TimeIntegrationOption m_timeIntegrationOption;
299  real64 m_maxForce = 0.0;
300  integer m_maxNumResolves;
301  integer m_strainTheory;
302 
307 
310 
311  real64 m_contactPenaltyStiffness;
312 
313 private:
314 
315  string m_contactRelationName;
316 
317  PhysicsSolverBase *m_surfaceGenerator;
318  string m_surfaceGeneratorName;
319 };
320 
322  "QuasiStatic",
323  "ImplicitDynamic",
324  "ExplicitDynamic" );
325 
326 //**********************************************************************************************************************
327 //**********************************************************************************************************************
328 //**********************************************************************************************************************
329 
330 
331 template< typename TYPE_LIST,
332  typename KERNEL_WRAPPER,
333  typename ... PARAMS >
334 real64 SolidMechanicsLagrangianFEM::assemblyLaunch( MeshLevel & mesh,
335  DofManager const & dofManager,
336  string_array const & regionNames,
337  string const & materialNamesString,
338  CRSMatrixView< real64, globalIndex const > const & localMatrix,
339  arrayView1d< real64 > const & localRhs,
340  real64 const dt,
341  PARAMS && ... params )
342 {
344 
345  NodeManager const & nodeManager = mesh.getNodeManager();
346 
347  string const dofKey = dofManager.getKey( fields::solidMechanics::totalDisplacement::key() );
348  arrayView1d< globalIndex const > const & dofNumber = nodeManager.getReference< globalIndex_array >( dofKey );
349 
350  real64 const gravityVectorData[3] = LVARRAY_TENSOROPS_INIT_LOCAL_3( gravityVector() );
351 
352  KERNEL_WRAPPER kernelWrapper( dofNumber,
353  dofManager.rankOffset(),
354  localMatrix,
355  localRhs,
356  dt,
357  gravityVectorData,
358  std::forward< PARAMS >( params )... );
359 
360  return finiteElement::
361  regionBasedKernelApplication< parallelDevicePolicy< >,
362  TYPE_LIST >( mesh,
363  regionNames,
364  this->getDiscretizationName(),
365  materialNamesString,
366  kernelWrapper );
367 
368 }
369 
370 } /* namespace geos */
371 
372 #endif /* GEOS_PHYSICSSOLVERS_SOLIDMECHANICS_SOLIDMECHANICSLAGRANGIANFEM_HPP_ */
#define GEOS_UNUSED_VAR(...)
Mark an unused variable and silence compiler warnings.
Definition: GeosxMacros.hpp:84
#define GEOS_MARK_FUNCTION
Mark function with both Caliper and NVTX if enabled.
The DoFManager is responsible for allocating global dofs, constructing sparsity patterns,...
Definition: DofManager.hpp:44
globalIndex rankOffset(string const &fieldName) const
string const & getKey(string const &fieldName) const
Return the key used to record the field in the DofManager.
Partition of the decomposed physical domain. It also manages the connexion information to its neighbo...
Class facilitating the representation of a multi-level discretization of a MeshBody.
Definition: MeshLevel.hpp:42
NodeManager const & getNodeManager() const
Get the node manager.
Definition: MeshLevel.hpp:155
The NodeManager class provides an interface to ObjectManagerBase in order to manage node data.
Definition: NodeManager.hpp:46
Base class for all physics solvers.
string getDiscretizationName() const
return the name of the discretization object
R1Tensor const gravityVector() const
return the value of the gravity vector specified in PhysicsSolverManager
bool m_performStressInitialization
Flag to indicate that the solver is going to perform stress initialization.
virtual void registerDataOnMesh(Group &meshBodies) override
Register wrappers that contain data on the mesh objects.
bool m_isFixedStressPoromechanicsUpdate
Flag to indicate that the solver is running in fixed stress (sequential) mode.
virtual ~SolidMechanicsLagrangianFEM() override
virtual void initializePreSubGroups() override
Called by Initialize() prior to initializing sub-Groups.
virtual void saveSequentialIterationState(DomainPartition &domain) override
Save the state of the solver for sequential iteration.
SolidMechanicsLagrangianFEM(const string &name, Group *const parent)
virtual void initializePostInitialConditionsPreSubGroups() override
Called by InitializePostInitialConditions() prior to initializing sub-Groups.
void applyDisplacementBCImplicit(real64 const time, DofManager const &dofManager, DomainPartition &domain, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs)
array1d< ParallelVector > m_rigidBodyModes
Rigid body modes.
virtual void setConstitutiveNamesCallSuper(ElementSubRegionBase &subRegion) const override
This function sets constitutive name fields on an ElementSubRegionBase, and calls the base function i...
virtual real64 scalingForSystemSolution(DomainPartition &domain, DofManager const &dofManager, arrayView1d< real64 const > const &localSolution) override
Function to determine if the solution vector should be scaled back in order to maintain a known const...
virtual void postInputInitialization() override
static string coupledSolverAttributePrefix()
String used to form the solverName used to register single-physics solvers in CoupledSolver.
GEOS_DECLTYPE_AUTO_RETURN getReference(LOOKUP_TYPE const &lookup) const
Look up a wrapper and get reference to wrapped object.
Definition: Group.hpp:1275
virtual void applySystemSolution(DofManager const &dofManager, arrayView1d< real64 const > const &localSolution, real64 const scalingFactor, real64 const dt, DomainPartition &domain) override
Function to apply the solution vector to the state.
virtual real64 calculateResidualNorm(real64 const &time_n, real64 const &dt, DomainPartition const &domain, DofManager const &dofManager, arrayView1d< real64 const > const &localRhs) override
calculate the norm of the global system residual
virtual void implicitStepComplete(real64 const &time, real64 const &dt, DomainPartition &domain) override
perform cleanup for implicit timestep
virtual real64 explicitStep(real64 const &time_n, real64 const &dt, integer const cycleNumber, DomainPartition &domain) override
Entry function for an explicit time integration step.
virtual void setupSystem(DomainPartition &domain, DofManager &dofManager, CRSMatrix< real64, globalIndex > &localMatrix, ParallelVector &rhs, ParallelVector &solution, bool const setSparsity=false) override
Set up the linear system (DOF indices and sparsity patterns)
virtual void assembleSystem(real64 const time, real64 const dt, DomainPartition &domain, DofManager const &dofManager, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs) override
function to assemble the linear system matrix and rhs
virtual void applyBoundaryConditions(real64 const time, real64 const dt, DomainPartition &domain, DofManager const &dofManager, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs) override
apply boundary condition to system
virtual void resetStateToBeginningOfStep(DomainPartition &domain) override
reset state of physics back to the beginning of the step.
virtual void updateState(DomainPartition &domain) override
Recompute all dependent quantities from primary variables (including constitutive models)
virtual real64 solverStep(real64 const &time_n, real64 const &dt, integer const cycleNumber, DomainPartition &domain) override
entry function to perform a solver step
virtual void setupDofs(DomainPartition const &domain, DofManager &dofManager) const override
Populate degree-of-freedom manager with fields relevant to this solver.
virtual void implicitStepSetup(real64 const &time_n, real64 const &dt, DomainPartition &domain) override
function to perform setup for implicit timestep
ArrayView< T, 1 > arrayView1d
Alias for 1D array view.
Definition: DataTypes.hpp:188
stdVector< string > string_array
A 1-dimensional array of geos::string types.
Definition: DataTypes.hpp:401
LvArray::CRSMatrixView< T, COL_INDEX, localIndex const, LvArray::ChaiBuffer > CRSMatrixView
Alias for CRS Matrix View.
Definition: DataTypes.hpp:318
array1d< globalIndex > globalIndex_array
A 1-dimensional array of geos::globalIndex types.
Definition: DataTypes.hpp:410
ENUM_STRINGS(LinearSolverParameters::SolverType, "direct", "cg", "gmres", "fgmres", "bicgstab", "preconditioner")
Declare strings associated with enumeration values.
std::string string
String type.
Definition: DataTypes.hpp:90
double real64
64-bit floating point type.
Definition: DataTypes.hpp:98
std::int32_t integer
Signed integer type.
Definition: DataTypes.hpp:81
LvArray::CRSMatrix< T, COL_INDEX, localIndex, LvArray::ChaiBuffer > CRSMatrix
Alias for CRS Matrix class.
Definition: DataTypes.hpp:314
LvArray::SortedArray< T, localIndex, LvArray::ChaiBuffer > SortedArray
A sorted array of local indices.
Definition: DataTypes.hpp:275
Array< T, 1 > array1d
Alias for 1D array.
Definition: DataTypes.hpp:184
LAInterface::ParallelVector ParallelVector
Alias for ParallelVector.
Structure to hold scoped key names.