GEOS
CompositionalMultiphaseBase.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_FLUIDFLOW_COMPOSITIONALMULTIPHASEBASE_HPP_
21 #define GEOS_PHYSICSSOLVERS_FLUIDFLOW_COMPOSITIONALMULTIPHASEBASE_HPP_
22 
24 #include "common/DataLayouts.hpp"
25 #include "constitutive/fluid/multifluid/Layouts.hpp"
26 #include "constitutive/relativePermeability/Layouts.hpp"
27 #include "constitutive/capillaryPressure/Layouts.hpp"
29 #include "constitutive/fluid/multifluid/MultiFluidBase.hpp"
30 #include "constitutive/solid/CoupledSolidBase.hpp"
31 #include "physicsSolvers/fluidFlow/kernels/compositional/AccumulationKernel.hpp"
34 
35 
36 namespace geos
37 {
38 
43 {
46 };
47 
52  "ComponentDensities",
53  "OverallComposition" );
54 
55 //START_SPHINX_INCLUDE_00
62 {
63 public:
64 
70  CompositionalMultiphaseBase( const string & name,
71  Group * const parent );
72 
75 
78 
81 
84 
87 
91  virtual ~CompositionalMultiphaseBase() override = default;
92 
93 //START_SPHINX_INCLUDE_01
94 
95  virtual void registerDataOnMesh( Group & meshBodies ) override;
96 
97  virtual void registerDataForCFL( Group & meshBodies ) { GEOS_UNUSED_VAR( meshBodies ); }
98 
106  virtual void
107  implicitStepSetup( real64 const & time_n,
108  real64 const & dt,
109  DomainPartition & domain ) override;
110 
111  virtual void
112  assembleSystem( real64 const time_n,
113  real64 const dt,
114  DomainPartition & domain,
115  DofManager const & dofManager,
116  CRSMatrixView< real64, globalIndex const > const & localMatrix,
117  arrayView1d< real64 > const & localRhs ) override;
118 
119  virtual void
121  real64 const dt,
122  DomainPartition & domain,
123  DofManager const & dofManager,
124  CRSMatrixView< real64, globalIndex const > const & localMatrix,
125  arrayView1d< real64 > const & localRhs ) override;
126 
136  template< typename SUBREGION_TYPE >
137  void accumulationAssemblyLaunch( DofManager const & dofManager,
138  SUBREGION_TYPE const & subRegion,
139  CRSMatrixView< real64, globalIndex const > const & localMatrix,
140  arrayView1d< real64 > const & localRhs ) const;
141 
142  virtual void
144 
145  virtual void
147  real64 const & dt,
148  DomainPartition & domain ) override;
149 
155 
161 
166  void updateFluidModel( ObjectManagerBase & dataGroup ) const;
167 
172  void updateRelPermModel( ObjectManagerBase & dataGroup ) const;
173 
178  void updateCapPressureModel( ObjectManagerBase & dataGroup ) const;
179 
184  void updateCompAmount( ElementSubRegionBase & subRegion ) const;
185 
190  void updateEnergy( ElementSubRegionBase & subRegion ) const;
191 
197 
202  virtual void updatePhaseMobility( ObjectManagerBase & dataGroup ) const = 0;
203 
204  real64 updateFluidState( ElementSubRegionBase & subRegion ) const;
205 
206  virtual void saveConvergedState( ElementSubRegionBase & subRegion ) const override final;
207 
208  virtual void saveSequentialIterationState( DomainPartition & domain ) override final;
209 
210  virtual void updateState( DomainPartition & domain ) override final;
211 
217 
222  integer numFluidPhases() const { return m_numPhases; }
223 
229 
233  virtual units::Unit getMassUnit() const override
234  { return m_useMass ? units::Unit::Mass : units::Unit::Mole; }
235 
246  DofManager const & dofManager,
247  CRSMatrixView< real64, globalIndex const > const & localMatrix,
248  arrayView1d< real64 > const & localRhs ) const;
249 
259  virtual void
261  DomainPartition const & domain,
262  DofManager const & dofManager,
263  CRSMatrixView< real64, globalIndex const > const & localMatrix,
264  arrayView1d< real64 > const & localRhs ) const = 0;
265 
274  virtual void
276  DomainPartition const & domain,
277  DofManager const & dofManager,
278  CRSMatrixView< real64, globalIndex const > const & localMatrix,
279  arrayView1d< real64 > const & localRhs ) const = 0;
280 
285  {
286  static constexpr char const * elemDofFieldString() { return "compositionalVariables"; }
287 
288  // inputs
289 
290  static constexpr char const * useMassFlagString() { return "useMass"; }
291  static constexpr char const * formulationTypeString() { return "formulationType"; }
292 
293  // time stepping controls
294 
295  static constexpr char const * solutionChangeScalingFactorString() { return "solutionChangeScalingFactor"; }
296  static constexpr char const * targetRelativePresChangeString() { return "targetRelativePressureChangeInTimeStep"; }
297  static constexpr char const * targetRelativeTempChangeString() { return "targetRelativeTemperatureChangeInTimeStep"; }
298  static constexpr char const * targetPhaseVolFracChangeString() { return "targetPhaseVolFractionChangeInTimeStep"; }
299  static constexpr char const * targetRelativeCompDensChangeString() { return "targetRelativeCompDensChangeInTimeStep"; }
300  static constexpr char const * targetCompFracChangeString() { return "targetCompFracChangeInTimeStep"; }
301  static constexpr char const * targetFlowCFLString() { return "targetFlowCFL"; }
302 
303 
304  // nonlinear solver parameters
305 
306  static constexpr char const * maxCompFracChangeString() { return "maxCompFractionChange"; }
307  static constexpr char const * maxRelativePresChangeString() { return "maxRelativePressureChange"; }
308  static constexpr char const * maxRelativeTempChangeString() { return "maxRelativeTemperatureChange"; }
309  static constexpr char const * maxRelativeCompDensChangeString() { return "maxRelativeCompDensChange"; }
310  static constexpr char const * allowLocalCompDensChoppingString() { return "allowLocalCompDensityChopping"; }
311  static constexpr char const * useTotalMassEquationString() { return "useTotalMassEquation"; }
312  static constexpr char const * useSimpleAccumulationString() { return "useSimpleAccumulation"; }
313  static constexpr char const * minCompDensString() { return "minCompDens"; }
314  static constexpr char const * minCompFracString() { return "minCompFrac"; }
315  static constexpr char const * maxSequentialCompDensChangeString() { return "maxSequentialCompDensChange"; }
316  static constexpr char const * minScalingFactorString() { return "minScalingFactor"; }
317 
318  static constexpr char const * relPermNamesString() { return "relPermNames"; }
319  static constexpr char const * capPressureNamesString() { return "capPressureNames"; }
320  static constexpr char const * diffusionNamesString() { return "diffusionNames"; }
321  static constexpr char const * dispersionNamesString() { return "dispersionNames"; }
322  };
323 
332  virtual void initializeFluidState( MeshLevel & mesh, string_array const & regionNames ) override;
333 
334  virtual void initializeThermalState( MeshLevel & mesh, string_array const & regionNames ) override;
335 
339  virtual void computeHydrostaticEquilibrium( DomainPartition & domain ) override;
340 
350  void applyDirichletBC( real64 const time,
351  real64 const dt,
352  DofManager const & dofManager,
353  DomainPartition & domain,
354  CRSMatrixView< real64, globalIndex const > const & localMatrix,
355  arrayView1d< real64 > const & localRhs ) const;
356 
366  void applySourceFluxBC( real64 const time,
367  real64 const dt,
368  DofManager const & dofManager,
369  DomainPartition & domain,
370  CRSMatrixView< real64, globalIndex const > const & localMatrix,
371  arrayView1d< real64 > const & localRhs ) const;
372 
382  virtual void applyAquiferBC( real64 const time,
383  real64 const dt,
384  DofManager const & dofManager,
385  DomainPartition & domain,
386  CRSMatrixView< real64, globalIndex const > const & localMatrix,
387  arrayView1d< real64 > const & localRhs ) const = 0;
388 
401  real64 const dt,
402  DofManager const & dofManager,
403  DomainPartition & domain,
404  CRSMatrixView< real64, globalIndex const > const & localMatrix,
405  arrayView1d< real64 > const & localRhs ) const;
406 
407 
413  void chopNegativeDensities( ElementSubRegionBase & subRegion );
414 
420 
421  virtual real64 setNextDtBasedOnStateChange( real64 const & currentDt,
422  DomainPartition & domain ) override;
423 
424 
425  virtual void computeCFLNumbers( DomainPartition & domain, real64 const & dt, real64 & maxPhaseCFL, real64 & maxCompCFL )
426  {
427  GEOS_UNUSED_VAR( domain, dt, maxPhaseCFL, maxCompCFL );
428  GEOS_ERROR( GEOS_FMT( "{}: computeCFLNumbers is not implemented for {}", getDataContext(), getCatalogName() ) );
429  }
430 
432 
433  integer useSimpleAccumulation() const { return m_useSimpleAccumulation; }
434 
435  integer useTotalMassEquation() const { return m_useTotalMassEquation; }
436 
437  virtual bool checkSequentialSolutionIncrements( DomainPartition & domain ) const override;
438 
439 protected:
440 
441  virtual void postInputInitialization() override;
442 
443  virtual void initializePreSubGroups() override;
444 
445 
452  void validateConstitutiveModels( DomainPartition const & domain ) const;
453 
458  void initializeAquiferBC( constitutive::ConstitutiveManager const & cm ) const;
459 
469  template< typename OBJECT_TYPE >
470  void applyFieldValue( real64 const & time_n,
471  real64 const & dt,
472  MeshLevel & mesh,
473  char const logMessage[],
474  string const fieldKey,
475  string const boundaryFieldKey ) const;
476 
479 
482 
485 
488 
491 
494 
497 
500 
503 
506 
509 
512 
515 
518 
521 
524 
527 
530 
533 
536 
539 
542 
545 
548 
551  real64 m_maxSequentialCompDensChange;
552 
553 private:
554 
560  bool validateDirichletBC( DomainPartition & domain,
561  real64 const time ) const;
562 
563  virtual void setConstitutiveNames( ElementSubRegionBase & subRegion ) const override;
564 
565 };
566 
567 template< typename OBJECT_TYPE >
569  real64 const & dt,
570  MeshLevel & mesh,
571  char const logMessage[],
572  string const fieldKey,
573  string const boundaryFieldKey ) const
574 {
576 
577  fsManager.apply< OBJECT_TYPE >( time_n + dt,
578  mesh,
579  fieldKey,
580  [&]( FieldSpecificationBase const & fs,
581  string const & setName,
583  OBJECT_TYPE & targetGroup,
584  string const & )
585  {
586  if( fs.getLogLevel() >= 1 && m_nonlinearSolverParameters.m_numNewtonIterations == 0 )
587  {
588  globalIndex const numTargetElems = MpiWrapper::sum< globalIndex >( lset.size() );
589  GEOS_LOG_RANK_0( GEOS_FMT( logMessage,
590  getName(), time_n+dt, fs.getCatalogName(), fs.getName(),
591  setName, targetGroup.getName(), fs.getScale(), numTargetElems ) );
592  }
593 
594  // Specify the bc value of the field
595  fs.applyFieldValue< FieldSpecificationEqual,
596  parallelDevicePolicy<> >( lset,
597  time_n + dt,
598  targetGroup,
599  boundaryFieldKey );
600  } );
601 }
602 
603 template< typename SUBREGION_TYPE >
605  SUBREGION_TYPE const & subRegion,
606  CRSMatrixView< real64, globalIndex const > const & localMatrix,
607  arrayView1d< real64 > const & localRhs ) const
608 {
609  constitutive::MultiFluidBase const & fluid =
610  getConstitutiveModel< constitutive::MultiFluidBase >( subRegion, subRegion.template getReference< string >( viewKeyStruct::fluidNamesString() ) );
611  constitutive::CoupledSolidBase const & solid =
612  getConstitutiveModel< constitutive::CoupledSolidBase >( subRegion, subRegion.template getReference< string >( viewKeyStruct::solidNamesString() ) );
613 
614  string const dofKey = dofManager.getKey( viewKeyStruct::elemDofFieldString() );
615 
616  using namespace isothermalCompositionalMultiphaseBaseKernels;
617 
618  BitFlags< KernelFlags > kernelFlags;
620  kernelFlags.set( KernelFlags::TotalMassEquation );
622  kernelFlags.set( KernelFlags::SimpleAccumulation );
623 
625  {
626  // isothermal for now
627  isothermalCompositionalMultiphaseBaseKernels::
628  AccumulationZFormulationKernelFactory::
629  createAndLaunch< parallelDevicePolicy<> >( m_numComponents,
630  m_numPhases,
631  dofManager.rankOffset(),
632  kernelFlags,
633  dofKey,
634  subRegion,
635  fluid,
636  solid,
637  localMatrix,
638  localRhs );
639  }
640  else
641  {
642  if( m_isThermal )
643  {
644  thermalCompositionalMultiphaseBaseKernels::
645  AccumulationKernelFactory::
646  createAndLaunch< parallelDevicePolicy<> >( m_numComponents,
647  m_numPhases,
648  dofManager.rankOffset(),
649  kernelFlags,
650  dofKey,
651  subRegion,
652  fluid,
653  solid,
654  localMatrix,
655  localRhs );
656  }
657  else
658  {
659  isothermalCompositionalMultiphaseBaseKernels::
660  AccumulationKernelFactory::
661  createAndLaunch< parallelDevicePolicy<> >( m_numComponents,
662  m_numPhases,
663  dofManager.rankOffset(),
664  kernelFlags,
665  dofKey,
666  subRegion,
667  fluid,
668  solid,
669  localMatrix,
670  localRhs );
671  }
672  }
673 }
674 
675 } // namespace geos
676 
677 #endif //GEOS_PHYSICSSOLVERS_FLUIDFLOW_COMPOSITIONALMULTIPHASEBASE_HPP_
#define GEOS_UNUSED_VAR(...)
Mark an unused variable and silence compiler warnings.
Definition: GeosxMacros.hpp:84
#define GEOS_ERROR(msg)
Raise a hard error and terminate the program.
Definition: Logger.hpp:157
#define GEOS_LOG_RANK_0(msg)
Log a message on screen on rank 0.
Definition: Logger.hpp:101
Unit
Enumerator of available unit types for given physical scales. Units are in SI by default.
Definition: Units.hpp:59
void initializeAquiferBC(constitutive::ConstitutiveManager const &cm) const
Initialize the aquifer boundary condition (gravity vector, water phase index)
bool m_hasCapPressure
flag to determine whether or not to apply capillary pressure
real64 m_solutionChangeScalingFactor
damping factor for solution change targets
CompositionalMultiphaseBase(CompositionalMultiphaseBase &&)=default
default move constructor
void applyFieldValue(real64 const &time_n, real64 const &dt, MeshLevel &mesh, char const logMessage[], string const fieldKey, string const boundaryFieldKey) const
Utility function that encapsulates the call to FieldSpecificationBase::applyFieldValue in BC applicat...
virtual void initializePreSubGroups() override
Called by Initialize() prior to initializing sub-Groups.
real64 m_targetPhaseVolFracChange
target (absolute) change in phase volume fraction in a time step
virtual void postInputInitialization() override
void keepVariablesConstantDuringInitStep(real64 const time, real64 const dt, DofManager const &dofManager, DomainPartition &domain, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs) const
Function to fix the initial state during the initialization step in coupled problems.
integer m_useMass
flag indicating whether mass or molar formulation should be used
integer m_useSimpleAccumulation
flag indicating whether simple accumulation form is used
real64 m_maxRelativeTempChange
maximum (relative) change in temperature in a Newton iteration
virtual bool checkSequentialSolutionIncrements(DomainPartition &domain) const override
Check if the solution increments are ok to use.
real64 m_targetCompFracChange
target (absolute) change in component fraction in a time step
virtual void initializeFluidState(MeshLevel &mesh, string_array const &regionNames) override
Initialize all variables from initial conditions.
real64 m_maxRelativePresChange
maximum (relative) change in pressure in a Newton iteration
real64 m_targetRelativePresChange
target (relative) change in pressure in a time step
real64 m_minCompFrac
minimum allowed global component fraction
void validateConstitutiveModels(DomainPartition const &domain) const
Utility function that checks the consistency of the constitutive models.
real64 m_minScalingFactor
minimum value of the scaling factor obtained by enforcing maxCompFracChange
virtual void applyAquiferBC(real64 const time, real64 const dt, DofManager const &dofManager, DomainPartition &domain, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs) const =0
Apply aquifer boundary conditions to the system.
bool m_hasDiffusion
flag to determine whether or not to apply diffusion
virtual void registerDataOnMesh(Group &meshBodies) override
Register wrappers that contain data on the mesh objects.
void applySourceFluxBC(real64 const time, real64 const dt, DofManager const &dofManager, DomainPartition &domain, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs) const
Apply source flux boundary conditions to the system.
virtual void initializePostInitialConditionsPreSubGroups() override
Called by InitializePostInitialConditions() prior to initializing sub-Groups.
real64 m_maxCompFracChange
maximum (absolute) change in a component fraction in a Newton iteration
real64 m_sequentialCompDensChange
maximum (absolute) component density change in a sequential iteration
real64 m_targetRelativeCompDensChange
target (relative) change in component density in a time step
virtual real64 setNextDtBasedOnStateChange(real64 const &currentDt, DomainPartition &domain) override
function to set the next dt based on state change
void chopNegativeCompFractions(DomainPartition &domain)
Sets all the negative component fractions (if any) to zero.
CompositionalMultiphaseBase()=delete
deleted default constructor
integer m_allowCompDensChopping
flag indicating whether local (cell-wise) chopping of negative compositions is allowed
string m_referenceFluidModelName
name of the fluid constitutive model used as a reference for component/phase description
integer m_numComponents
the number of fluid components
CompositionalMultiphaseBase(CompositionalMultiphaseBase const &)=delete
deleted copy constructor
integer m_useTotalMassEquation
flag indicating whether total mass equation is used
virtual ~CompositionalMultiphaseBase() override=default
default destructor
void chopNegativeDensities(DomainPartition &domain)
Sets all the negative component densities (if any) to zero.
real64 m_maxRelativeCompDensChange
maximum (relative) change in component density in a Newton iteration
virtual void computeHydrostaticEquilibrium(DomainPartition &domain) override
Compute the hydrostatic equilibrium using the compositions and temperature input tables.
real64 m_minCompDens
minimum allowed global component density
CompositionalMultiphaseFormulationType m_formulationType
formulation type
CompositionalMultiphaseBase(const string &name, Group *const parent)
main constructor for Group Objects
integer m_numPhases
the max number of fluid phases
bool m_hasDispersion
flag to determine whether or not to apply dispersion
real64 m_targetRelativeTempChange
target (relative) change in temperature in a time step
void applyDirichletBC(real64 const time, real64 const dt, DofManager const &dofManager, DomainPartition &domain, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs) const
Function to perform the Application of Dirichlet type BC's.
CompositionalMultiphaseBase & operator=(CompositionalMultiphaseBase const &)=delete
deleted assignment operator
CompositionalMultiphaseBase & operator=(CompositionalMultiphaseBase &&)=delete
deleted move operator
The DoFManager is responsible for allocating global dofs, constructing sparsity patterns,...
Definition: DofManager.hpp:45
globalIndex rankOffset(string const &fieldName) const
string const & getKey(string const &fieldName) const
Return the key used to record the field in the DofManager.
Partition of the decomposed physical domain. It also manages the connexion information to its neighbo...
void apply(real64 const time, MeshLevel &mesh, string const &fieldName, LAMBDA &&lambda) const
This function is the main driver for the field applications.
static FieldSpecificationManager & getInstance()
integer m_isThermal
flag to determine whether or not this is a thermal simulation
Class facilitating the representation of a multi-level discretization of a MeshBody.
Definition: MeshLevel.hpp:42
integer m_numNewtonIterations
The number of nonlinear iterations that have been exectued.
The ObjectManagerBase is the base object of all object managers in the mesh data hierachy.
virtual string getCatalogName() const =0
NonlinearSolverParameters m_nonlinearSolverParameters
Nonlinear solver parameters.
DataContext const & getDataContext() const
Definition: Group.hpp:1345
string const & getName() const
Get group name.
Definition: Group.hpp:1331
virtual void resetStateToBeginningOfStep(DomainPartition &domain) override
reset state of physics back to the beginning of the step.
virtual void implicitStepSetup(real64 const &time_n, real64 const &dt, DomainPartition &domain) override
function to perform setup for implicit timestep
real64 updatePhaseVolumeFraction(ObjectManagerBase &dataGroup) const
Recompute phase volume fractions (saturations) from constitutive and primary variables.
virtual units::Unit getMassUnit() const override
virtual void assembleSystem(real64 const time_n, 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 assembleStabilizedFluxTerms(real64 const dt, DomainPartition const &domain, DofManager const &dofManager, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs) const =0
assembles the flux terms for all cells with pressure jump stabilization
virtual void saveSequentialIterationState(DomainPartition &domain) override final
Utility function to save the iteration state (useful for sequential simulations)
virtual void implicitStepComplete(real64 const &time, real64 const &dt, DomainPartition &domain) override
perform cleanup for implicit timestep
integer numFluidComponents() const
Getter for the number of fluid components (species)
virtual void applyBoundaryConditions(real64 const time_n, real64 const dt, DomainPartition &domain, DofManager const &dofManager, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs) override
apply boundary condition to system
virtual void assembleFluxTerms(real64 const dt, DomainPartition const &domain, DofManager const &dofManager, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs) const =0
assembles the flux terms for all cells
void updateFluidModel(ObjectManagerBase &dataGroup) const
Update all relevant fluid models using current values of pressure and composition.
virtual void updatePhaseMobility(ObjectManagerBase &dataGroup) const =0
Recompute phase mobility from constitutive and primary variables.
void updateRelPermModel(ObjectManagerBase &dataGroup) const
Update all relevant relperm models using current values of phase volume fraction.
virtual void updateState(DomainPartition &domain) override final
Recompute all dependent quantities from primary variables (including constitutive models)
void assembleLocalTerms(DomainPartition &domain, DofManager const &dofManager, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs) const
assembles the accumulation other local terms for all cells
void updateCapPressureModel(ObjectManagerBase &dataGroup) const
Update all relevant capillary pressure models using current values of phase volume fraction.
void updateEnergy(ElementSubRegionBase &subRegion) const
Update energy.
void updateSolidInternalEnergyModel(ObjectManagerBase &dataGroup) const
Update all relevant solid internal energy models using current values of temperature.
string referenceFluidModelName() const
Getter for the name of the reference fluid model name.
virtual void saveConvergedState(ElementSubRegionBase &subRegion) const override final
Utility function to save the converged state.
void updateCompAmount(ElementSubRegionBase &subRegion) const
Update components mass/moles.
void accumulationAssemblyLaunch(DofManager const &dofManager, SUBREGION_TYPE const &subRegion, CRSMatrixView< real64, globalIndex const > const &localMatrix, arrayView1d< real64 > const &localRhs) const
assembles the accumulation terms for all cells of a spcefici subRegion.
integer numFluidPhases() const
Getter for the number of fluid phases.
void updateGlobalComponentFraction(ObjectManagerBase &dataGroup) const
Recompute global component fractions from primary variables (component densities)
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
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
LvArray::SortedArrayView< T, localIndex, LvArray::ChaiBuffer > SortedArrayView
A sorted array view of local indices.
Definition: DataTypes.hpp:270
CompositionalMultiphaseFormulationType
Options for flow formulation.
@ OverallComposition
use overall composition (z_c) as primary variables
@ ComponentDensities
use component densities as primary variables
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.