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

This commit is contained in:
JamesEMcclure 2019-09-06 11:42:57 -04:00
commit 3dcd8c6d1e
2 changed files with 42 additions and 13 deletions

View File

@ -175,6 +175,23 @@ void SubPhase::Basic(){
double nA,nB;
double count_w = 0.0;
double count_n = 0.0;
for (k=0; k<Nz; k++){
for (j=0; j<Ny; j++){
for (i=0; i<Nx; i++){
n = k*Nx*Ny + j*Nx + i;
// Compute volume averages
if ( Dm->id[n] > 0 ){
// compute density
double nA = Rho_n(n);
double nB = Rho_w(n);
double phi = (nA-nB)/(nA+nB);
Phi(n) = phi;
}
}
}
}
for (k=kmin; k<kmax; k++){
for (j=jmin; j<Ny-1; j++){
for (i=imin; i<Nx-1; i++){
@ -184,10 +201,9 @@ void SubPhase::Basic(){
// compute density
double nA = Rho_n(n);
double nB = Rho_w(n);
double phi = (nA-nB)/(nA+nB);
Phi(n) = phi;
double phi = (nA-nB)/(nA+nB);
if ( phi > 0.0 ){
nA = 1.0;
nb.V += 1.0;
nb.M += nA*rho_n;
// velocity
@ -196,6 +212,7 @@ void SubPhase::Basic(){
nb.Pz += rho_n*nA*Vel_z(n);
}
else{
nB = 1.0;
wb.M += nB*rho_w;
wb.V += 1.0;
@ -522,8 +539,7 @@ void SubPhase::Full(){
double ux = Vel_x(n);
double uy = Vel_y(n);
double uz = Vel_z(n);
Phi(n) = phi;
if (DelPhi(n) > 1e-3){
// interface region
double nx = 0.5*(Phi(i+1,j,k)-Phi(i-1,j,k));
@ -554,6 +570,7 @@ void SubPhase::Full(){
}
if ( phi > 0.0){
if (morph_n->label(i,j,k) > 0 ){
nA = 1.0;
nd.M += nA*rho_n;
nd.Px += nA*rho_n*ux;
nd.Py += nA*rho_n*uy;
@ -561,6 +578,7 @@ void SubPhase::Full(){
nd.K += nA*rho_n*(ux*ux + uy*uy + uz*uz);
}
else{
nA = 1.0;
nc.M += nA*rho_n;
nc.Px += nA*rho_n*ux;
nc.Py += nA*rho_n*uy;
@ -571,6 +589,7 @@ void SubPhase::Full(){
else{
// water region
if (morph_w->label(i,j,k) > 0 ){
nB = 1.0;
wd.M += nB*rho_w;
wd.Px += nB*rho_w*ux;
wd.Py += nB*rho_w*uy;
@ -578,6 +597,7 @@ void SubPhase::Full(){
wd.K += nB*rho_w*(ux*ux + uy*uy + uz*uz);
}
else{
nB = 1.0;
wc.M += nB*rho_w;
wc.Px += nB*rho_w*ux;
wc.Py += nB*rho_w*uy;
@ -626,6 +646,11 @@ void SubPhase::Full(){
giwn.Kw=sumReduce( Dm->Comm, iwn.Kw);
// pressure averaging
gnc.p=sumReduce( Dm->Comm, nc.p);
gnd.p=sumReduce( Dm->Comm, nd.p);
gwc.p=sumReduce( Dm->Comm, wc.p);
gwd.p=sumReduce( Dm->Comm, wd.p);
if (vol_wc_bulk > 0.0)
wc.p = wc.p /vol_wc_bulk;
if (vol_nc_bulk > 0.0)

View File

@ -819,14 +819,18 @@ void ScaLBL_ColorModel::Run(){
double flow_rate_A_connected = Vol_nc*(vAc_x*dir_x + vAc_y*dir_y + vAc_z*dir_z);
double flow_rate_B_connected = Vol_wc*(vBc_x*dir_x + vBc_y*dir_y + vBc_z*dir_z);
double flow_rate_A_disconnected = Vol_nd*(vAd_x*dir_x + vAd_y*dir_y + vAd_z*dir_z);
double flow_rate_B_disconnected = Vol_wd*(vBd_x*dir_x + vBd_y*dir_y + vBd_z*dir_z);
double flow_rate_A_disconnected = (Vol_nd)*(vAd_x*dir_x + vAd_y*dir_y + vAd_z*dir_z);
double flow_rate_B_disconnected = (Vol_wd)*(vBd_x*dir_x + vBd_y*dir_y + vBd_z*dir_z);
double kAeff_connected = h*h*muA*flow_rate_A_connected/(rhoA*force_mag);
double kBeff_connected = h*h*muB*flow_rate_B_connected/(rhoB*force_mag);
double kAeff_connected = h*h*muA*flow_rate_A_connected/(force_mag);
double kBeff_connected = h*h*muB*flow_rate_B_connected/(force_mag);
double kAeff = h*h*muA*(flow_rate_A_connected+flow_rate_A_disconnected)/(rhoA*force_mag);
double kBeff = h*h*muB*(flow_rate_B_connected+flow_rate_B_disconnected)/(rhoB*force_mag);
double kAeff_disconnected = h*h*muA*flow_rate_A_disconnected/(rhoA*force_mag);
double kBeff_disconnected = h*h*muB*flow_rate_B_disconnected/(rhoB*force_mag);
double kAeff = h*h*muA*(flow_rate_A)/(rhoA*force_mag);
double kBeff = h*h*muB*(flow_rate_B)/(rhoB*force_mag);
double viscous_pressure_drop = (rhoA*volA + rhoB*volB)*force_mag;
double Mobility = muA/muB;
@ -838,9 +842,9 @@ void ScaLBL_ColorModel::Run(){
WriteHeader=true;
kr_log_file = fopen("relperm.csv","a");
if (WriteHeader)
fprintf(kr_log_file,"timesteps sat.water eff.perm.oil eff.perm.water eff.perm.oil.connected eff.perm.water.connected cap.pressure cap.pressure.connected pressure.drop Ca M\n",CURRENT_STEADY_TIMESTEPS,current_saturation,kAeff,kBeff,pAB,viscous_pressure_drop,Ca,Mobility);
fprintf(kr_log_file,"timesteps sat.water eff.perm.oil eff.perm.water eff.perm.oil.connected eff.perm.water.connected eff.perm.oil.disconnected eff.perm.water.disconnected cap.pressure cap.pressure.connected 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 %.5g %.5g %.5g\n",CURRENT_STEADY_TIMESTEPS,current_saturation,kAeff,kBeff,kAeff_connected,kBeff_connected,pAB,pAB_connected,viscous_pressure_drop,Ca,Mobility);
fprintf(kr_log_file,"%i %.5g %.5g %.5g %.5g %.5g %.5g %.5g %.5g %.5g %.5g %.5g %.5g %.5g\n",CURRENT_STEADY_TIMESTEPS,current_saturation,kAeff,kBeff,kAeff_connected,kBeff_connected,kAeff_disconnected,kBeff_disconnected,pAB,pAB_connected,viscous_pressure_drop,Ca,Mobility);
fclose(kr_log_file);
printf(" Measured capillary number %f \n ",Ca);