GEOS
SolidMechanicsLagrangeContact.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 
21 #ifndef GEOS_PHYSICSSOLVERS_CONTACT_SOLIDMECHANICSLAGRANGECONTACT_HPP_
22 #define GEOS_PHYSICSSOLVERS_CONTACT_SOLIDMECHANICSLAGRANGECONTACT_HPP_
23 
24 #include "physicsSolvers/solidMechanics/contact/ContactSolverBase.hpp"
25 
26 namespace geos
27 {
28 
29 class NumericalMethodsManager;
30 
32 {
33 public:
34  SolidMechanicsLagrangeContact( const string & name,
35  Group * const parent );
36 
38 
43  static string catalogName()
44  {
45  return "SolidMechanicsLagrangeContact";
46  }
50  string getCatalogName() const override { return catalogName(); }
51 
52  virtual void initializePreSubGroups() override;
53 
54  virtual void registerDataOnMesh( Group & MeshBodies ) override final;
55 
56  virtual void
57  setupDofs( DomainPartition const & domain,
58  DofManager & dofManager ) const override;
59 
60  virtual void
62  DofManager & dofManager,
64  ParallelVector & rhs,
65  ParallelVector & solution,
66  bool const setSparsity = true ) override final;
67 
68  virtual std::unique_ptr< PreconditionerBase< LAInterface > >
69  createPreconditioner( DomainPartition & domain ) const override;
70 
71  virtual void
72  implicitStepSetup( real64 const & time_n,
73  real64 const & dt,
74  DomainPartition & domain ) override final;
75 
76  virtual void
77  implicitStepComplete( real64 const & time_n,
78  real64 const & dt,
79  DomainPartition & domain ) override final;
80 
81  virtual void
82  assembleSystem( real64 const time,
83  real64 const dt,
84  DomainPartition & domain,
85  DofManager const & dofManager,
86  CRSMatrixView< real64, globalIndex const > const & localMatrix,
87  arrayView1d< real64 > const & localRhs ) override;
88 
89  virtual real64
91  real64 const & dt,
92  DomainPartition const & domain,
93  DofManager const & dofManager,
94  arrayView1d< real64 const > const & localRhs ) override;
95 
96  virtual void
97  applySystemSolution( DofManager const & dofManager,
98  arrayView1d< real64 const > const & localSolution,
99  real64 const scalingFactor,
100  real64 const dt,
101  DomainPartition & domain ) override;
102 
103  virtual void
105 
106  void updateState( DomainPartition & domain ) override final;
107 
108  void assembleContact( DomainPartition & domain,
109  DofManager const & dofManager,
110  CRSMatrixView< real64, globalIndex const > const & localMatrix,
111  arrayView1d< real64 > const & localRhs );
112 
113  void assembleForceResidualDerivativeWrtTraction( MeshLevel const & mesh,
114  string_array const & regionNames,
115  DofManager const & dofManager,
116  CRSMatrixView< real64, globalIndex const > const & localMatrix,
117  arrayView1d< real64 > const & localRhs );
118 
119  void assembleTractionResidualDerivativeWrtDisplacementAndTraction( MeshLevel const & mesh,
120  string_array const & regionNames,
121  DofManager const & dofManager,
122  CRSMatrixView< real64, globalIndex const > const & localMatrix,
123  arrayView1d< real64 > const & localRhs );
124 
125  void assembleForceResidualPressureContribution( MeshLevel const & mesh,
126  string_array const & regionNames,
127  DofManager const & dofManager,
128  CRSMatrixView< real64, globalIndex const > const & localMatrix,
129  arrayView1d< real64 > const & localRhs );
130 
131  void assembleStabilization( MeshLevel const & mesh,
132  NumericalMethodsManager const & numericalMethodManager,
133  DofManager const & dofManager,
134  CRSMatrixView< real64, globalIndex const > const & localMatrix,
135  arrayView1d< real64 > const & localRhs );
136 
137  bool resetConfigurationToDefault( DomainPartition & domain ) const override final;
138 
140  integer const configurationLoopIter ) override final;
141 
142  bool isFractureAllInStickCondition( DomainPartition const & domain ) const;
143 
144  void computeRotationMatrices( DomainPartition & domain ) const;
145 
146  void computeTolerances( DomainPartition & domain ) const;
147 
148  void computeFaceNodalArea( localIndex const kf0,
149  arrayView2d< real64 const, nodes::REFERENCE_POSITION_USD > const & nodePosition,
150  ArrayOfArraysView< localIndex const > const & faceToNodeMap,
151  ArrayOfArraysView< localIndex const > const & faceToEdgeMap,
152  arrayView2d< localIndex const > const & edgeToNodeMap,
153  arrayView2d< real64 const > const faceCenters,
154  arrayView2d< real64 const > const faceNormals,
155  arrayView1d< real64 const > const faceAreas,
156  stackArray1d< real64, FaceManager::maxFaceNodes() > & nodalArea ) const;
157 
158  void computeFaceIntegrals( arrayView2d< real64 const, nodes::REFERENCE_POSITION_USD > const & nodesCoords,
159  localIndex const (&faceToNodes)[11],
160  localIndex const (&faceToEdges)[11],
161  localIndex const & numFaceVertices,
162  real64 const & faceArea,
163  real64 const (&faceCenter)[3],
164  real64 const (&faceNormal)[3],
165  arrayView2d< localIndex const > const & edgeToNodes,
166  real64 const & invCellDiameter,
167  real64 const (&cellCenter)[3],
168  stackArray1d< real64, FaceManager::maxFaceNodes() > & basisIntegrals,
169  real64 ( &threeDMonomialIntegrals )[3] ) const;
170 
171  real64 const machinePrecision = std::numeric_limits< real64 >::epsilon();
172 
173  string getStabilizationName() const { return m_stabilizationName; }
174 
175 protected:
176 
177  real64 calculateContactResidualNorm( DomainPartition const & domain,
178  DofManager const & dofManager,
179  arrayView1d< real64 const > const & localRhs );
180 
181  virtual void postInputInitialization() override final;
182 
183  void setMGRStrategy();
184 
185 private:
186  string m_stabilizationName;
187 
188  real64 const m_slidingCheckTolerance = 0.05;
189 
190  real64 m_stabilizationScalingCoefficient = 1.0;
191 
192  static const localIndex m_maxFaceNodes; // Maximum number of nodes on a contact face
193 
194  void computeFaceDisplacementJump( DomainPartition & domain );
195 
196  struct viewKeyStruct : ContactSolverBase::viewKeyStruct
197  {
198  constexpr static char const * stabilizationNameString() { return "stabilizationName"; }
199 
200  constexpr static char const * rotationMatrixString() { return "rotationMatrix"; }
201 
202  constexpr static char const * normalDisplacementToleranceString() { return "normalDisplacementTolerance"; }
203 
204  constexpr static char const * normalTractionToleranceString() { return "normalTractionTolerance"; }
205 
206  constexpr static char const * slidingToleranceString() { return "slidingTolerance"; }
207 
208  constexpr static char const * transMultiplierString() { return "penaltyStiffnessTransMultiplier"; }
209 
210  constexpr static char const * stabilizationScalingCoefficientString() { return "stabilizationScalingCoefficient"; }
211 
212  constexpr static char const * useLocalYieldAccelerationString() { return "useLocalYieldAcceleration"; }
213 
214  constexpr static char const * localYieldAccelerationBufferString() { return "localYieldAccelerationBuffer"; }
215  };
216 
218 
219  integer m_useLocalYieldAcceleration; // flag for applying acceleration to yield
220  real64 m_localYieldAccelerationBuffer; // buffer to control the acceleration
221 
222  void tryLocalYieldAcceleration( FaceElementSubRegion & subRegion,
223  integer configurationLoopIter,
224  localIndex kfe,
225  real64 currentTau_unscaled,
226  real64 limitTau,
227  real64 currentTau,
228  integer & fractureState );
229 
230 };
231 
232 } /* namespace geos */
233 
234 #endif /* GEOS_PHYSICSSOLVERS_CONTACT_SOLIDMECHANICSLAGRANGECONTACT_HPP_ */
The DoFManager is responsible for allocating global dofs, constructing sparsity patterns,...
Definition: DofManager.hpp:45
Partition of the decomposed physical domain. It also manages the connexion information to its neighbo...
The FaceManager class provides an interface to ObjectManagerBase in order to manage face data.
Definition: FaceManager.hpp:44
Class facilitating the representation of a multi-level discretization of a MeshBody.
Definition: MeshLevel.hpp:42
Common interface for preconditioning operators.
virtual void registerDataOnMesh(Group &MeshBodies) override final
Register wrappers that contain data on the mesh objects.
bool updateConfiguration(DomainPartition &domain, integer const configurationLoopIter) override final
updates the configuration (if needed) based on the state after a converged Newton loop.
virtual void setupDofs(DomainPartition const &domain, DofManager &dofManager) const override
Populate degree-of-freedom manager with fields relevant to this solver.
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 implicitStepSetup(real64 const &time_n, real64 const &dt, DomainPartition &domain) override final
function to perform setup for implicit timestep
virtual void initializePreSubGroups() override
Called by Initialize() prior to initializing sub-Groups.
bool resetConfigurationToDefault(DomainPartition &domain) const override final
resets the configuration to the default value.
virtual void resetStateToBeginningOfStep(DomainPartition &domain) override
reset state of physics back to the beginning of the step.
virtual std::unique_ptr< PreconditionerBase< LAInterface > > createPreconditioner(DomainPartition &domain) const override
Create a preconditioner for this solver's linear system.
virtual void postInputInitialization() override final
virtual real64 calculateResidualNorm(real64 const &time, real64 const &dt, DomainPartition const &domain, DofManager const &dofManager, arrayView1d< real64 const > const &localRhs) override
calculate the norm of the global system residual
void updateState(DomainPartition &domain) override final
Recompute all dependent quantities from primary variables (including constitutive models)
virtual void setupSystem(DomainPartition &domain, DofManager &dofManager, CRSMatrix< real64, globalIndex > &localMatrix, ParallelVector &rhs, ParallelVector &solution, bool const setSparsity=true) override final
Set up the linear system (DOF indices and sparsity patterns)
static string catalogName()
name of the node manager in the object catalog
virtual void implicitStepComplete(real64 const &time_n, real64 const &dt, DomainPartition &domain) override final
perform cleanup for implicit timestep
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.
ArrayView< T, 1 > arrayView1d
Alias for 1D array view.
Definition: DataTypes.hpp:179
stdVector< string > string_array
A 1-dimensional array of geos::string types.
Definition: DataTypes.hpp:361
GEOS_GLOBALINDEX_TYPE globalIndex
Global index type (for indexing objects across MPI partitions).
Definition: DataTypes.hpp:87
LvArray::CRSMatrix< T, COL_INDEX, INDEX_TYPE, LvArray::ChaiBuffer > CRSMatrix
Alias for CRS Matrix class.
Definition: DataTypes.hpp:305
StackArray< T, 1, MAXSIZE > stackArray1d
Alias for 1D stack array.
Definition: DataTypes.hpp:187
LvArray::ArrayOfArraysView< T, INDEX_TYPE const, CONST_SIZES, LvArray::ChaiBuffer > ArrayOfArraysView
View of array of variable-sized arrays. See LvArray::ArrayOfArraysView for details.
Definition: DataTypes.hpp:285
double real64
64-bit floating point type.
Definition: DataTypes.hpp:98
GEOS_LOCALINDEX_TYPE localIndex
Local index type (for indexing objects within an MPI partition).
Definition: DataTypes.hpp:84
GEOS_CONCAT(GEOS_LA_INTERFACE, Interface) LAInterface
Alias for current interface.
LvArray::CRSMatrixView< T, COL_INDEX, INDEX_TYPE const, LvArray::ChaiBuffer > CRSMatrixView
Alias for CRS Matrix View.
Definition: DataTypes.hpp:309
ArrayView< T, 2, USD > arrayView2d
Alias for 2D array view.
Definition: DataTypes.hpp:195
int integer
Signed integer type.
Definition: DataTypes.hpp:81
LAInterface::ParallelVector ParallelVector
Alias for ParallelVector.