From 88665f12b19e62a47c1e06c29eacefd12f438c11 Mon Sep 17 00:00:00 2001 From: Tom Byer <149726499+tjb-ltk@users.noreply.github.com> Date: Fri, 20 Feb 2026 09:05:07 -0800 Subject: [PATCH 1/4] call flash with P&T --- .../fluidFlow/wells/SinglePhaseWell.cpp | 32 ++++++++++++++++--- 1 file changed, 27 insertions(+), 5 deletions(-) diff --git a/src/coreComponents/physicsSolvers/fluidFlow/wells/SinglePhaseWell.cpp b/src/coreComponents/physicsSolvers/fluidFlow/wells/SinglePhaseWell.cpp index 6414693bf22..d4999d44d2e 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/wells/SinglePhaseWell.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/wells/SinglePhaseWell.cpp @@ -266,6 +266,8 @@ void SinglePhaseWell::updateVolRateForConstraint( ElementRegionManager const & e arrayView1d< real64 const > const pres = subRegion.getField< well::pressure >(); + arrayView1d< real64 const > const & temp = subRegion.getField< well::temperature >(); + arrayView1d< real64 const > const & connRate = subRegion.getField< well::connectionRate >(); @@ -282,11 +284,12 @@ void SinglePhaseWell::updateVolRateForConstraint( ElementRegionManager const & e string const wellControlsName = wellControls.getName(); bool const logSurfaceCondition = isLogLevelActive< logInfo::WellControl >( wellControls.getLogLevel()); integer const useSurfaceConditions = wellControls.useSurfaceConditions(); - real64 flashPressure; + real64 flashPressure, flashTemperature; if( useSurfaceConditions ) { // use surface conditions flashPressure = wellControls.getSurfacePressure(); + flashTemperature = wellControls.getSurfaceTemperature(); } else { @@ -309,10 +312,12 @@ void SinglePhaseWell::updateVolRateForConstraint( ElementRegionManager const & e } // use region conditions flashPressure = wellControls.getRegionAveragePressure(); + flashTemperature = wellControls.getRegionAverageTemperature(); if( flashPressure < 0.0 ) { // use segment conditions flashPressure = pres[iwelemRef]; + flashTemperature = temp[iwelemRef]; } } real64 & currentVolRate = @@ -332,12 +337,14 @@ void SinglePhaseWell::updateVolRateForConstraint( ElementRegionManager const & e // bring everything back to host, capture the scalars by reference forAll< serialPolicy >( 1, [fluidWrapper, pres, + temp, connRate, dens, dDens, logSurfaceCondition, &useSurfaceConditions, &flashPressure, + &flashTemperature, ¤tVolRate, dCurrentVolRate, &iwelemRef, @@ -350,12 +357,19 @@ void SinglePhaseWell::updateVolRateForConstraint( ElementRegionManager const & e if( useSurfaceConditions ) { // we need to compute the surface density - fluidWrapper.update( iwelemRef, 0, flashPressure ); + if constexpr ( IS_THERMAL ) + { + fluidWrapper.update( iwelemRef, 0, flashPressure, flashTemperature ); + } + else + { + fluidWrapper.update( iwelemRef, 0, flashPressure ); + } if( logSurfaceCondition ) { - GEOS_LOG_RANK( GEOS_FMT( "{}: surface density computed with P_surface = {} Pa", - wellControlsName, flashPressure ) ); + GEOS_LOG_RANK( GEOS_FMT( "{}: surface density computed with P_surface = {} Pa, T_surface = {} K", + wellControlsName, flashPressure, flashTemperature ) ); } #ifdef GEOS_USE_HIP @@ -366,7 +380,15 @@ void SinglePhaseWell::updateVolRateForConstraint( ElementRegionManager const & e else { real64 const refPres = pres[iwelemRef]; - fluidWrapper.update( iwelemRef, 0, refPres ); + real64 const refTemp = temp[iwelemRef]; + if constexpr ( IS_THERMAL ) + { + fluidWrapper.update( iwelemRef, 0, refPres, refTemp ); + } + else + { + fluidWrapper.update( iwelemRef, 0, refPres, refTemp ); + } } real64 const densInv = 1.0 / dens[iwelemRef][0]; From c051dd7a4744c0f32e809adf175f6c5af46c2325 Mon Sep 17 00:00:00 2001 From: Tom Byer <149726499+tjb-ltk@users.noreply.github.com> Date: Tue, 24 Feb 2026 10:51:38 -0800 Subject: [PATCH 2/4] pr comment requests --- .../fluidFlow/wells/SinglePhaseWell.cpp | 22 +++++-------------- .../fluidFlow/wells/WellControls.cpp | 15 +++++++++++-- 2 files changed, 19 insertions(+), 18 deletions(-) diff --git a/src/coreComponents/physicsSolvers/fluidFlow/wells/SinglePhaseWell.cpp b/src/coreComponents/physicsSolvers/fluidFlow/wells/SinglePhaseWell.cpp index d4999d44d2e..0de3205d8ea 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/wells/SinglePhaseWell.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/wells/SinglePhaseWell.cpp @@ -357,14 +357,9 @@ void SinglePhaseWell::updateVolRateForConstraint( ElementRegionManager const & e if( useSurfaceConditions ) { // we need to compute the surface density - if constexpr ( IS_THERMAL ) - { - fluidWrapper.update( iwelemRef, 0, flashPressure, flashTemperature ); - } - else - { - fluidWrapper.update( iwelemRef, 0, flashPressure ); - } + + fluidWrapper.update( iwelemRef, 0, flashPressure, flashTemperature ); + if( logSurfaceCondition ) { @@ -381,14 +376,9 @@ void SinglePhaseWell::updateVolRateForConstraint( ElementRegionManager const & e { real64 const refPres = pres[iwelemRef]; real64 const refTemp = temp[iwelemRef]; - if constexpr ( IS_THERMAL ) - { - fluidWrapper.update( iwelemRef, 0, refPres, refTemp ); - } - else - { - fluidWrapper.update( iwelemRef, 0, refPres, refTemp ); - } + + fluidWrapper.update( iwelemRef, 0, refPres, refTemp ); + } real64 const densInv = 1.0 / dens[iwelemRef][0]; diff --git a/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.cpp b/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.cpp index 684306cce46..0a5535a46b5 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.cpp @@ -41,8 +41,8 @@ WellControls::WellControls( string const & name, Group * const parent ) m_targetPhaseRate( 0.0 ), m_targetMassRate( 0.0 ), m_useSurfaceConditions( 0 ), - m_surfacePres( 0.0 ), - m_surfaceTemp( 0.0 ), + m_surfacePres( -1.0 ), + m_surfaceTemp( -1.0 ), m_isCrossflowEnabled( 1 ), m_initialPressureCoefficient( 0.1 ), m_rateSign( -1.0 ), @@ -305,6 +305,17 @@ void WellControls::postInputInitialization() "The flag to select surface/reservoir conditions must be equal to 0 or 1", InputError, getWrapperDataContext( viewKeyStruct::useSurfaceConditionsString() ) ); + if( m_useSurfaceConditions ) + { + // 3a) check the flag for surface condition P&T + GEOS_THROW_IF( m_surfacePres <= 0.0 , + "Surface conditions set to 1 but surface pressure is not set", + InputError, getWrapperDataContext( viewKeyStruct::useSurfaceConditionsString() ) ); + GEOS_THROW_IF( m_surfaceTemp <= 0.0 , + "Surface conditions set to 1 but surface temperature is not set", + InputError, getWrapperDataContext( viewKeyStruct::useSurfaceConditionsString() ) ); + + } // 4) check that at least one rate constraint has been defined GEOS_THROW_IF( ((m_targetPhaseRate <= 0.0 && m_targetPhaseRateTableName.empty()) && (m_targetMassRate <= 0.0 && m_targetMassRateTableName.empty()) && From 7e469c75f6e6dff7c20fbf91da09469b22b6ebd8 Mon Sep 17 00:00:00 2001 From: Tom Byer <149726499+tjb-ltk@users.noreply.github.com> Date: Tue, 24 Feb 2026 12:04:31 -0800 Subject: [PATCH 3/4] uncrust --- .../physicsSolvers/fluidFlow/wells/WellControls.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.cpp b/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.cpp index 0a5535a46b5..120a28352bc 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/wells/WellControls.cpp @@ -308,10 +308,10 @@ void WellControls::postInputInitialization() if( m_useSurfaceConditions ) { // 3a) check the flag for surface condition P&T - GEOS_THROW_IF( m_surfacePres <= 0.0 , + GEOS_THROW_IF( m_surfacePres <= 0.0, "Surface conditions set to 1 but surface pressure is not set", InputError, getWrapperDataContext( viewKeyStruct::useSurfaceConditionsString() ) ); - GEOS_THROW_IF( m_surfaceTemp <= 0.0 , + GEOS_THROW_IF( m_surfaceTemp <= 0.0, "Surface conditions set to 1 but surface temperature is not set", InputError, getWrapperDataContext( viewKeyStruct::useSurfaceConditionsString() ) ); From f2dd4c070f4c7297e309445b31bd495a168ae00d Mon Sep 17 00:00:00 2001 From: Tom Byer <149726499+tjb-ltk@users.noreply.github.com> Date: Wed, 25 Feb 2026 04:42:30 -0800 Subject: [PATCH 4/4] add temp to sep def for wells --- .../poromechanics/PoroElastic_staircase_singlephase_3d_fim.xml | 1 + .../PoroElastic_staircase_singlephase_3d_sequential.xml | 1 + .../singlePhaseWell/staircase_single_phase_wells_hybrid_3d.xml | 1 + 3 files changed, 3 insertions(+) diff --git a/inputFiles/poromechanics/PoroElastic_staircase_singlephase_3d_fim.xml b/inputFiles/poromechanics/PoroElastic_staircase_singlephase_3d_fim.xml index 95f2090ada1..395292ba46d 100644 --- a/inputFiles/poromechanics/PoroElastic_staircase_singlephase_3d_fim.xml +++ b/inputFiles/poromechanics/PoroElastic_staircase_singlephase_3d_fim.xml @@ -63,6 +63,7 @@ targetBHP="5e7" targetTotalRateTableName="injectorTotalRateTable" useSurfaceConditions="1" + surfaceTemperature="300" surfacePressure="101325"/> diff --git a/inputFiles/poromechanics/PoroElastic_staircase_singlephase_3d_sequential.xml b/inputFiles/poromechanics/PoroElastic_staircase_singlephase_3d_sequential.xml index 0c952c6934f..af1007c6c6a 100755 --- a/inputFiles/poromechanics/PoroElastic_staircase_singlephase_3d_sequential.xml +++ b/inputFiles/poromechanics/PoroElastic_staircase_singlephase_3d_sequential.xml @@ -72,6 +72,7 @@ targetBHP="5e7" targetTotalRateTableName="injectorTotalRateTable" useSurfaceConditions="1" + surfaceTemperature="300" surfacePressure="101325"/> diff --git a/inputFiles/singlePhaseWell/staircase_single_phase_wells_hybrid_3d.xml b/inputFiles/singlePhaseWell/staircase_single_phase_wells_hybrid_3d.xml index 2298a32e483..52f243e506d 100644 --- a/inputFiles/singlePhaseWell/staircase_single_phase_wells_hybrid_3d.xml +++ b/inputFiles/singlePhaseWell/staircase_single_phase_wells_hybrid_3d.xml @@ -38,6 +38,7 @@ control="totalVolRate" useSurfaceConditions="1" surfacePressure="101325" + surfaceTemperature="300" referenceElevation="12" targetBHP="1e7" targetTotalRate="1e-6"/>