diff --git a/.integrated_tests.yaml b/.integrated_tests.yaml
index 6009de53427..746f223b649 100644
--- a/.integrated_tests.yaml
+++ b/.integrated_tests.yaml
@@ -1,6 +1,6 @@
baselines:
bucket: geosx
- baseline: integratedTests/baseline_integratedTests-pr3384-9542-4e0f693
+ baseline: integratedTests/baseline_integratedTests-pr3495-9552-257c87e
allow_fail:
all: ''
streak: ''
diff --git a/BASELINE_NOTES.md b/BASELINE_NOTES.md
index 26bf74a58a5..c5617ba8ed9 100644
--- a/BASELINE_NOTES.md
+++ b/BASELINE_NOTES.md
@@ -6,6 +6,10 @@ This file is designed to track changes to the integrated test baselines.
Any developer who updates the baseline ID in the .integrated_tests.yaml file is expected to create an entry in this file with the pull request number, date, and their justification for rebaselining.
These notes should be in reverse-chronological order, and use the following time format: (YYYY-MM-DD).
+PR #3495 (2024-01-08)
+=====================
+Add missing logic to support switching from fixed mass rate injection rate constraint to max injection pressure.
+
PR #3384 (2024-01-07)
=====================
Added plastic strain output.
diff --git a/inputFiles/compositionalMultiphaseWell/compositionalMultiphaseWell.ats b/inputFiles/compositionalMultiphaseWell/compositionalMultiphaseWell.ats
index ea1d810c064..16de5c4d90c 100644
--- a/inputFiles/compositionalMultiphaseWell/compositionalMultiphaseWell.ats
+++ b/inputFiles/compositionalMultiphaseWell/compositionalMultiphaseWell.ats
@@ -76,6 +76,22 @@ decks = [
partitions=[(1, 1, 1), (2, 2, 1)],
restart_step=12,
check_step=25,
+ restartcheck_params=RestartcheckParameters(**restartcheck_params)),
+ TestDeck(
+ name="isothm_mass_inj_table",
+ description=
+ "Isothermal mass injection constraint testing control switching",
+ partitions=[(1, 1, 1), (1, 2, 1)],
+ restart_step=24,
+ check_step=28,
+ restartcheck_params=RestartcheckParameters(**restartcheck_params)),
+ TestDeck(
+ name="isothm_vol_inj_table",
+ description=
+ "Isothermal vol injection constraint testing control switching, should have same results as isothm_mass_inj_table.xml",
+ partitions=[(1, 1, 1), (1, 2, 1)],
+ restart_step=24,
+ check_step=28,
restartcheck_params=RestartcheckParameters(**restartcheck_params))
]
diff --git a/inputFiles/compositionalMultiphaseWell/isothm_mass_inj_table.xml b/inputFiles/compositionalMultiphaseWell/isothm_mass_inj_table.xml
new file mode 100644
index 00000000000..14076fab917
--- /dev/null
+++ b/inputFiles/compositionalMultiphaseWell/isothm_mass_inj_table.xml
@@ -0,0 +1,229 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/inputFiles/compositionalMultiphaseWell/isothm_vol_inj_table.xml b/inputFiles/compositionalMultiphaseWell/isothm_vol_inj_table.xml
new file mode 100644
index 00000000000..d39421e36fc
--- /dev/null
+++ b/inputFiles/compositionalMultiphaseWell/isothm_vol_inj_table.xml
@@ -0,0 +1,234 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/wells/CompositionalMultiphaseWell.cpp b/src/coreComponents/physicsSolvers/fluidFlow/wells/CompositionalMultiphaseWell.cpp
index 90f3a8677ec..b891467863d 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/wells/CompositionalMultiphaseWell.cpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/wells/CompositionalMultiphaseWell.cpp
@@ -260,6 +260,8 @@ void CompositionalMultiphaseWell::registerDataOnMesh( Group & meshBodies )
wellControls.registerWrapper< real64 >( viewKeyStruct::massDensityString() );
+ wellControls.registerWrapper< real64 >( viewKeyStruct::currentMassRateString() );
+
// write rates output header
// the rank that owns the reference well element is responsible
if( m_writeCSV > 0 && subRegion.isLocallyOwned() )
@@ -452,6 +454,7 @@ void CompositionalMultiphaseWell::validateWellConstraints( real64 const & time_n
WellControls::Control const currentControl = wellControls.getControl();
real64 const & targetTotalRate = wellControls.getTargetTotalRate( time_n + dt );
real64 const & targetPhaseRate = wellControls.getTargetPhaseRate( time_n + dt );
+ real64 const & targetMassRate = wellControls.getTargetMassRate( time_n + dt );
GEOS_THROW_IF( wellControls.isInjector() && currentControl == WellControls::Control::PHASEVOLRATE,
"WellControls " << wellControls.getDataContext() <<
@@ -470,6 +473,18 @@ void CompositionalMultiphaseWell::validateWellConstraints( real64 const & time_n
"WellControls " << wellControls.getDataContext() <<
": Target phase rate cannot be used for injectors",
InputError );
+ GEOS_THROW_IF( wellControls.isProducer() && !isZero( targetTotalRate ),
+ "WellControls " << wellControls.getDataContext() <<
+ ": Target total rate cannot be used for producers",
+ InputError );
+ GEOS_THROW_IF( wellControls.isProducer() && !isZero( targetMassRate ),
+ "WellControls " << wellControls.getDataContext() <<
+ ": Target mass rate cannot be used for producers",
+ InputError );
+ GEOS_THROW_IF( !m_useMass && !isZero( targetMassRate ),
+ "WellControls " << wellControls.getDataContext() <<
+ ": Target mass rate cannot with useMass=0",
+ InputError );
// The user always provides positive rates, but these rates are later multiplied by -1 internally for producers
GEOS_THROW_IF( wellControls.isProducer() && targetPhaseRate > 0.0,
@@ -686,6 +701,10 @@ void CompositionalMultiphaseWell::updateVolRatesForConstraint( WellElementSubReg
real64 & currentTotalVolRate =
wellControls.getReference< real64 >( CompositionalMultiphaseWell::viewKeyStruct::currentTotalVolRateString() );
+
+ real64 & currentMassRate =
+ wellControls.getReference< real64 >( CompositionalMultiphaseWell::viewKeyStruct::currentMassRateString() );
+
arrayView1d< real64 > const & dCurrentTotalVolRate =
wellControls.getReference< array1d< real64 > >( CompositionalMultiphaseWell::viewKeyStruct::dCurrentTotalVolRateString() );
@@ -721,6 +740,7 @@ void CompositionalMultiphaseWell::updateVolRatesForConstraint( WellElementSubReg
dCurrentTotalVolRate,
currentPhaseVolRate,
dCurrentPhaseVolRate,
+ ¤tMassRate,
&iwelemRef,
&logLevel,
&wellControlsName,
@@ -757,7 +777,8 @@ void CompositionalMultiphaseWell::updateVolRatesForConstraint( WellElementSubReg
// Step 2: update the total volume rate
real64 const currentTotalRate = connRate[iwelemRef];
-
+ // Assumes useMass is true
+ currentMassRate = currentTotalRate;
// Step 2.1: compute the inverse of the total density and derivatives
massDensity = totalDens[iwelemRef][0];
real64 const totalDensInv = 1.0 / totalDens[iwelemRef][0];
@@ -1959,6 +1980,12 @@ void CompositionalMultiphaseWell::assemblePressureRelations( real64 const & time
GEOS_LOG_LEVEL_INFO_RANK_0( logInfo::WellControl,
GEOS_FMT( "Control switch for well {} from BHP constraint to phase volumetric rate constraint", subRegion.getName() ) );
}
+ else if( wellControls.getInputControl() == WellControls::Control::MASSRATE )
+ {
+ wellControls.switchToMassRateControl( wellControls.getTargetMassRate( timeAtEndOfStep ) );
+ GEOS_LOG_LEVEL_INFO_RANK_0( logInfo::WellControl,
+ GEOS_FMT( "Control switch for well {} from BHP constraint to mass rate constraint", subRegion.getName()) );
+ }
else
{
wellControls.switchToTotalRateControl( wellControls.getTargetTotalRate( timeAtEndOfStep ) );
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/wells/CompositionalMultiphaseWell.hpp b/src/coreComponents/physicsSolvers/fluidFlow/wells/CompositionalMultiphaseWell.hpp
index aa70d4f4507..72614596376 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/wells/CompositionalMultiphaseWell.hpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/wells/CompositionalMultiphaseWell.hpp
@@ -296,6 +296,8 @@ class CompositionalMultiphaseWell : public WellSolverBase
static constexpr char const * currentTotalVolRateString() { return "currentTotalVolumetricRate"; }
static constexpr char const * dCurrentTotalVolRateString() { return "dCurrentTotalVolumetricRate"; }
+ static constexpr char const * currentMassRateString() { return "currentMassRate"; }
+
static constexpr char const * dCurrentTotalVolRate_dPresString() { return "dCurrentTotalVolumetricRate_dPres"; }
static constexpr char const * dCurrentTotalVolRate_dCompDensString() { return "dCurrentTotalVolumetricRate_dCompDens"; }
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.cpp b/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.cpp
index 7ce65b6cbec..47877476c2f 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.cpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.cpp
@@ -71,27 +71,32 @@ WellControls::WellControls( string const & name, Group * const parent )
registerWrapper( viewKeyStruct::targetBHPString(), &m_targetBHP ).
setDefaultValue( 0.0 ).
setInputFlag( InputFlags::OPTIONAL ).
+ setRestartFlags( RestartFlags::WRITE_AND_READ ).
setDescription( "Target bottom-hole pressure [Pa]" );
registerWrapper( viewKeyStruct::targetTotalRateString(), &m_targetTotalRate ).
setDefaultValue( 0.0 ).
setInputFlag( InputFlags::OPTIONAL ).
+ setRestartFlags( RestartFlags::WRITE_AND_READ ).
setDescription( "Target total volumetric rate (if useSurfaceConditions: [surface m^3/s]; else [reservoir m^3/s])" );
registerWrapper( viewKeyStruct::targetPhaseRateString(), &m_targetPhaseRate ).
setDefaultValue( 0.0 ).
setInputFlag( InputFlags::OPTIONAL ).
+ setRestartFlags( RestartFlags::WRITE_AND_READ ).
setDescription( "Target phase volumetric rate (if useSurfaceConditions: [surface m^3/s]; else [reservoir m^3/s])" );
registerWrapper( viewKeyStruct::targetMassRateString(), &m_targetMassRate ).
setDefaultValue( 0.0 ).
setInputFlag( InputFlags::OPTIONAL ).
+ setRestartFlags( RestartFlags::WRITE_AND_READ ).
setDescription( "Target Mass Rate rate ( [kg^3/s])" );
registerWrapper( viewKeyStruct::targetPhaseNameString(), &m_targetPhaseName ).
setRTTypeName( rtTypes::CustomTypes::groupNameRef ).
setDefaultValue( "" ).
setInputFlag( InputFlags::OPTIONAL ).
+ setRestartFlags( RestartFlags::WRITE_AND_READ ).
setDescription( "Name of the target phase" );
registerWrapper( viewKeyStruct::refElevString(), &m_refElevation ).
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.hpp b/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.hpp
index 3335b7d98dc..9e202096385 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.hpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.hpp
@@ -149,6 +149,12 @@ class WellControls : public dataRepository::Group
*/
Control getControl() const { return m_currentControl; }
+ /**
+ * @brief Get the input control type for the well.
+ * @return the Control enum enforced at the well
+ */
+ Control getInputControl() const { return m_inputControl; }
+
/**
* @brief Getter for the reference elevation where the BHP control is enforced
* @return the reference elevation
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/wells/kernels/CompositionalMultiphaseWellKernels.cpp b/src/coreComponents/physicsSolvers/fluidFlow/wells/kernels/CompositionalMultiphaseWellKernels.cpp
index 42ec6911080..01fb94d2e72 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/wells/kernels/CompositionalMultiphaseWellKernels.cpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/wells/kernels/CompositionalMultiphaseWellKernels.cpp
@@ -37,6 +37,7 @@ inline
void
ControlEquationHelper::
switchControl( bool const isProducer,
+ WellControls::Control const & inputControl,
WellControls::Control const & currentControl,
integer const phasePhaseIndex,
real64 const & targetBHP,
@@ -46,6 +47,7 @@ ControlEquationHelper::
real64 const & currentBHP,
arrayView1d< real64 const > const & currentPhaseVolRate,
real64 const & currentTotalVolRate,
+ real64 const & currentMassRate,
WellControls::Control & newControl )
{
// if isViable is true at the end of the following checks, no need to switch
@@ -58,7 +60,7 @@ ControlEquationHelper::
// Currently, the available constraints are:
// - Producer: BHP, PHASEVOLRATE
- // - Injector: BHP, TOTALVOLRATE
+ // - Injector: BHP, TOTALVOLRATE, MASSRATE
// TODO: support GRAT, WRAT, LIQUID for producers and check if any of the active constraint is violated
@@ -71,6 +73,10 @@ ControlEquationHelper::
controlIsViable = ( LvArray::math::abs( currentPhaseVolRate[phasePhaseIndex] ) <= LvArray::math::abs( targetPhaseRate ) );
}
// the control is viable if the reference total rate is below the max rate for injectors
+ else if( inputControl == WellControls::Control::MASSRATE )
+ {
+ controlIsViable = ( LvArray::math::abs( currentMassRate ) <= LvArray::math::abs( targetMassRate ) );
+ }
else
{
controlIsViable = ( LvArray::math::abs( currentTotalVolRate ) <= LvArray::math::abs( targetTotalRate ) );
@@ -219,17 +225,6 @@ ControlEquationHelper::
if constexpr ( IS_THERMAL )
dControlEqn[COFFSET_WJ::dT] = massDensity*dCurrentTotalVolRate[COFFSET_WJ::dT];
}
- // Total mass rate control
- else if( currentControl == WellControls::Control::MASSRATE )
- {
- controlEqn = massDensity*currentTotalVolRate - targetMassRate;
- dControlEqn[COFFSET_WJ::dP] = massDensity*dCurrentTotalVolRate[COFFSET_WJ::dP];
- dControlEqn[COFFSET_WJ::dQ] = massDensity*dCurrentTotalVolRate[COFFSET_WJ::dQ];
- for( integer ic = 0; ic < NC; ++ic )
- {
- dControlEqn[COFFSET_WJ::dC+ic] = massDensity*dCurrentTotalVolRate[COFFSET_WJ::dC+ic];
- }
- }
else
{
GEOS_ERROR( "This constraint is not supported in CompositionalMultiphaseWell" );
@@ -323,6 +318,7 @@ PressureRelationKernel::
// static well control data
bool const isProducer = wellControls.isProducer();
WellControls::Control const currentControl = wellControls.getControl();
+ WellControls::Control const inputControl = wellControls.getInputControl();
real64 const targetBHP = wellControls.getTargetBHP( timeAtEndOfStep );
real64 const targetTotalRate = wellControls.getTargetTotalRate( timeAtEndOfStep );
real64 const targetPhaseRate = wellControls.getTargetPhaseRate( timeAtEndOfStep );
@@ -344,6 +340,9 @@ PressureRelationKernel::
arrayView1d< real64 const > const & dCurrentTotalVolRate =
wellControls.getReference< array1d< real64 > >( CompositionalMultiphaseWell::viewKeyStruct::dCurrentTotalVolRateString() );
+ real64 const & currentMassRate =
+ wellControls.getReference< real64 >( CompositionalMultiphaseWell::viewKeyStruct::currentMassRateString() );
+
real64 const & massDensity =
wellControls.getReference< real64 >( CompositionalMultiphaseWell::viewKeyStruct::massDensityString() );
@@ -358,6 +357,7 @@ PressureRelationKernel::
{
WellControls::Control newControl = currentControl;
ControlEquationHelper::switchControl( isProducer,
+ inputControl,
currentControl,
targetPhaseIndex,
targetBHP,
@@ -367,6 +367,7 @@ PressureRelationKernel::
currentBHP,
currentPhaseVolRate,
currentTotalVolRate,
+ currentMassRate,
newControl );
if( currentControl != newControl )
{
@@ -737,7 +738,6 @@ RateInitializationKernel::
else if( control == WellControls::Control::MASSRATE )
{
connRate[iwelem] = targetMassRate;
- connRate[iwelem] = targetMassRate* totalDens[iwelem][0];
}
else
{
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/wells/kernels/CompositionalMultiphaseWellKernels.hpp b/src/coreComponents/physicsSolvers/fluidFlow/wells/kernels/CompositionalMultiphaseWellKernels.hpp
index aae86fcf48c..c2e9078ba47 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/wells/kernels/CompositionalMultiphaseWellKernels.hpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/wells/kernels/CompositionalMultiphaseWellKernels.hpp
@@ -134,6 +134,7 @@ struct ControlEquationHelper
static
void
switchControl( bool const isProducer,
+ WellControls::Control const & inputControl,
WellControls::Control const & currentControl,
integer const phasePhaseIndex,
real64 const & targetBHP,
@@ -143,6 +144,7 @@ struct ControlEquationHelper
real64 const & currentBHP,
arrayView1d< real64 const > const & currentPhaseVolRate,
real64 const & currentTotalVolRate,
+ real64 const & currentMassRate,
WellControls::Control & newControl );
template< integer NC, integer IS_THERMAL >