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 setSparsityPattern( DomainPartition & domain,
61  DofManager & dofManager,
63  SparsityPattern< globalIndex > & pattern ) override;
64 
65  virtual std::unique_ptr< PreconditionerBase< LAInterface > >
66  createPreconditioner( DomainPartition & domain ) const override;
67 
68  virtual void
69  implicitStepSetup( real64 const & time_n,
70  real64 const & dt,
71  DomainPartition & domain ) override final;
72 
73  virtual void
74  implicitStepComplete( real64 const & time_n,
75  real64 const & dt,
76  DomainPartition & domain ) override final;
77 
78  virtual void
79  assembleSystem( real64 const time,
80  real64 const dt,
81  DomainPartition & domain,
82  DofManager const & dofManager,
84  arrayView1d< real64 > const & localRhs ) override;
85 
86  virtual real64
88  real64 const & dt,
89  DomainPartition const & domain,
90  DofManager const & dofManager,
91  arrayView1d< real64 const > const & localRhs ) override;
92 
93  virtual void
94  applySystemSolution( DofManager const & dofManager,
95  arrayView1d< real64 const > const & localSolution,
96  real64 const scalingFactor,
97  real64 const dt,
98  DomainPartition & domain ) override;
99 
100  virtual void
102 
103  void updateState( DomainPartition & domain ) override final;
104 
105  void assembleContact( DomainPartition & domain,
106  DofManager const & dofManager,
107  CRSMatrixView< real64, globalIndex const > const & localMatrix,
108  arrayView1d< real64 > const & localRhs );
109 
110  void assembleForceResidualDerivativeWrtTraction( MeshLevel const & mesh,
111  string_array const & regionNames,
112  DofManager const & dofManager,
113  CRSMatrixView< real64, globalIndex const > const & localMatrix,
114  arrayView1d< real64 > const & localRhs );
115 
116  void assembleTractionResidualDerivativeWrtDisplacementAndTraction( MeshLevel const & mesh,
117  string_array const & regionNames,
118  DofManager const & dofManager,
119  CRSMatrixView< real64, globalIndex const > const & localMatrix,
120  arrayView1d< real64 > const & localRhs );
121 
122  void assembleForceResidualPressureContribution( MeshLevel const & mesh,
123  string_array const & regionNames,
124  DofManager const & dofManager,
125  CRSMatrixView< real64, globalIndex const > const & localMatrix,
126  arrayView1d< real64 > const & localRhs );
127 
128  void assembleStabilization( MeshLevel const & mesh,
129  NumericalMethodsManager const & numericalMethodManager,
130  DofManager const & dofManager,
131  CRSMatrixView< real64, globalIndex const > const & localMatrix,
132  arrayView1d< real64 > const & localRhs );
133 
134  bool resetConfigurationToDefault( DomainPartition & domain ) const override final;
135 
137  integer const configurationLoopIter ) override final;
138 
139  bool isFractureAllInStickCondition( DomainPartition const & domain ) const;
140 
141  void computeRotationMatrices( DomainPartition & domain ) const;
142 
143  void computeTolerances( DomainPartition & domain ) const;
144 
145  void computeFaceNodalArea( localIndex const kf0,
147  ArrayOfArraysView< localIndex const > const & faceToNodeMap,
148  ArrayOfArraysView< localIndex const > const & faceToEdgeMap,
149  arrayView2d< localIndex const > const & edgeToNodeMap,
150  arrayView2d< real64 const > const faceCenters,
151  arrayView2d< real64 const > const faceNormals,
152  arrayView1d< real64 const > const faceAreas,
153  stackArray1d< real64, FaceManager::maxFaceNodes() > & nodalArea ) const;
154 
155  void computeFaceIntegrals( arrayView2d< real64 const, nodes::REFERENCE_POSITION_USD > const & nodesCoords,
156  localIndex const (&faceToNodes)[11],
157  localIndex const (&faceToEdges)[11],
158  localIndex const & numFaceVertices,
159  real64 const & faceArea,
160  real64 const (&faceCenter)[3],
161  real64 const (&faceNormal)[3],
162  arrayView2d< localIndex const > const & edgeToNodes,
163  real64 const & invCellDiameter,
164  real64 const (&cellCenter)[3],
165  stackArray1d< real64, FaceManager::maxFaceNodes() > & basisIntegrals,
166  real64 ( &threeDMonomialIntegrals )[3] ) const;
167 
168  real64 const machinePrecision = std::numeric_limits< real64 >::epsilon();
169 
170  string getStabilizationName() const { return m_stabilizationName; }
171 
172 protected:
173 
174  real64 calculateContactResidualNorm( DomainPartition const & domain,
175  DofManager const & dofManager,
176  arrayView1d< real64 const > const & localRhs );
177 
178  virtual void postInputInitialization() override final;
179 
180  void setMGRStrategy();
181 
182 private:
183  string m_stabilizationName;
184 
185  real64 const m_slidingCheckTolerance = 0.05;
186 
187  real64 m_stabilizationScalingCoefficient = 1.0;
188 
189  static const localIndex m_maxFaceNodes; // Maximum number of nodes on a contact face
190 
191  void computeFaceDisplacementJump( DomainPartition & domain );
192 
193  struct viewKeyStruct : ContactSolverBase::viewKeyStruct
194  {
195  constexpr static char const * stabilizationNameString() { return "stabilizationName"; }
196 
197  constexpr static char const * rotationMatrixString() { return "rotationMatrix"; }
198 
199  constexpr static char const * normalDisplacementToleranceString() { return "normalDisplacementTolerance"; }
200 
201  constexpr static char const * normalTractionToleranceString() { return "normalTractionTolerance"; }
202 
203  constexpr static char const * slidingToleranceString() { return "slidingTolerance"; }
204 
205  constexpr static char const * transMultiplierString() { return "penaltyStiffnessTransMultiplier"; }
206 
207  constexpr static char const * stabilizationScalingCoefficientString() { return "stabilizationScalingCoefficient"; }
208 
209  constexpr static char const * useLocalYieldAccelerationString() { return "useLocalYieldAcceleration"; }
210 
211  constexpr static char const * localYieldAccelerationBufferString() { return "localYieldAccelerationBuffer"; }
212  };
213 
215 
216  integer m_useLocalYieldAcceleration; // flag for applying acceleration to yield
217  real64 m_localYieldAccelerationBuffer; // buffer to control the acceleration
218 
219  void tryLocalYieldAcceleration( FaceElementSubRegion & subRegion,
220  integer configurationLoopIter,
221  localIndex kfe,
222  real64 currentTau_unscaled,
223  real64 limitTau,
224  real64 currentTau,
225  integer & fractureState );
226 
227 };
228 
229 } /* namespace geos */
230 
231 #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...
constexpr static int maxFaceNodes()
}@
Class facilitating the representation of a multi-level discretization of a MeshBody.
Definition: MeshLevel.hpp:42
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)
static string catalogName()
name of the node manager in the object catalog
virtual void setSparsityPattern(DomainPartition &domain, DofManager &dofManager, CRSMatrix< real64, globalIndex > &localMatrix, SparsityPattern< globalIndex > &pattern) override
Set the sparsity pattern of the linear system matrix.
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
LvArray::CRSMatrix< T, COL_INDEX, INDEX_TYPE, LvArray::ChaiBuffer > CRSMatrix
Alias for CRS Matrix class.
Definition: DataTypes.hpp:305
LvArray::SparsityPattern< COL_INDEX, INDEX_TYPE, LvArray::ChaiBuffer > SparsityPattern
Alias for Sparsity pattern class.
Definition: DataTypes.hpp:297
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
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