21 #ifndef GEOS_PHYSICSSOLVERS_MULTIPHYSICS_COUPLEDRESERVOIRANDWELLSBASE_HPP_
22 #define GEOS_PHYSICSSOLVERS_MULTIPHYSICS_COUPLEDRESERVOIRANDWELLSBASE_HPP_
27 #include "constitutive/permeability/PermeabilityFields.hpp"
28 #include "constitutive/permeability/PermeabilityBase.hpp"
31 #include "physicsSolvers/fluidFlow/wells/WellControls.hpp"
37 namespace coupledReservoirAndWellsInternal
57 string const & resElemDofName,
58 string const & wellElemDofName );
73 template<
typename RESERVOIR_SOLVER,
typename WELL_SOLVER >
102 :
Base( name, parent ),
105 this->
template getWrapper< string >( Base::viewKeyStruct::discretizationString() ).
123 reservoirSolver()->setSparsityPattern( domain, dofManager, localMatrix, patternDiag );
127 for(
localIndex localRow = 0; localRow < patternDiag.numRows(); ++localRow )
129 rowLengths[localRow] = patternDiag.numNonZeros( localRow );
136 pattern.resizeFromRowCapacities< parallelHostPolicy >( patternDiag.numRows(), patternDiag.numColumns(), rowLengths.data());
139 for(
localIndex localRow = 0; localRow < patternDiag.numRows(); ++localRow )
141 globalIndex const * cols = patternDiag.getColumns( localRow ).dataIfContiguous();
142 pattern.insertNonZeros( localRow, cols, cols + patternDiag.numNonZeros( localRow ));
170 DomainPartition & domain = this->
template getGroupByPath< DomainPartition >(
"/Problem/domain" );
173 if( !validateWellPerforations( domain ))
175 GEOS_ERROR( GEOS_FMT(
"{}: well perforations validation failed, bad perforations found", this->
getName()));
200 computeWellTransmissibility( domain );
208 assembleFluxTerms(
real64 const dt,
209 DomainPartition
const & domain,
210 DofManager
const & dofManager,
211 CRSMatrixView< real64, globalIndex const >
const & localMatrix,
212 arrayView1d< real64 >
const & localRhs )
const
213 {
reservoirSolver()->assembleFluxTerms( dt, domain, dofManager, localMatrix, localRhs ); }
216 assembleStabilizedFluxTerms(
real64 const dt,
217 DomainPartition
const & domain,
218 DofManager
const & dofManager,
219 CRSMatrixView< real64, globalIndex const >
const & localMatrix,
220 arrayView1d< real64 >
const & localRhs )
const
221 {
reservoirSolver()->assembleStabilizedFluxTerms( dt, domain, dofManager, localMatrix, localRhs ); }
223 real64 updateFluidState( ElementSubRegionBase & subRegion )
const
225 void updatePorosityAndPermeability( CellElementSubRegion & subRegion )
const
227 void updateSolidInternalEnergyModel( ObjectManagerBase & dataGroup )
const
232 void enableJumpStabilization()
235 void enableFixedStressPoromechanicsUpdate()
238 void setKeepVariablesConstantDuringInitStep(
bool const keepVariablesConstantDuringInitStep )
240 reservoirSolver()->setKeepVariablesConstantDuringInitStep( keepVariablesConstantDuringInitStep );
241 wellSolver()->setKeepVariablesConstantDuringInitStep( keepVariablesConstantDuringInitStep );
260 coupledReservoirAndWellsInternal::
261 addCouplingNumNonzeros(
this,
282 virtual void setMGRStrategy()
309 this->
template forDiscretizationOnMeshTargets<>( domain.
getMeshBodies(), [&] (
string const &,
313 ElementRegionManager & elemManager = meshLevel.getElemManager();
315 ElementRegionManager::ElementViewAccessor< arrayView2d< real64 > > const elemCenter =
316 elemManager.constructViewAccessor< array2d< real64 >, arrayView2d< real64 > >( ElementSubRegionBase::viewKeyStruct::elementCenterString() );
319 elemManager.forElementSubRegions< WellElementSubRegion >( regionNames, [&]( localIndex const,
320 WellElementSubRegion & subRegion )
322 array1d< array1d< arrayView3d< real64 const > > > const permeability =
323 elemManager.constructMaterialFieldAccessor< constitutive::PermeabilityBase,
324 fields::permeability::permeability >();
326 PerforationData & perforationData = *subRegion.getPerforationData();
327 WellControls const & wellControls = wellSolver()->getWellControls( subRegion );
330 perforationData.computeWellTransmissibility( meshLevel, subRegion, permeability );
333 if( wellControls.getLogLevel() >= 2 )
335 arrayView2d< real64 const > const perfLocation =
336 perforationData.getField< fields::perforation::location >();
337 arrayView1d< real64 const > const perfTrans =
338 perforationData.getField< fields::perforation::wellTransmissibility >();
341 arrayView1d< localIndex const > const resElemRegion =
342 perforationData.getField< fields::perforation::reservoirElementRegion >();
343 arrayView1d< localIndex const > const resElemSubRegion =
344 perforationData.getField< fields::perforation::reservoirElementSubRegion >();
345 arrayView1d< localIndex const > const resElemIndex =
346 perforationData.getField< fields::perforation::reservoirElementIndex >();
348 GEOS_UNUSED_VAR( perfLocation );
349 GEOS_UNUSED_VAR( perfTrans );
350 GEOS_UNUSED_VAR( resElemRegion );
351 GEOS_UNUSED_VAR( resElemSubRegion );
352 GEOS_UNUSED_VAR( resElemIndex );
354 forAll< serialPolicy >( perforationData.size(), [=] ( localIndex const iperf )
356 GEOS_UNUSED_VAR( iperf );
357 GEOS_LOG_RANK( GEOS_FMT(
"{}: perforation at ({},{},{}), perforated element center = ({},{},{}), transmissibility = {} [{}]",
358 this->getName(), perfLocation[iperf][0], perfLocation[iperf][1], perfLocation[iperf][2],
359 elemCenter[resElemRegion[iperf]][resElemSubRegion[iperf]][resElemIndex[iperf]][0],
360 elemCenter[resElemRegion[iperf]][resElemSubRegion[iperf]][resElemIndex[iperf]][1],
361 elemCenter[resElemRegion[iperf]][resElemSubRegion[iperf]][resElemIndex[iperf]][2],
362 perfTrans[iperf], getSymbol( units::Transmissibility ) ) );
bool validateWellPerforations(PhysicsSolverBase const *const reservoirSolver, WellSolverBase const *const wellSolver, DomainPartition const &domain)
Validate the well perforations ensuring that each perforation is located in a reservoir region that i...
void addCouplingNumNonzeros(PhysicsSolverBase const *const solver, DomainPartition &domain, DofManager &dofManager, arrayView1d< localIndex > const &rowLengths, integer const resNumDof, integer const wellNumDof, string const &resElemDofName, string const &wellElemDofName)
Utility function for the implementation details of addCouplingNumZeros.
#define GEOS_ERROR(msg)
Raise a hard error and terminate the program.
bool m_isWellTransmissibilityComputed
Flag to determine whether the well transmissibility needs to be computed.
virtual void initializePostInitialConditionsPreSubGroups() override
Called by InitializePostInitialConditions() prior to initializing sub-Groups.
virtual void addCouplingSparsityPattern(DomainPartition const &domain, DofManager const &dofManager, SparsityPatternView< globalIndex > const &pattern) const =0
static string coupledSolverAttributePrefix()
String used to form the solverName used to register solvers in CoupledSolver.
CoupledReservoirAndWellsBase(const string &name, dataRepository::Group *const parent)
main constructor for ManagedGroup Objects
virtual void postInputInitialization() override
WELL_SOLVER * wellSolver() const
accessor for the pointer to the well solver
void addCouplingNumNonzeros(DomainPartition &domain, DofManager &dofManager, arrayView1d< localIndex > const &rowLengths) const
virtual void saveSequentialIterationState(DomainPartition &domain) override
Save the state of the solver for sequential iteration.
virtual void implicitStepSetup(real64 const &time_n, real64 const &dt, DomainPartition &domain) override
function to perform setup for implicit timestep
RESERVOIR_SOLVER * reservoirSolver() const
accessor for the pointer to the reservoir solver
virtual void postInputInitialization() override
std::array< string, sizeof...(SOLVERS) > m_names
Names of the single-physics solvers.
std::tuple< SOLVERS *... > m_solvers
Pointers of the single-physics solvers.
The DoFManager is responsible for allocating global dofs, constructing sparsity patterns,...
Partition of the decomposed physical domain. It also manages the connexion information to its neighbo...
Group const & getMeshBodies() const
Get the mesh bodies, const version.
Class facilitating the representation of a multi-level discretization of a MeshBody.
Base class for all physics solvers.
virtual string getCatalogName() const =0
CRSMatrix< real64, globalIndex > m_localMatrix
Local system matrix and rhs.
DofManager m_dofManager
Data structure to handle degrees of freedom.
ParallelVector m_solution
System solution vector.
ParallelVector m_rhs
System right-hand side vector.
LinearSolverParametersInput m_linearSolverParameters
Linear solver parameters.
string const & getName() const
Get group name.
virtual void initializePostInitialConditionsPreSubGroups()
Called by InitializePostInitialConditions() prior to initializing sub-Groups.
virtual void implicitStepSetup(real64 const &time_n, real64 const &dt, DomainPartition &domain) override
virtual void setSparsityPattern(DomainPartition &domain, DofManager &dofManager, CRSMatrix< real64, globalIndex > &localMatrix, SparsityPattern< globalIndex > &pattern) override
Set the sparsity pattern of the linear system matrix.
@ FALSE
Not read from input.
ArrayView< T, 1 > arrayView1d
Alias for 1D array view.
stdVector< string > string_array
A 1-dimensional array of geos::string types.
GEOS_GLOBALINDEX_TYPE globalIndex
Global index type (for indexing objects across MPI partitions).
LvArray::CRSMatrix< T, COL_INDEX, INDEX_TYPE, LvArray::ChaiBuffer > CRSMatrix
Alias for CRS Matrix class.
LvArray::SparsityPatternView< COL_INDEX, INDEX_TYPE const, LvArray::ChaiBuffer > SparsityPatternView
Alias for Sparsity pattern View.
LvArray::SparsityPattern< COL_INDEX, INDEX_TYPE, LvArray::ChaiBuffer > SparsityPattern
Alias for Sparsity pattern class.
double real64
64-bit floating point type.
GEOS_LOCALINDEX_TYPE localIndex
Local index type (for indexing objects within an MPI partition).
int integer
Signed integer type.
Array< T, 1 > array1d
Alias for 1D array.
@ mgr
Multigrid reduction (Hypre only)
PreconditionerType preconditionerType
Preconditioner type.