Skip to content

fix: SinglePhaseWell Separator flash #3978

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
tjb-ltk wants to merge 1 commit into
base: develop
Choose a base branch
from
Open

fix: SinglePhaseWell Separator flash #3978

Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
32 changes: 27 additions & 5 deletions src/coreComponents/physicsSolvers/fluidFlow/wells/SinglePhaseWell.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -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 >();

Expand All @@ -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
{
Expand All @@ -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 =
Expand All @@ -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,
&currentVolRate,
dCurrentVolRate,
&iwelemRef,
Expand All @@ -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
Expand All @@ -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];
Expand Down
Loading