Merge branch 'morphLBM' of github.com:JamesEMcClure/LBPM-WIA into morphLBM

This commit is contained in:
James E McClure 2019-04-26 11:00:52 -04:00
commit 088f50c6b1
4 changed files with 98 additions and 24 deletions

View File

@ -249,11 +249,12 @@ void SubPhase::Basic(){
double not_water_flow_rate=gnb.V*(gnb.Px*dir_x + gnb.Py*dir_y + gnb.Pz*dir_z)/gnb.M/ Dm->Volume;
double total_flow_rate = water_flow_rate + not_water_flow_rate;
double fractional_flow= water_flow_rate / total_flow_rate;
double krn = nu_n*not_water_flow_rate / force_mag ;
double krw = nu_w*water_flow_rate / force_mag;
double h = Dm->voxel_length;
double krn = h*h*nu_n*not_water_flow_rate / force_mag ;
double krw = h*h*nu_w*water_flow_rate / force_mag;
//printf(" water saturation = %f, fractional flow =%f \n",saturation,fractional_flow);
fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %.5g %.5g %.5g\n",saturation,krw,krn,water_flow_rate,not_water_flow_rate, gwb.p, gnb.p);
fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %.5g %.5g %.5g\n",saturation,krw,krn,h*h*h*water_flow_rate,h*h*h*not_water_flow_rate, gwb.p, gnb.p);
fflush(TIMELOG);
}

View File

@ -0,0 +1,49 @@
Color {
tauA = 0.7;
tauB = 0.7;
rhoA = 1.0;
rhoB = 1.0;
alpha = 1e-3;
beta = 0.95;
F = 0, 0, 0
Restart = false
pBC = 0
din = 1.0
dout = 1.0
timestepMax = 3000
interval = 1000
tol = 1e-5;
das = 0.1
dbs = 0.9
flux = 0.0
ComponentLabels = 0, -1 // labels for immobile components
ComponentAffinity = -1.0, 0.5 // wetting condition for immobile components (-1 = water wet / +1 oil wet)
}
Domain {
Filename = "CornerFlow.raw"
nproc = 1, 1, 2 // Number of processors (Npx,Npy,Npz)
n = 64, 64, 64 // Size of local domain (Nx,Ny,Nz)
N = 64, 64, 128 // size of the input image
n_spheres = 1 // Number of spheres
L = 1, 1, 1 // Length of domain (x,y,z)
BC = 0 // Boundary condition type
ReadType = "8bit"
ReadValues = 0, 1, 2 // labels within the original input image
WriteValues = 0, 1, 2 // labels to write (if they aren't the same)
ComponentLabels = 0 // labels that are immobile
HistoryLabels = -1 // new labels to assign to each immobile component based on fluid history
Sw = 0.3 // target saturation for morphological routines
}
Analysis {
blobid_interval = 1000 // Frequency to perform blob identification
analysis_interval = 1000 // Frequency to perform analysis
restart_interval = 1000 // Frequency to write restart data
visualization_interval = 2000 // Frequency to write visualization data
restart_file = "Restart" // Filename to use for restart file (will append rank)
N_threads = 4 // Number of threads to use
load_balance = "independent" // Load balance method to use: "none", "default", "independent"
}

View File

@ -626,6 +626,8 @@ void ScaLBL_ColorModel::Run(){
if ( morph_timesteps > morph_interval ){
double volB = Averages->gwb.V;
double volA = Averages->gnb.V;
volA /= double((Nx-2)*(Ny-2)*(Nz-2)*nprocs);
volB /= double((Nx-2)*(Ny-2)*(Nz-2)*nprocs);
double vA_x = Averages->gnb.Px/Averages->gnb.M;
double vA_y = Averages->gnb.Py/Averages->gnb.M;
double vA_z = Averages->gnb.Pz/Averages->gnb.M;
@ -639,12 +641,18 @@ void ScaLBL_ColorModel::Run(){
double dir_x = Fx/force_mag;
double dir_y = Fy/force_mag;
double dir_z = Fz/force_mag;
double flow_rate_A = (vA_x*dir_x + vA_y*dir_y + vA_z*dir_z);
double flow_rate_B = (vB_x*dir_x + vB_y*dir_y + vB_z*dir_z);
if (force_mag == 0.0){
// default to z direction
dir_x = 0.0;
dir_y = 0.0;
dir_z = 1.0;
force_mag = 1.0;
}
double current_saturation = volB/(volA+volB);
double Ca = fabs(volA*muA*flow_rate_A + volB*muB*flow_rate_B)/(5.796*alpha*double((Nx-2)*(Ny-2)*(Nz-2)*nprocs));
double force_magnitude = sqrt(Fx*Fx + Fy*Fy + Fz*Fz);
double flow_rate_A = volA*(vA_x*dir_x + vA_y*dir_y + vA_z*dir_z);
double flow_rate_B = volB*(vB_x*dir_x + vB_y*dir_y + vB_z*dir_z);
double Ca = fabs(muA*flow_rate_A + muB*flow_rate_B)/(5.796*alpha);
bool isSteady = false;
if ( (fabs((Ca - Ca_previous)/Ca) < tolerance && CURRENT_STEADY_TIMESTEPS > MIN_STEADY_TIMESTEPS))
@ -664,11 +672,27 @@ void ScaLBL_ColorModel::Run(){
if (rank==0){
printf("** WRITE STEADY POINT *** ");
printf("Ca = %f, (previous = %f) \n",Ca,Ca_previous);
volA /= double((Nx-2)*(Ny-2)*(Nz-2)*nprocs);
volB /= double((Nx-2)*(Ny-2)*(Nz-2)*nprocs);
FILE * kr_log_file = fopen("relperm.csv","a");
fprintf(kr_log_file,"%i %.5g %.5g %.5g %.5g %.5g %.5g ",CURRENT_STEADY_TIMESTEPS,muA,muB,5.796*alpha,Fx,Fy,Fz);
fprintf(kr_log_file,"%.5g %.5g %.5g %.5g %.5g %.5g %.5g %.5g\n",volA,volB,vA_x,vA_y,vA_z,vB_x,vB_y,vB_z);
double pA = Averages->gnb.p;
double pB = Averages->gwb.p;
double h = Dm->voxel_length;
double kAeff = h*h*muA*flow_rate_A/(rhoA*force_mag);
double kBeff = h*h*muB*flow_rate_B/(rhoB*force_mag);
double pAB = (pA-pB)/(h*5.796*alpha);
double viscous_pressure_drop = (rhoA*volA + rhoB*volB)*force_mag;
double Mobility = muA/muB;
bool WriteHeader=false;
FILE * kr_log_file = fopen("relperm.csv","r");
if (kr_log_file != NULL)
fclose(kr_log_file);
else
WriteHeader=true;
kr_log_file = fopen("relperm.csv","a");
if (WriteHeader)
fprintf(kr_log_file,"timesteps sat.water eff.perm.oil eff.perm.water cap.pressure pressure.drop Ca M\n",CURRENT_STEADY_TIMESTEPS,current_saturation,kAeff,kBeff,pAB,viscous_pressure_drop,Ca,Mobility);
fprintf(kr_log_file,"%i %.5g %.5g %.5g %.5g %.5g %.5g %.5g\n",CURRENT_STEADY_TIMESTEPS,current_saturation,kAeff,kBeff,pAB,viscous_pressure_drop,Ca,Mobility);
fclose(kr_log_file);
printf(" Measured capillary number %f \n ",Ca);
@ -679,15 +703,15 @@ void ScaLBL_ColorModel::Run(){
Fy *= capillary_number / Ca;
Fz *= capillary_number / Ca;
if (force_magnitude > 1e-3){
Fx *= 1e-3/force_magnitude; // impose ceiling for stability
Fy *= 1e-3/force_magnitude;
Fz *= 1e-3/force_magnitude;
if (force_mag > 1e-3){
Fx *= 1e-3/force_mag; // impose ceiling for stability
Fy *= 1e-3/force_mag;
Fz *= 1e-3/force_mag;
}
if (force_magnitude < 1e-7){
Fx *= 1e-7/force_magnitude; // impose floor
Fy *= 1e-7/force_magnitude;
Fz *= 1e-7/force_magnitude;
if (force_mag < 1e-7){
Fx *= 1e-7/force_mag; // impose floor
Fy *= 1e-7/force_mag;
Fz *= 1e-7/force_mag;
}
if (rank == 0) printf(" -- adjust force by factor %f \n ",capillary_number / Ca);
Averages->SetParams(rhoA,rhoB,tauA,tauB,Fx,Fy,Fz,alpha,beta);
@ -702,7 +726,7 @@ void ScaLBL_ColorModel::Run(){
morph_timesteps=0;
Ca_previous = Ca;
}
if (MORPH_ADAPT ){
if (MORPH_ADAPT && morph_delta != 0.0){
CURRENT_MORPH_TIMESTEPS += analysis_interval;
if (rank==0) printf("***Morphological step with target volume change %f ***\n", delta_volume_target);
//double delta_volume_target = volB - (volA + volB)*TARGET_SATURATION; // change in volume to A

View File

@ -279,7 +279,7 @@ void ScaLBL_MRTModel::Run(){
printf(" %f\n",absperm);
FILE * log_file = fopen("Permeability.csv","a");
fprintf(log_file,"%i %.8g %.8g %.8g %.8g %.8g %.8g %.8g %.8g %.8g %.8g %.8g %.8g\n",timestep, Fx, Fy, Fz, mu,
Vs,As,Hs,Xs,vax,vay,vaz, absperm);
h*h*h*Vs,h*h*As,h*Hs,Xs,vax,vay,vaz, absperm);
fclose(log_file);
}
}