GEOS
HydrofractureSolver.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_MULTIPHYSICS_HYDROFRACTURESOLVER_HPP_
21 #define GEOS_PHYSICSSOLVERS_MULTIPHYSICS_HYDROFRACTURESOLVER_HPP_
22 
26 
27 namespace geos
28 {
29 
30 using dataRepository::Group;
31 
32 template< typename POROMECHANICS_SOLVER = SinglePhasePoromechanics<> >
33 class HydrofractureSolver : public POROMECHANICS_SOLVER
34 {
35 public:
36 
37  using Base = POROMECHANICS_SOLVER;
38  using Base::m_solvers;
39  using Base::m_names;
40  using Base::m_dofManager;
41  using Base::m_localMatrix;
42  using Base::m_rhs;
43  using Base::m_solution;
44  using Base::m_linearSolverParameters;
45 
46  using Base::registerWrapper;
47  using Base::forDiscretizationOnMeshTargets;
48  using Base::getMeshModificationTimestamp;
49  using Base::getSystemSetupTimestamp;
50  using Base::nonlinearImplicitStep;
51  using Base::implicitStepComplete;
52  using Base::getLogLevel;
53  using Base::setSystemSetupTimestamp;
54  using Base::setupDofs;
55  using Base::flowSolver;
56  using Base::solidMechanicsSolver;
57  using Base::assembleElementBasedTerms;
58  using Base::resetStateToBeginningOfStep;
59 
60 
66  HydrofractureSolver( const string & name,
67  Group * const parent );
68 
70  ~HydrofractureSolver() override {}
71 
72  static string catalogName()
73  {
74  // single phase
75  if constexpr ( std::is_same_v< POROMECHANICS_SOLVER, SinglePhasePoromechanics< SinglePhaseBase > > )
76  {
77  return "Hydrofracture";
78  }
79 // // multi phase (TODO)
80 // else if constexpr ( std::is_same_v< POROMECHANICS_SOLVER, MultiphasePoromechanics< CompositionalMultiphaseBase > > )
81 // {
82 // return "MultiphaseHydrofracture";
83 // }
84  }
88  string getCatalogName() const override { return catalogName(); }
89 
91  static string coupledSolverAttributePrefix() { return "poromechanics"; }
92 
100  virtual void registerDataOnMesh( Group & MeshBodies ) override final;
101 
102  virtual void setupCoupling( DomainPartition const & domain,
103  DofManager & dofManager ) const override final;
104 
105  virtual void setSparsityPattern( DomainPartition & domain,
106  DofManager & dofManager,
107  CRSMatrix< real64, globalIndex > & localMatrix,
108  SparsityPattern< globalIndex > & pattern ) override;
109 
110  virtual void implicitStepSetup( real64 const & time_n,
111  real64 const & dt,
112  DomainPartition & domain ) override final;
113 
114  virtual void assembleSystem( real64 const time,
115  real64 const dt,
116  DomainPartition & domain,
117  DofManager const & dofManager,
118  CRSMatrixView< real64, globalIndex const > const & localMatrix,
119  arrayView1d< real64 > const & localRhs ) override;
120 
121  virtual real64 setNextDt( real64 const & currentTime,
122  real64 const & currentDt,
123  DomainPartition & domain ) override;
124 
125  virtual void updateState( DomainPartition & domain ) override final;
126 
127  virtual void implicitStepComplete( real64 const & time_n,
128  real64 const & dt,
129  DomainPartition & domain ) override final;
130 
131  virtual void resetStateToBeginningOfStep( DomainPartition & domain ) override final;
132 
135  void updateHydraulicApertureAndFracturePermeability( DomainPartition & domain );
136 
137  void assembleForceResidualDerivativeWrtPressure( DomainPartition & domain,
138  CRSMatrixView< real64, globalIndex const > const & localMatrix,
139  arrayView1d< real64 > const & localRhs );
140 
141  void assembleFluidMassResidualDerivativeWrtDisplacement( DomainPartition const & domain,
142  CRSMatrixView< real64, globalIndex const > const & localMatrix );
143 
144  std::unique_ptr< CRSMatrix< real64, localIndex > > & getRefDerivativeFluxResidual_dAperture()
145  {
146  return m_derivativeFluxResidual_dAperture;
147  }
148 
149  CRSMatrixView< real64, localIndex const > getDerivativeFluxResidual_dNormalJump()
150  {
151  return m_derivativeFluxResidual_dAperture->toViewConstSizes();
152  }
153 
154  CRSMatrixView< real64 const, localIndex const > getDerivativeFluxResidual_dNormalJump() const
155  {
156  return m_derivativeFluxResidual_dAperture->toViewConst();
157  }
158 
159  enum class InitializationType : integer
160  {
161  Pressure,
162  Displacement,
163  };
164 
165  struct viewKeyStruct : Base::viewKeyStruct
166  {
167  constexpr static char const * surfaceGeneratorNameString() { return "surfaceGeneratorName"; }
168 
169  constexpr static char const * maxNumResolvesString() { return "maxNumResolves"; }
170 
171  constexpr static char const * isMatrixPoroelasticString() { return "isMatrixPoroelastic"; }
172 
173  constexpr static char const * newFractureInitializationTypeString() { return "newFractureInitializationType"; }
174 
175  constexpr static char const * useQuasiNewtonString() { return "useQuasiNewton"; }
176 
177  constexpr static char const * isLaggingFractureStencilWeightsUpdateString() { return "isLaggingFractureStencilWeightsUpdate"; }
178 
179  constexpr static char const * leakoffConstString() {return "leakoffCoefficient"; }
180 
181  constexpr static char const * fractureCreationTimeString() {return "fractureCreationTime"; }
182 
183 #ifdef GEOS_USE_SEPARATION_COEFFICIENT
184  constexpr static char const * separationCoeff0String() { return "separationCoeff0"; }
185  constexpr static char const * apertureAtFailureString() { return "apertureAtFailure"; }
186 #endif
187  };
188 
189 protected:
190 
191  virtual void postInputInitialization() override final;
192 
200  DofManager & dofManager,
201  arrayView1d< localIndex > const & rowLengths ) const;
202 
203 
211  DofManager & dofManager,
212  SparsityPatternView< globalIndex > const & pattern ) const;
213 
214 
215  void setUpDflux_dApertureMatrix( DomainPartition & domain,
216  DofManager const & dofManager,
217  CRSMatrix< real64, globalIndex > & localMatrix );
218 
219  virtual void setMGRStrategy() override;
220 
221 private:
222 
223  virtual real64 fullyCoupledSolverStep( real64 const & time_n,
224  real64 const & dt,
225  integer const cycleNumber,
226  DomainPartition & domain ) override final;
227 
228  void checkRockOnlyMatrix( dataRepository::Group & meshBodies );
229 
230  void assembleFluidLeakSource( real64 time,
231  real64 dt,
232  DomainPartition & domain,
233  DofManager const & dofManager,
234  CRSMatrixView< real64, globalIndex const > const & localMatrix,
235  arrayView1d< real64 > const & localRhs ) const;
240  void initializeNewFractureFields( real64 time, DomainPartition & domain );
241 
242  // name of the contact relation
243  string m_contactRelationName;
244 
246  string m_surfaceGeneratorName;
247 
249  SurfaceGenerator * m_surfaceGenerator;
250 
251  // it is only important for this case.
252  std::unique_ptr< CRSMatrix< real64, localIndex > > m_derivativeFluxResidual_dAperture;
253 
254  integer m_maxNumResolves;
255  integer m_numResolves[2];
256 
257  integer m_isMatrixPoroelastic;
258 
259  // flag to determine which initialization type to use for the new fracture cell
260  InitializationType m_newFractureInitializationType;
261 
262  integer m_useQuasiNewton; // use Quasi-Newton (see https://arxiv.org/abs/2111.00264)
263 
264  // flag to determine whether or not to apply lagging update for the fracture stencil weights
265  integer m_isLaggingFractureStencilWeightsUpdate;
266 
267  // analytical leakoff coefficient
268  real64 m_leakoffCoefficient;
269 };
270 
272  "Pressure",
273  "Displacement" );
274 
275 
276 } /* namespace geos */
277 
278 #endif /* GEOS_PHYSICSSOLVERS_MULTIPHYSICS_HYDROFRACTURESOLVER_HPP_ */
@ Pressure
Pressure in Pascal.
Definition: Units.hpp:67
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...
~HydrofractureSolver() override
Destructor for the class.
HydrofractureSolver(const string &name, Group *const parent)
main constructor for HydrofractureSolver objects
static string coupledSolverAttributePrefix()
String used to form the solverName used to register solvers in CoupledSolver.
void addFluxApertureCouplingSparsityPattern(DomainPartition &domain, DofManager &dofManager, SparsityPatternView< globalIndex > const &pattern) const
virtual void postInputInitialization() override final
string getCatalogName() const override
void addFluxApertureCouplingNNZ(DomainPartition &domain, DofManager &dofManager, arrayView1d< localIndex > const &rowLengths) const
virtual void setupCoupling(DomainPartition const &domain, DofManager &dofManager) const override final
Utility function to set the coupling between degrees of freedom.
virtual real64 setNextDt(real64 const &currentTime, real64 const &currentDt, DomainPartition &domain) override
function to set the next time step size
virtual void resetStateToBeginningOfStep(DomainPartition &domain) override final
reset state of physics back to the beginning of the step.
virtual void implicitStepComplete(real64 const &time_n, real64 const &dt, DomainPartition &domain) override final
perform cleanup for implicit timestep
virtual void registerDataOnMesh(Group &MeshBodies) override final
Register wrappers that contain data on the mesh objects.
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 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 updateState(DomainPartition &domain) override final
Recompute all dependent quantities from primary variables (including constitutive models)
virtual void implicitStepSetup(real64 const &time_n, real64 const &dt, DomainPartition &domain) override final
function to perform setup for implicit timestep
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
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
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
int integer
Signed integer type.
Definition: DataTypes.hpp:81
ENUM_STRINGS(LinearSolverParameters::SolverType, "direct", "cg", "gmres", "fgmres", "bicgstab", "richardson", "preconditioner")
Declare strings associated with enumeration values.