OilfieldToolsNet/LBPM
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

fractional flow protocl with endpoint termination

Browse Source
This commit is contained in:
James McClure 2021-06-13 10:08:32 -04:00
parent 9d968ee2ca
commit 8327ff9d38
1 changed files with 25 additions and 12 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

37
tests/lbpm_color_simulator.cpp
View File

@ -9,8 +9,6 @@
#include "common/Utilities.h" #include "common/Utilities.h"
#include "models/ColorModel.h" #include "models/ColorModel.h"
//#define WRE_SURFACES
/* /*
* Simulator for two-phase flow in porous media * Simulator for two-phase flow in porous media
* James E. McClure 2013-2014 * James E. McClure 2013-2014
@ -71,29 +69,33 @@ int main( int argc, char **argv )
if (SimulationMode == "development"){ if (SimulationMode == "development"){
double MLUPS=0.0; double MLUPS=0.0;
int timestep = 0; int timestep = 0;
bool ContinueSimulation = true;
int ANALYSIS_INTERVAL = ColorModel.timestepMax;
/* flow adaptor keys to control */ /* flow adaptor keys to control */
int SKIP_TIMESTEPS = 0; int SKIP_TIMESTEPS = 0;
int MAX_STEADY_TIME = 1000000; int MAX_STEADY_TIME = 1000000;
double ENDPOINT_THRESHOLD = 0.1;
double FRACTIONAL_FLOW_INCREMENT = 0.05; double FRACTIONAL_FLOW_INCREMENT = 0.05;
if (ColorModel.db->keyExists( "FlowAdaptor" )){ if (ColorModel.db->keyExists( "FlowAdaptor" )){
auto flow_db = ColorModel.db->getDatabase( "FlowAdaptor" ); auto flow_db = ColorModel.db->getDatabase( "FlowAdaptor" );
MAX_STEADY_TIME = flow_db->getWithDefault<int>( "max_steady_timesteps", 1000000 ); MAX_STEADY_TIME = flow_db->getWithDefault<int>( "max_steady_timesteps", 1000000 );
SKIP_TIMESTEPS = flow_db->getWithDefault<int>( "skip_timesteps", 100000 ); SKIP_TIMESTEPS = flow_db->getWithDefault<int>( "skip_timesteps", 100000 );
FRACTIONAL_FLOW_INCREMENT = flow_db->getWithDefault<double>( "fractional_flow_increment", 0.05); FRACTIONAL_FLOW_INCREMENT = flow_db->getWithDefault<double>( "fractional_flow_increment", 0.05);
ENDPOINT_THRESHOLD = flow_db->getWithDefault<double>( "endpoint_threshold", 0.1);
} }
int ANALYSIS_INTERVAL = ColorModel.timestepMax;
if (ColorModel.analysis_db->keyExists( "analysis_interval" )){ if (ColorModel.analysis_db->keyExists( "analysis_interval" )){
ANALYSIS_INTERVAL = ColorModel.analysis_db->getScalar<int>( "analysis_interval" ); ANALYSIS_INTERVAL = ColorModel.analysis_db->getScalar<int>( "analysis_interval" );
} }
/* Launch the simulation */ /* Launch the simulation */
FlowAdaptor Adapt(ColorModel); FlowAdaptor Adapt(ColorModel);
runAnalysis analysis(ColorModel); runAnalysis analysis(ColorModel);
while (ColorModel.timestep < ColorModel.timestepMax){ if (ContinueSimulation){
/* this will run steady points */ /* this will run steady points */
timestep += MAX_STEADY_TIME; timestep += MAX_STEADY_TIME;
MLUPS = ColorModel.Run(timestep); MLUPS = ColorModel.Run(timestep);
if (rank==0) printf("Lattice update rate (per MPI process)= %f MLUPS \n", MLUPS); if (rank==0) printf("Lattice update rate (per MPI process)= %f MLUPS \n", MLUPS);
if (ColorModel.timestep > ColorModel.timestepMax) ContinueSimulation = false;
/* update the fractional flow by adding mass */ /* update the fractional flow by adding mass */
int skip_time = 0; int skip_time = 0;
timestep = ColorModel.timestep; timestep = ColorModel.timestep;
@ -101,6 +103,17 @@ int main( int argc, char **argv )
double volB = ColorModel.Averages->gwb.V; double volB = ColorModel.Averages->gwb.V;
double volA = ColorModel.Averages->gnb.V; double volA = ColorModel.Averages->gnb.V;
double initialSaturation = volB/(volA + volB); double initialSaturation = volB/(volA + volB);
double vA_x = ColorModel.Averages->gnb.Px/ColorModel.Averages->gnb.M;
double vA_y = ColorModel.Averages->gnb.Py/ColorModel.Averages->gnb.M;
double vA_z = ColorModel.Averages->gnb.Pz/ColorModel.Averages->gnb.M;
double vB_x = ColorModel.Averages->gwb.Px/ColorModel.Averages->gwb.M;
double vB_y = ColorModel.Averages->gwb.Py/ColorModel.Averages->gwb.M;
double vB_z = ColorModel.Averages->gwb.Pz/ColorModel.Averages->gwb.M;
double speedA = sqrt(vA_x*vA_x + vA_y*vA_y + vA_z*vA_z);
double speedB = sqrt(vB_x*vB_x + vB_y*vB_y + vB_z*vB_z);
if (volA*speedA < ENDPOINT_THRESHOLD*volB*speedB) ContinueSimulation = false;
if (ContinueSimulation){
while (skip_time < SKIP_TIMESTEPS && fabs(SaturationChange) < fabs(FRACTIONAL_FLOW_INCREMENT) ){ while (skip_time < SKIP_TIMESTEPS && fabs(SaturationChange) < fabs(FRACTIONAL_FLOW_INCREMENT) ){
timestep += ANALYSIS_INTERVAL; timestep += ANALYSIS_INTERVAL;
Adapt.UpdateFractionalFlow(ColorModel); Adapt.UpdateFractionalFlow(ColorModel);
@ -109,11 +122,11 @@ int main( int argc, char **argv )
double volA = ColorModel.Averages->gnb.V; double volA = ColorModel.Averages->gnb.V;
SaturationChange = volB/(volA + volB) - initialSaturation; SaturationChange = volB/(volA + volB) - initialSaturation;
skip_time += ANALYSIS_INTERVAL; skip_time += ANALYSIS_INTERVAL;
} }
if (rank==0) printf(" ********************************************************************* \n"); if (rank==0) printf(" ********************************************************************* \n");
if (rank==0) printf(" Updated fractional flow with saturation change = %f \n", SaturationChange); if (rank==0) printf(" Updated fractional flow with saturation change = %f \n", SaturationChange);
if (rank==0) printf(" ********************************************************************* \n"); if (rank==0) printf(" ********************************************************************* \n");
}
/* apply timestep skipping algorithm to accelerate steady-state */ /* apply timestep skipping algorithm to accelerate steady-state */
/* skip_time = 0; /* skip_time = 0;
timestep = ColorModel.timestep; timestep = ColorModel.timestep;
@ -125,7 +138,7 @@ int main( int argc, char **argv )
} }
*/ */
} }
ColorModel.WriteDebug(); //ColorModel.WriteDebug();
} }
else else