GEOS
SolidMechanicsAugmentedLagrangianContact.hpp
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 
16 /*
17  * SolidMechanicsAugmentedLagrangianContact.hpp
18  *
19  */
20 
21 #ifndef GEOS_PHYSICSSOLVERS_CONTACT_SOLIDMECHANICSAUGMENTEDLAGRANGIANCONTACT_HPP_
22 #define GEOS_PHYSICSSOLVERS_CONTACT_SOLIDMECHANICSAUGMENTEDLAGRANGIANCONTACT_HPP_
23 
24 #include "physicsSolvers/solidMechanics/contact/ContactSolverBase.hpp"
25 
26 namespace geos
27 {
28 
30 {
31 public:
32  SolidMechanicsAugmentedLagrangianContact( const string & name,
33  Group * const parent );
34 
36 
41  static string catalogName()
42  {
43  return "SolidMechanicsAugmentedLagrangianContact";
44  }
48  string getCatalogName() const override { return catalogName(); }
49 
50  virtual void registerDataOnMesh( dataRepository::Group & meshBodies ) override final;
51 
52  virtual void setupDofs( DomainPartition const & domain,
53  DofManager & dofManager ) const override;
54 
55  virtual void setupSystem( DomainPartition & domain,
56  DofManager & dofManager,
58  ParallelVector & rhs,
59  ParallelVector & solution,
60  bool const setSparsity = true ) override final;
61 
62  virtual void setSparsityPattern( DomainPartition & domain,
63  DofManager & dofManager,
64  CRSMatrix< real64, globalIndex > & localMatrix,
65  SparsityPattern< globalIndex > & pattern ) override final;
66 
67  virtual void implicitStepSetup( real64 const & time_n,
68  real64 const & dt,
69  DomainPartition & domain ) override final;
70 
71  virtual void implicitStepComplete( real64 const & time_n,
72  real64 const & dt,
73  DomainPartition & domain ) override final;
74 
75  virtual void assembleSystem( real64 const time,
76  real64 const dt,
77  DomainPartition & domain,
78  DofManager const & dofManager,
79  CRSMatrixView< real64, globalIndex const > const & localMatrix,
80  arrayView1d< real64 > const & localRhs ) override;
81 
82  virtual real64 calculateResidualNorm( real64 const & time_n,
83  real64 const & dt,
84  DomainPartition const & domain,
85  DofManager const & dofManager,
86  arrayView1d< real64 const > const & localRhs ) override;
87 
88  virtual void applySystemSolution( DofManager const & dofManager,
89  arrayView1d< real64 const > const & localSolution,
90  real64 const scalingFactor,
91  real64 const dt,
92  DomainPartition & domain ) override;
93 
94  void updateState( DomainPartition & domain ) override final;
95 
96  virtual bool updateConfiguration( DomainPartition & domain,
97  integer configurationLoopIter ) override final;
98 
99 
107  template< typename LAMBDA >
108  void forFiniteElementOnFractureSubRegions( string const & meshName, LAMBDA && lambda ) const
109  {
110 
111  std::map< string,
112  array1d< localIndex > > const & faceTypesToFaceElements = m_faceTypesToFaceElements.at( meshName );
113 
114  for( const auto & [finiteElementName, faceElementList] : faceTypesToFaceElements )
115  {
116  arrayView1d< localIndex const > const faceElemList = faceElementList.toViewConst();
117 
118  finiteElement::FiniteElementBase const & subRegionFE = *(m_faceTypeToFiniteElements.at( finiteElementName ));
119 
120  lambda( finiteElementName, subRegionFE, faceElemList );
121  }
122 
123  }
124 
132  template< typename LAMBDA >
133  void forFiniteElementOnStickFractureSubRegions( string const & meshName, LAMBDA && lambda ) const
134  {
135 
136  bool const isStickState = true;
137 
138  std::map< string, array1d< localIndex > > const &
139  faceTypesToFaceElements = m_faceTypesToFaceElementsStick.at( meshName );
140 
141  for( const auto & [finiteElementName, faceElementList] : faceTypesToFaceElements )
142  {
143  arrayView1d< localIndex const > const faceElemList = faceElementList.toViewConst();
144 
145  finiteElement::FiniteElementBase const & subRegionFE = *(m_faceTypeToFiniteElements.at( finiteElementName ));
146 
147  lambda( finiteElementName, subRegionFE, faceElemList, isStickState );
148  }
149 
150  }
151 
159  template< typename LAMBDA >
160  void forFiniteElementOnSlipFractureSubRegions( string const & meshName, LAMBDA && lambda ) const
161  {
162 
163  bool const isStickState = false;
164 
165  std::map< string, array1d< localIndex > > const &
166  faceTypesToFaceElements = m_faceTypesToFaceElementsSlip.at( meshName );
167 
168  for( const auto & [finiteElementName, faceElementList] : faceTypesToFaceElements )
169  {
170  arrayView1d< localIndex const > const faceElemList = faceElementList.toViewConst();
171 
172  finiteElement::FiniteElementBase const & subRegionFE = *(m_faceTypeToFiniteElements.at( finiteElementName ));
173 
174  lambda( finiteElementName, subRegionFE, faceElemList, isStickState );
175  }
176 
177  }
178 
185  void updateStickSlipList( DomainPartition const & domain );
186 
192  void createFaceTypeList( DomainPartition const & domain );
193 
199  void createBubbleCellList( DomainPartition & domain ) const;
200 
201 private:
202 
210  void addCouplingNumNonzeros( DomainPartition & domain,
211  DofManager & dofManager,
212  arrayView1d< localIndex > const & rowLengths ) const;
213 
220  void addCouplingSparsityPattern( DomainPartition const & domain,
221  DofManager const & dofManager,
222  SparsityPatternView< globalIndex > const & pattern ) const;
223 
224  void computeTolerances( DomainPartition & domain ) const;
225 
227  std::map< string, std::map< string, array1d< localIndex > > > m_faceTypesToFaceElements;
228 
230  std::map< string, std::map< string, array1d< localIndex > > > m_faceTypesToFaceElementsStick;
231 
233  std::map< string, std::map< string, array1d< localIndex > > > m_faceTypesToFaceElementsSlip;
234 
236  std::map< string, std::unique_ptr< geos::finiteElement::FiniteElementBase > > m_faceTypeToFiniteElements;
237 
238  struct viewKeyStruct : ContactSolverBase::viewKeyStruct
239  {
240 
241  constexpr static char const * normalDisplacementToleranceString() { return "normalDisplacementTolerance"; }
242 
243  constexpr static char const * normalTractionToleranceString() { return "normalTractionTolerance"; }
244 
245  constexpr static char const * slidingToleranceString() { return "slidingTolerance"; }
246 
247  constexpr static char const * dispJumpUpdPenaltyString() { return "dispJumpUpdPenalty"; }
248 
249  constexpr static char const * simultaneousString() { return "simultaneous"; }
250 
251  constexpr static char const * symmetricString() { return "symmetric"; }
252 
253  constexpr static char const * iterativePenaltyNFacString() { return "iterPenaltyN"; }
254 
255  constexpr static char const * iterativePenaltyTFacString() { return "iterPenaltyT"; }
256 
257  constexpr static char const * tolJumpDispNFacString() { return "tolJumpN"; }
258 
259  constexpr static char const * tolJumpDispTFacString() { return "tolJumpT"; }
260 
261  constexpr static char const * tolNormalTracFacString() { return "tolNormalTrac"; }
262 
263  constexpr static char const * tolTauLimitString() { return "tolTauLimit"; }
264 
265  };
266 
269  real64 m_slidingCheckTolerance = 5.e-02;
270 
272  int m_simultaneous = 1;
273 
276  int m_symmetric = 1;
277 
279  real64 m_iterPenaltyNFac = 10.0;
280 
282  real64 m_iterPenaltyTFac = 0.1;
283 
285  real64 m_tolJumpDispNFac = 1.e-07;
286 
288  real64 m_tolJumpDispTFac = 1.e-05;
289 
291  real64 m_tolNormalTracFac = 0.5;
292 
293 };
294 
295 } /* namespace geos */
296 
297 #endif /* GEOS_PHYSICSSOLVERS_CONTACT_SOLIDMECHANICSAUGMENTEDLAGRANGIANCONTACT_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...
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
void createFaceTypeList(DomainPartition const &domain)
Create the list of finite elements of the same type for each FaceElementSubRegion (Triangle or Quadri...
void createBubbleCellList(DomainPartition &domain) const
Create the list of elements belonging to CellElementSubRegion that are enriched with the bubble basis...
virtual void setSparsityPattern(DomainPartition &domain, DofManager &dofManager, CRSMatrix< real64, globalIndex > &localMatrix, SparsityPattern< globalIndex > &pattern) override final
Set the sparsity pattern of the linear system matrix.
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 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 setupDofs(DomainPartition const &domain, DofManager &dofManager) const override
Populate degree-of-freedom manager with fields relevant to this solver.
virtual void implicitStepComplete(real64 const &time_n, real64 const &dt, DomainPartition &domain) override final
perform cleanup for implicit timestep
virtual void registerDataOnMesh(dataRepository::Group &meshBodies) override final
Register wrappers that contain data on the mesh objects.
virtual bool updateConfiguration(DomainPartition &domain, integer configurationLoopIter) override final
updates the configuration (if needed) based on the state after a converged Newton loop.
void forFiniteElementOnFractureSubRegions(string const &meshName, LAMBDA &&lambda) const
Loop over the finite element type on the fracture subregions of meshName and apply callback.
void forFiniteElementOnSlipFractureSubRegions(string const &meshName, LAMBDA &&lambda) const
Loop over the finite element type on the slip fracture subregions of meshName and apply callback.
void updateState(DomainPartition &domain) override final
Recompute all dependent quantities from primary variables (including constitutive models)
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 void implicitStepSetup(real64 const &time_n, real64 const &dt, DomainPartition &domain) override final
function to perform setup for implicit timestep
void updateStickSlipList(DomainPartition const &domain)
Create the list of finite elements of the same type for each FaceElementSubRegion (Triangle or Quadri...
void forFiniteElementOnStickFractureSubRegions(string const &meshName, LAMBDA &&lambda) const
Loop over the finite element type on the stick fracture subregions of meshName and apply callback.
Base class for FEM element implementations.
ArrayView< T, 1 > arrayView1d
Alias for 1D array view.
Definition: DataTypes.hpp:179
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
std::string string
String type.
Definition: DataTypes.hpp:90
LvArray::SparsityPatternView< COL_INDEX, INDEX_TYPE const, LvArray::ChaiBuffer > SparsityPatternView
Alias for Sparsity pattern View.
Definition: DataTypes.hpp:301
LvArray::SparsityPattern< COL_INDEX, INDEX_TYPE, LvArray::ChaiBuffer > SparsityPattern
Alias for Sparsity pattern class.
Definition: DataTypes.hpp:297
double real64
64-bit floating point type.
Definition: DataTypes.hpp:98
LvArray::CRSMatrixView< T, COL_INDEX, INDEX_TYPE const, LvArray::ChaiBuffer > CRSMatrixView
Alias for CRS Matrix View.
Definition: DataTypes.hpp:309
mapBase< TKEY, TVAL, std::integral_constant< bool, true > > map
Ordered map type.
Definition: DataTypes.hpp:339
int integer
Signed integer type.
Definition: DataTypes.hpp:81
Array< T, 1 > array1d
Alias for 1D array.
Definition: DataTypes.hpp:175
LAInterface::ParallelVector ParallelVector
Alias for ParallelVector.