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

Updated morphdrain preprocessor to take a target saturation as input

Browse Source
This commit is contained in:
James E McClure 2016-11-07 14:28:39 -05:00
parent 55df5b0c4d
commit 802dbb461e
1 changed files with 146 additions and 124 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

38
tests/lbpm_morphdrain_pp.cpp
View File

@ -53,10 +53,11 @@ int main(int argc, char **argv)
MPI_Comm_rank(comm,&rank);
MPI_Comm_size(comm,&nprocs);
double Rcrit;
Rcrit=strtod(argv[1],NULL);
//double Rcrit;
double SW=strtod(argv[1],NULL);
if (rank==0){
printf("Critical radius %f (voxels)\n",Rcrit);
//printf("Critical radius %f (voxels)\n",Rcrit);
printf("Target saturation %f \n",SW);
}
// }
@ -192,7 +193,8 @@ int main(int argc, char **argv)
int kproc = Dm.kproc;
// Generate the NWP configuration
if (rank==0) printf("Performing morphological drainage with critical radius %f \n", Rcrit);
//if (rank==0) printf("Performing morphological drainage with critical radius %f \n", Rcrit);
if (rank==0) printf("Performing morphological drainage with target saturation %f \n", SW);
// GenerateResidual(id,nx,ny,nz,Saturation);
// Communication buffers
@ -250,10 +252,11 @@ int main(int argc, char **argv)
int Nx = nx;
int Ny = ny;
int Nz = nz;
float sat = 0.f;
double sw=1.f;
int GlobalNumber = 1;
int Window=int(Rcrit);
double radius,Rcrit;
radius = 0.0;
// Layer the inlet with NWP
if (kproc == 0){
for(j=0; j<Ny; j++){
@ -261,11 +264,29 @@ int main(int argc, char **argv)
n = j*nx+i;
// n = nx*ny + j*nx+i;
id[n]=1;
if (SignDist(i,j,k) > radius){
radius=SignDist(i,j,k);
}
}
}
}
MPI_Allreduce(&radius,&Rcrit,1,MPI_DOUBLE,MPI_MAX,comm);
int Window=int(Rcrit);
if (rank==0) printf("Starting morhpological drainage with critical radius = %f \n",Rcrit);
int imin,jmin,kmin,imax,jmax,kmax;
// Decrease the critical radius until the target saturation is met
double deltaR=0.05; // amount to change the radius in voxel units
while (sw > SW){
// decrease critical radius
Rcrit -= deltaR;
Window=int(Rcrit);
GlobalNumber = 1;
while (GlobalNumber != 0){
if (rank==0) printf("GlobalNumber=%i \n",GlobalNumber);
@ -416,8 +437,9 @@ int main(int argc, char **argv)
}
}
MPI_Allreduce(&count,&countGlobal,1,MPI_INT,MPI_SUM,comm);
sat = float(countGlobal)/totalGlobal;
if (rank==0) printf("Final saturation=%f\n",sat);
sw= double(countGlobal)/totalGlobal;
if (rank==0) printf("Final saturation=%f\n",sw);
}
sprintf(LocalRankFilename,"ID.%05i",rank);
FILE *ID = fopen(LocalRankFilename,"wb");