drop velocity seeding alg

This commit is contained in:
James McClure 2020-04-06 06:07:46 -04:00
parent b81d4199a1
commit 6d4eaebf47

View File

@ -56,8 +56,6 @@ void ScaLBL_ColorModel::ReadCheckpoint(char *FILENAME, double *cPhi, double *cfq
File.close();
}
*/
void ScaLBL_ColorModel::ReadParams(string filename){
// read the input database
db = std::make_shared<Database>( filename );
@ -1176,190 +1174,81 @@ double ScaLBL_ColorModel::MorphOpenConnected(double target_volume_change){
}
return(volume_change);
}
double ScaLBL_ColorModel::SeedPhaseField(const double seed_water_in_oil){
srand(time(NULL));
double mass_loss =0.f;
double count =0.f;
double *Aq_tmp, *Bq_tmp;
double *Vel_tmp;
Aq_tmp = new double [7*Np];
Bq_tmp = new double [7*Np];
Vel_tmp = new double [3*Np];
srand(time(NULL));
double mass_loss =0.f;
double count =0.f;
double *Aq_tmp, *Bq_tmp;
Aq_tmp = new double [7*Np];
Bq_tmp = new double [7*Np];
ScaLBL_CopyToHost(Aq_tmp, Aq, 7*Np*sizeof(double));
ScaLBL_CopyToHost(Bq_tmp, Bq, 7*Np*sizeof(double));
ScaLBL_CopyToHost(Vel_tmp, Velocity, 3*Np*sizeof(double));
//Extract averged velocity
double vx_glb = (Averages->gnb.Px+Averages->gwb.Px)/(Averages->gnb.M+Averages->gwb.M);
double vy_glb = (Averages->gnb.Py+Averages->gwb.Py)/(Averages->gnb.M+Averages->gwb.M);
double vz_glb = (Averages->gnb.Pz+Averages->gwb.Pz)/(Averages->gnb.M+Averages->gwb.M);
double v_mag_glb = sqrt(vx_glb*vx_glb+vy_glb*vy_glb+vz_glb*vz_glb);
ScaLBL_CopyToHost(Aq_tmp, Aq, 7*Np*sizeof(double));
ScaLBL_CopyToHost(Bq_tmp, Bq, 7*Np*sizeof(double));
for (int n=0; n < ScaLBL_Comm->LastExterior(); n++){
double random_value = seed_water_in_oil*double(rand())/ RAND_MAX;
double dA = Aq_tmp[n] + Aq_tmp[n+Np] + Aq_tmp[n+2*Np] + Aq_tmp[n+3*Np] + Aq_tmp[n+4*Np] + Aq_tmp[n+5*Np] + Aq_tmp[n+6*Np];
double dB = Bq_tmp[n] + Bq_tmp[n+Np] + Bq_tmp[n+2*Np] + Bq_tmp[n+3*Np] + Bq_tmp[n+4*Np] + Bq_tmp[n+5*Np] + Bq_tmp[n+6*Np];
double phase_id = (dA - dB) / (dA + dB);
if (phase_id > 0.0){
Aq_tmp[n] -= 0.3333333333333333*random_value;
Aq_tmp[n+Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+2*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+3*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+4*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+5*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+6*Np] -= 0.1111111111111111*random_value;
Bq_tmp[n] += 0.3333333333333333*random_value;
Bq_tmp[n+Np] += 0.1111111111111111*random_value;
Bq_tmp[n+2*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+3*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+4*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+5*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+6*Np] += 0.1111111111111111*random_value;
}
mass_loss += random_value*seed_water_in_oil;
}
for (int n=0; n < ScaLBL_Comm->LastExterior(); n++){
double v_mag_local = sqrt(Vel_tmp[n]*Vel_tmp[n]+Vel_tmp[n+1*Np]*Vel_tmp[n+1*Np]+Vel_tmp[n+2*Np]*Vel_tmp[n+2*Np]);
double weight = (v_mag_local<v_mag_glb) ? v_mag_local/v_mag_glb : 1.0;
double random_value = weight*seed_water_in_oil*double(rand())/ RAND_MAX;
double dA = Aq_tmp[n] + Aq_tmp[n+Np] + Aq_tmp[n+2*Np] + Aq_tmp[n+3*Np] + Aq_tmp[n+4*Np] + Aq_tmp[n+5*Np] + Aq_tmp[n+6*Np];
double dB = Bq_tmp[n] + Bq_tmp[n+Np] + Bq_tmp[n+2*Np] + Bq_tmp[n+3*Np] + Bq_tmp[n+4*Np] + Bq_tmp[n+5*Np] + Bq_tmp[n+6*Np];
double phase_id = (dA - dB) / (dA + dB);
if (phase_id > 0.0){
Aq_tmp[n] -= 0.3333333333333333*random_value;
Aq_tmp[n+Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+2*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+3*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+4*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+5*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+6*Np] -= 0.1111111111111111*random_value;
Bq_tmp[n] += 0.3333333333333333*random_value;
Bq_tmp[n+Np] += 0.1111111111111111*random_value;
Bq_tmp[n+2*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+3*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+4*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+5*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+6*Np] += 0.1111111111111111*random_value;
count += 1.0;
}
mass_loss += random_value*seed_water_in_oil;
}
for (int n=ScaLBL_Comm->FirstInterior(); n < ScaLBL_Comm->LastInterior(); n++){
double random_value = seed_water_in_oil*double(rand())/ RAND_MAX;
double dA = Aq_tmp[n] + Aq_tmp[n+Np] + Aq_tmp[n+2*Np] + Aq_tmp[n+3*Np] + Aq_tmp[n+4*Np] + Aq_tmp[n+5*Np] + Aq_tmp[n+6*Np];
double dB = Bq_tmp[n] + Bq_tmp[n+Np] + Bq_tmp[n+2*Np] + Bq_tmp[n+3*Np] + Bq_tmp[n+4*Np] + Bq_tmp[n+5*Np] + Bq_tmp[n+6*Np];
double phase_id = (dA - dB) / (dA + dB);
if (phase_id > 0.0){
Aq_tmp[n] -= 0.3333333333333333*random_value;
Aq_tmp[n+Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+2*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+3*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+4*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+5*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+6*Np] -= 0.1111111111111111*random_value;
Bq_tmp[n] += 0.3333333333333333*random_value;
Bq_tmp[n+Np] += 0.1111111111111111*random_value;
Bq_tmp[n+2*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+3*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+4*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+5*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+6*Np] += 0.1111111111111111*random_value;
}
mass_loss += random_value*seed_water_in_oil;
}
for (int n=ScaLBL_Comm->FirstInterior(); n < ScaLBL_Comm->LastInterior(); n++){
double v_mag_local = sqrt(Vel_tmp[n]*Vel_tmp[n]+Vel_tmp[n+1*Np]*Vel_tmp[n+1*Np]+Vel_tmp[n+2*Np]*Vel_tmp[n+2*Np]);
double weight = (v_mag_local<v_mag_glb) ? v_mag_local/v_mag_glb : 1.0;
double random_value = weight*seed_water_in_oil*double(rand())/ RAND_MAX;
double dA = Aq_tmp[n] + Aq_tmp[n+Np] + Aq_tmp[n+2*Np] + Aq_tmp[n+3*Np] + Aq_tmp[n+4*Np] + Aq_tmp[n+5*Np] + Aq_tmp[n+6*Np];
double dB = Bq_tmp[n] + Bq_tmp[n+Np] + Bq_tmp[n+2*Np] + Bq_tmp[n+3*Np] + Bq_tmp[n+4*Np] + Bq_tmp[n+5*Np] + Bq_tmp[n+6*Np];
double phase_id = (dA - dB) / (dA + dB);
if (phase_id > 0.0){
Aq_tmp[n] -= 0.3333333333333333*random_value;
Aq_tmp[n+Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+2*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+3*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+4*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+5*Np] -= 0.1111111111111111*random_value;
Aq_tmp[n+6*Np] -= 0.1111111111111111*random_value;
Bq_tmp[n] += 0.3333333333333333*random_value;
Bq_tmp[n+Np] += 0.1111111111111111*random_value;
Bq_tmp[n+2*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+3*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+4*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+5*Np] += 0.1111111111111111*random_value;
Bq_tmp[n+6*Np] += 0.1111111111111111*random_value;
count += 1.0;
}
mass_loss += random_value*seed_water_in_oil;
}
count= sumReduce( Dm->Comm, count);
mass_loss= sumReduce( Dm->Comm, mass_loss);
if (rank == 0) printf("Remove mass %f from %f voxels \n",mass_loss,count);
count= sumReduce( Dm->Comm, count);
mass_loss= sumReduce( Dm->Comm, mass_loss);
if (rank == 0) printf("Remove mass %f from %f voxels \n",mass_loss,count);
// Need to initialize Aq, Bq, Den, Phi directly
//ScaLBL_CopyToDevice(Phi,phase.data(),7*Np*sizeof(double));
ScaLBL_CopyToDevice(Aq, Aq_tmp, 7*Np*sizeof(double));
ScaLBL_CopyToDevice(Bq, Bq_tmp, 7*Np*sizeof(double));
// Need to initialize Aq, Bq, Den, Phi directly
//ScaLBL_CopyToDevice(Phi,phase.data(),7*Np*sizeof(double));
ScaLBL_CopyToDevice(Aq, Aq_tmp, 7*Np*sizeof(double));
ScaLBL_CopyToDevice(Bq, Bq_tmp, 7*Np*sizeof(double));
return(mass_loss);
return(mass_loss);
}
//double ScaLBL_ColorModel::SeedPhaseField(const double seed_water_in_oil){
// srand(time(NULL));
// double mass_loss =0.f;
// double count =0.f;
// double *Aq_tmp, *Bq_tmp;
//
// Aq_tmp = new double [7*Np];
// Bq_tmp = new double [7*Np];
//
// ScaLBL_CopyToHost(Aq_tmp, Aq, 7*Np*sizeof(double));
// ScaLBL_CopyToHost(Bq_tmp, Bq, 7*Np*sizeof(double));
//
///* for (int k=1; k<Nz-1; k++){
// for (int j=1; j<Ny-1; j++){
// for (int i=1; i<Nx-1; i++){
// double random_value = double(rand())/ RAND_MAX;
//
// if (Averages->SDs(i,j,k) < 0.f){
// // skip
// }
// else if (phase(i,j,k) > 0.f ){
// phase(i,j,k) -= random_value*seed_water_in_oil;
// mass_loss += random_value*seed_water_in_oil;
// count++;
// }
// else {
//
// }
// }
// }
// }
// */
// for (int n=0; n < ScaLBL_Comm->LastExterior(); n++){
// double random_value = seed_water_in_oil*double(rand())/ RAND_MAX;
// double dA = Aq_tmp[n] + Aq_tmp[n+Np] + Aq_tmp[n+2*Np] + Aq_tmp[n+3*Np] + Aq_tmp[n+4*Np] + Aq_tmp[n+5*Np] + Aq_tmp[n+6*Np];
// double dB = Bq_tmp[n] + Bq_tmp[n+Np] + Bq_tmp[n+2*Np] + Bq_tmp[n+3*Np] + Bq_tmp[n+4*Np] + Bq_tmp[n+5*Np] + Bq_tmp[n+6*Np];
// double phase_id = (dA - dB) / (dA + dB);
// if (phase_id > 0.0){
// Aq_tmp[n] -= 0.3333333333333333*random_value;
// Aq_tmp[n+Np] -= 0.1111111111111111*random_value;
// Aq_tmp[n+2*Np] -= 0.1111111111111111*random_value;
// Aq_tmp[n+3*Np] -= 0.1111111111111111*random_value;
// Aq_tmp[n+4*Np] -= 0.1111111111111111*random_value;
// Aq_tmp[n+5*Np] -= 0.1111111111111111*random_value;
// Aq_tmp[n+6*Np] -= 0.1111111111111111*random_value;
//
// Bq_tmp[n] += 0.3333333333333333*random_value;
// Bq_tmp[n+Np] += 0.1111111111111111*random_value;
// Bq_tmp[n+2*Np] += 0.1111111111111111*random_value;
// Bq_tmp[n+3*Np] += 0.1111111111111111*random_value;
// Bq_tmp[n+4*Np] += 0.1111111111111111*random_value;
// Bq_tmp[n+5*Np] += 0.1111111111111111*random_value;
// Bq_tmp[n+6*Np] += 0.1111111111111111*random_value;
// }
// mass_loss += random_value*seed_water_in_oil;
// }
//
// for (int n=ScaLBL_Comm->FirstInterior(); n < ScaLBL_Comm->LastInterior(); n++){
// double random_value = seed_water_in_oil*double(rand())/ RAND_MAX;
// double dA = Aq_tmp[n] + Aq_tmp[n+Np] + Aq_tmp[n+2*Np] + Aq_tmp[n+3*Np] + Aq_tmp[n+4*Np] + Aq_tmp[n+5*Np] + Aq_tmp[n+6*Np];
// double dB = Bq_tmp[n] + Bq_tmp[n+Np] + Bq_tmp[n+2*Np] + Bq_tmp[n+3*Np] + Bq_tmp[n+4*Np] + Bq_tmp[n+5*Np] + Bq_tmp[n+6*Np];
// double phase_id = (dA - dB) / (dA + dB);
// if (phase_id > 0.0){
// Aq_tmp[n] -= 0.3333333333333333*random_value;
// Aq_tmp[n+Np] -= 0.1111111111111111*random_value;
// Aq_tmp[n+2*Np] -= 0.1111111111111111*random_value;
// Aq_tmp[n+3*Np] -= 0.1111111111111111*random_value;
// Aq_tmp[n+4*Np] -= 0.1111111111111111*random_value;
// Aq_tmp[n+5*Np] -= 0.1111111111111111*random_value;
// Aq_tmp[n+6*Np] -= 0.1111111111111111*random_value;
//
// Bq_tmp[n] += 0.3333333333333333*random_value;
// Bq_tmp[n+Np] += 0.1111111111111111*random_value;
// Bq_tmp[n+2*Np] += 0.1111111111111111*random_value;
// Bq_tmp[n+3*Np] += 0.1111111111111111*random_value;
// Bq_tmp[n+4*Np] += 0.1111111111111111*random_value;
// Bq_tmp[n+5*Np] += 0.1111111111111111*random_value;
// Bq_tmp[n+6*Np] += 0.1111111111111111*random_value;
// }
// mass_loss += random_value*seed_water_in_oil;
// }
//
// count = Dm->Comm.sumReduce( count );
// mass_loss = Dm->Comm.sumReduce( mass_loss );
// if (rank == 0) printf("Remove mass %f from %f voxels \n",mass_loss,count);
//
// // Need to initialize Aq, Bq, Den, Phi directly
// //ScaLBL_CopyToDevice(Phi,phase.data(),7*Np*sizeof(double));
// ScaLBL_CopyToDevice(Aq, Aq_tmp, 7*Np*sizeof(double));
// ScaLBL_CopyToDevice(Bq, Bq_tmp, 7*Np*sizeof(double));
//
// return(mass_loss);
//}
double ScaLBL_ColorModel::MorphInit(const double beta, const double target_delta_volume){
const RankInfoStruct rank_info(rank,nprocx,nprocy,nprocz);