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

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

Browse Source
This commit is contained in:
JamesEMcclure
2019-10-01 14:27:07 -04:00
parent 017f06039e c52ff20519
commit cc0f04a3d3
5 changed files with 136 additions and 92 deletions
Download Patch File Download Diff File Expand all files Collapse all files

81
analysis/SubPhase.cpp
View File

@@ -689,30 +689,11 @@ void SubPhase::AggregateLabels(char *FILENAME){
int ny = Dm->Ny;
int nz = Dm->Nz;
int npx = Dm->nprocx();
int npy = Dm->nprocy();
int npz = Dm->nprocz();
int ipx = Dm->iproc();
int ipy = Dm->jproc();
int ipz = Dm->kproc();
int nprocs = Dm->nprocx()*Dm->nprocy()*Dm->nprocz();
int full_nx = npx*(nx-2);
int full_ny = npy*(ny-2);
int full_nz = npz*(nz-2);
int local_size = (nx-2)*(ny-2)*(nz-2);
long int full_size = long(full_nx)*long(full_ny)*long(full_nz);
signed char *LocalID;
LocalID = new signed char [local_size];
//printf("aggregate labels: local size=%i, global size = %i",local_size, full_size);
// assign the ID for the local sub-region
for (int k=1; k<nz-1; k++){
for (int j=1; j<ny-1; j++){
for (int i=1; i<nx-1; i++){
// assign the ID from the phase indicator field
for (int k=0; k<nz; k++){
for (int j=0; j<ny; j++){
for (int i=0; i<nx; i++){
int n = k*nx*ny+j*nx+i;
signed char local_id_val = Dm->id[n];
if (local_id_val > 0){
@@ -720,63 +701,13 @@ void SubPhase::AggregateLabels(char *FILENAME){
if (value > 0.0) local_id_val = 1;
else local_id_val = 2;
}
LocalID[(k-1)*(nx-2)*(ny-2) + (j-1)*(nx-2) + i-1] = local_id_val;
Dm->id[n] = local_id_val;
}
}
}
MPI_Barrier(Dm->Comm);
// populate the FullID
if (Dm->rank() == 0){
signed char *FullID;
FullID = new signed char [full_size];
// first handle local ID for rank 0
for (int k=1; k<nz-1; k++){
for (int j=1; j<ny-1; j++){
for (int i=1; i<nx-1; i++){
int x = i-1;
int y = j-1;
int z = k-1;
int n_local = (k-1)*(nx-2)*(ny-2) + (j-1)*(nx-2) + i-1;
int n_full = z*full_nx*full_ny + y*full_nx + x;
FullID[n_full] = LocalID[n_local];
}
}
}
// next get the local ID from the other ranks
for (int rnk = 1; rnk<nprocs; rnk++){
ipz = rnk / (npx*npy);
ipy = (rnk - ipz*npx*npy) / npx;
ipx = (rnk - ipz*npx*npy - ipy*npx);
//printf("ipx=%i ipy=%i ipz=%i\n", ipx, ipy, ipz);
int tag = 15+rnk;
MPI_Recv(LocalID,local_size,MPI_CHAR,rnk,tag,Dm->Comm,MPI_STATUS_IGNORE);
for (int k=1; k<nz-1; k++){
for (int j=1; j<ny-1; j++){
for (int i=1; i<nx-1; i++){
int x = i-1 + ipx*(nx-2);
int y = j-1 + ipy*(ny-2);
int z = k-1 + ipz*(nz-2);
int n_local = (k-1)*(nx-2)*(ny-2) + (j-1)*(nx-2) + i-1;
int n_full = z*full_nx*full_ny + y*full_nx + x;
FullID[n_full] = LocalID[n_local];
}
}
}
}
// write the output
FILE *OUTFILE;
OUTFILE = fopen(FILENAME,"wb");
fwrite(FullID,1,full_size,OUTFILE);
fclose(OUTFILE);
}
else{
// send LocalID to rank=0
int tag = 15+ Dm->rank();
int dstrank = 0;
MPI_Send(LocalID,local_size,MPI_CHAR,dstrank,tag,Dm->Comm);
}
MPI_Barrier(Dm->Comm);
Dm->AggregateLabels(FILENAME);
}

100
common/Domain.cpp
View File

@@ -377,10 +377,10 @@ void Domain::Decomp(std::string Filename)
// relabel the data
std::vector<long int> LabelCount(ReadValues.size(),0);
for (int k = 0; k<Nz; k++){
for (int j = 0; j<Ny; j++){
for (int i = 0; i<Nx; i++){
n = k*Nx*Ny+j*Nx+i;
for (int k = 0; k<global_Nz; k++){
for (int j = 0; j<global_Ny; j++){
for (int i = 0; i<global_Nx; i++){
n = k*global_Nx*global_Ny+j*global_Nx+i;
//char locval = loc_id[n];
char locval = SegData[n];
for (int idx=0; idx<ReadValues.size(); idx++){
@@ -595,6 +595,98 @@ void Domain::Decomp(std::string Filename)
MPI_Barrier(Comm);
}
void Domain::AggregateLabels(char *FILENAME){
int nx = Nx;
int ny = Ny;
int nz = Nz;
int npx = nprocx();
int npy = nprocy();
int npz = nprocz();
int ipx = iproc();
int ipy = jproc();
int ipz = kproc();
int nprocs = nprocx()*nprocy()*nprocz();
int full_nx = npx*(nx-2);
int full_ny = npy*(ny-2);
int full_nz = npz*(nz-2);
int local_size = (nx-2)*(ny-2)*(nz-2);
long int full_size = long(full_nx)*long(full_ny)*long(full_nz);
signed char *LocalID;
LocalID = new signed char [local_size];
//printf("aggregate labels: local size=%i, global size = %i",local_size, full_size);
// assign the ID for the local sub-region
for (int k=1; k<nz-1; k++){
for (int j=1; j<ny-1; j++){
for (int i=1; i<nx-1; i++){
int n = k*nx*ny+j*nx+i;
signed char local_id_val = id[n];
LocalID[(k-1)*(nx-2)*(ny-2) + (j-1)*(nx-2) + i-1] = local_id_val;
}
}
}
MPI_Barrier(Comm);
// populate the FullID
if (rank() == 0){
signed char *FullID;
FullID = new signed char [full_size];
// first handle local ID for rank 0
for (int k=1; k<nz-1; k++){
for (int j=1; j<ny-1; j++){
for (int i=1; i<nx-1; i++){
int x = i-1;
int y = j-1;
int z = k-1;
int n_local = (k-1)*(nx-2)*(ny-2) + (j-1)*(nx-2) + i-1;
int n_full = z*full_nx*full_ny + y*full_nx + x;
FullID[n_full] = LocalID[n_local];
}
}
}
// next get the local ID from the other ranks
for (int rnk = 1; rnk<nprocs; rnk++){
ipz = rnk / (npx*npy);
ipy = (rnk - ipz*npx*npy) / npx;
ipx = (rnk - ipz*npx*npy - ipy*npx);
//printf("ipx=%i ipy=%i ipz=%i\n", ipx, ipy, ipz);
int tag = 15+rnk;
MPI_Recv(LocalID,local_size,MPI_CHAR,rnk,tag,Comm,MPI_STATUS_IGNORE);
for (int k=1; k<nz-1; k++){
for (int j=1; j<ny-1; j++){
for (int i=1; i<nx-1; i++){
int x = i-1 + ipx*(nx-2);
int y = j-1 + ipy*(ny-2);
int z = k-1 + ipz*(nz-2);
int n_local = (k-1)*(nx-2)*(ny-2) + (j-1)*(nx-2) + i-1;
int n_full = z*full_nx*full_ny + y*full_nx + x;
FullID[n_full] = LocalID[n_local];
}
}
}
}
// write the output
FILE *OUTFILE;
OUTFILE = fopen(FILENAME,"wb");
fwrite(FullID,1,full_size,OUTFILE);
fclose(OUTFILE);
}
else{
// send LocalID to rank=0
int tag = 15+ rank();
int dstrank = 0;
MPI_Send(LocalID,local_size,MPI_CHAR,dstrank,tag,Comm);
}
MPI_Barrier(Comm);
}
/********************************************************
* Initialize communication *
********************************************************/

1
common/Domain.h
View File

@@ -197,6 +197,7 @@ public: // Public variables (need to create accessors instead)
void CommunicateMeshHalo(DoubleArray &Mesh);
void CommInit();
int PoreCount();
void AggregateLabels(char *FILENAME);
private:

44
tests/lbpm_morphdrain_pp.cpp
View File

@@ -32,19 +32,13 @@ int main(int argc, char **argv)
//.......................................................................
// Reading the domain information file
//.......................................................................
int nprocx, nprocy, nprocz, nx, ny, nz, nspheres;
double Lx, Ly, Lz;
int i,j,k,n;
int BC=0;
// char fluidValue,solidValue;
int MAXTIME=1000;
int READ_FROM_BLOCK=0;
int n, nprocx, nprocy, nprocz, nx, ny, nz;
char LocalRankString[8];
char LocalRankFilename[40];
char FILENAME[128];
string filename;
double Rcrit_new, SW;
double SW,Rcrit_new;
if (argc > 1){
filename=argv[1];
Rcrit_new=0.f;
@@ -61,6 +55,7 @@ int main(int argc, char **argv)
auto ReadValues = domain_db->getVector<int>( "ReadValues" );
auto WriteValues = domain_db->getVector<int>( "WriteValues" );
SW = domain_db->getScalar<double>("Sw");
auto READFILE = domain_db->getScalar<std::string>( "Filename" );
// Generate the NWP configuration
//if (rank==0) printf("Initializing morphological distribution with critical radius %f \n", Rcrit);
@@ -77,19 +72,24 @@ int main(int argc, char **argv)
int N = (nx+2)*(ny+2)*(nz+2);
std::shared_ptr<Domain> Dm (new Domain(domain_db,comm));
std::shared_ptr<Domain> Mask (new Domain(domain_db,comm));
// std::shared_ptr<Domain> Dm (new Domain(nx,ny,nz,rank,nprocx,nprocy,nprocz,Lx,Ly,Lz,BC));
for (n=0; n<N; n++) Dm->id[n]=1;
Dm->CommInit();
signed char *id;
id = new signed char [N];
sprintf(LocalRankFilename,"ID.%05i",rank);
Mask->Decomp(READFILE);
Mask->CommInit();
/*sprintf(LocalRankFilename,"ID.%05i",rank);
size_t readID;
FILE *IDFILE = fopen(LocalRankFilename,"rb");
if (IDFILE==NULL) ERROR("Error opening file: ID.xxxxx");
readID=fread(id,1,N,IDFILE);
if (readID != size_t(N)) printf("lbpm_morphopen_pp: Error reading ID (rank=%i) \n",rank);
fclose(IDFILE);
*/
nx+=2; ny+=2; nz+=2;
// Generate the signed distance map
@@ -102,9 +102,13 @@ int main(int argc, char **argv)
for (int j=0;j<ny;j++){
for (int i=0;i<nx;i++){
int n = k*nx*ny+j*nx+i;
id[n] = Mask->id[n];
// Initialize the solid phase
if (id[n] > 0) id_solid(i,j,k) = 1;
else id_solid(i,j,k) = 0;
if (Mask->id[n] > 0){
id_solid(i,j,k) = 1;
}
else
id_solid(i,j,k) = 0;
}
}
}
@@ -188,6 +192,22 @@ int main(int argc, char **argv)
FILE *ID = fopen(LocalRankFilename,"wb");
fwrite(id,1,N,ID);
fclose(ID);
// write the geometry to a single file
for (int k=0;k<nz;k++){
for (int j=0;j<ny;j++){
for (int i=0;i<nx;i++){
int n = k*nx*ny+j*nx+i;
Mask->id[n] = id[n];
}
}
}
MPI_Barrier(comm);
sprintf(FILENAME,READFILE.c_str());
sprintf(FILENAME+strlen(FILENAME),".morphdrain.raw");
if (rank==0) printf("Writing file to: %s \n", FILENAME);
Mask->AggregateLabels(FILENAME);
}
MPI_Barrier(comm);