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

refactoring MRT model

Browse Source
This commit is contained in:
James E McClure
2018-07-28 16:01:24 -04:00
parent 1e3f3cd466
commit 29442f44f3
6 changed files with 83 additions and 23 deletions
Download Patch File Download Diff File Expand all files Collapse all files

16
analysis/Minkowski.cpp
View File

@@ -67,6 +67,22 @@ void Minkowski::Initialize()
vol_n = euler = Jn = An = Kn = 0.0; vol_n = euler = Jn = An = Kn = 0.0;
} }
double Minkowski::V(){
return vol_n_global;
}
double Minkowski::A(){
return An_global;
}
double Minkowski::J(){
return Jn_global;
}
double Minkowski::X(){
return euler_global_global;
}
void Minkowski::UpdateMeshValues() void Minkowski::UpdateMeshValues()
{ {

6
analysis/Minkowski.h
View File

@@ -69,6 +69,11 @@ public:
DoubleArray SDn_y; DoubleArray SDn_y;
DoubleArray SDn_z; DoubleArray SDn_z;
double V();
double A();
double J();
double X();
//........................................................................... //...........................................................................
Minkowski(std::shared_ptr <Domain> Dm); Minkowski(std::shared_ptr <Domain> Dm);
~Minkowski(); ~Minkowski();
@@ -80,6 +85,7 @@ public:
void PrintAll(); void PrintAll();
int GetCubeLabel(int i, int j, int k, IntArray &BlobLabel); int GetCubeLabel(int i, int j, int k, IntArray &BlobLabel);
void SortBlobs(); void SortBlobs();
}; };

2
common/ScaLBL.cpp
View File

@@ -1121,14 +1121,12 @@ void ScaLBL_Communicator::RecvGrad(double *phi, double *grad){
ScaLBL_DeviceBarrier(); ScaLBL_DeviceBarrier();
//................................................................................... //...................................................................................
// NOTE: AA Routine writes to opposite
// Unpack the gradributions on the device // Unpack the gradributions on the device
//................................................................................... //...................................................................................
//...Unpacking for x face(2,8,10,12,14)................................ //...Unpacking for x face(2,8,10,12,14)................................
ScaLBL_Gradient_Unpack(1.0,-1,0,0,dvcRecvDist_x,0,recvCount_x,recvbuf_x,phi,grad,N); ScaLBL_Gradient_Unpack(1.0,-1,0,0,dvcRecvDist_x,0,recvCount_x,recvbuf_x,phi,grad,N);
ScaLBL_Gradient_Unpack(0.5,-1,-1,0,dvcRecvDist_x,recvCount_x,recvCount_x,recvbuf_x,phi,grad,N); ScaLBL_Gradient_Unpack(0.5,-1,-1,0,dvcRecvDist_x,recvCount_x,recvCount_x,recvbuf_x,phi,grad,N);
ScaLBL_Gradient_Unpack(0.5,-1,1,0,dvcRecvDist_x,2*recvCount_x,recvCount_x,recvbuf_x,phi,grad,N); ScaLBL_Gradient_Unpack(0.5,-1,1,0,dvcRecvDist_x,2*recvCount_x,recvCount_x,recvbuf_x,phi,grad,N);
ScaLBL_Gradient_Unpack(0.5,-1,0,-1,dvcRecvDist_x,3*recvCount_x,recvCount_x,recvbuf_x,phi,grad,N);
ScaLBL_Gradient_Unpack(0.5,-1,0,1,dvcRecvDist_x,4*recvCount_x,recvCount_x,recvbuf_x,phi,grad,N); ScaLBL_Gradient_Unpack(0.5,-1,0,1,dvcRecvDist_x,4*recvCount_x,recvCount_x,recvbuf_x,phi,grad,N);
//................................................................................... //...................................................................................
//...Packing for X face(1,7,9,11,13)................................ //...Packing for X face(1,7,9,11,13)................................

24
models/ColorModel.cpp
View File

@@ -246,8 +246,8 @@ void ScaLBL_ColorModel::Create(){
} }
} }
// check that TmpMap is valid // check that TmpMap is valid
for (int idx=0; idx<ScaLBL_Comm->last_interior; idx++){ for (int idx=0; idx<ScaLBL_Comm->LastInterior(); idx++){
if (idx == ScaLBL_Comm->next) idx = ScaLBL_Comm->first_interior; if (idx == ScaLBL_Comm->LastExterior()) idx = ScaLBL_Comm->FirstInterior();
int n = TmpMap[idx]; int n = TmpMap[idx];
if (n > Nx*Ny*Nz){ if (n > Nx*Ny*Nz){
printf("Bad value! idx=%i \n"); printf("Bad value! idx=%i \n");
@@ -316,8 +316,8 @@ void ScaLBL_ColorModel::Initialize(){
if (rank==0) printf ("Initializing distributions \n"); if (rank==0) printf ("Initializing distributions \n");
ScaLBL_D3Q19_Init(fq, Np); ScaLBL_D3Q19_Init(fq, Np);
if (rank==0) printf ("Initializing phase field \n"); if (rank==0) printf ("Initializing phase field \n");
ScaLBL_PhaseField_Init(dvcMap, Phi, Den, Aq, Bq, 0, ScaLBL_Comm->next, Np); ScaLBL_PhaseField_Init(dvcMap, Phi, Den, Aq, Bq, 0, ScaLBL_Comm->LastExterior(), Np);
ScaLBL_PhaseField_Init(dvcMap, Phi, Den, Aq, Bq, ScaLBL_Comm->first_interior, ScaLBL_Comm->last_interior, Np); ScaLBL_PhaseField_Init(dvcMap, Phi, Den, Aq, Bq, ScaLBL_Comm->FirstInterior(), ScaLBL_Comm->LastInterior(), Np);
if (BoundaryCondition >0 ){ if (BoundaryCondition >0 ){
if (Dm->kproc()==0){ if (Dm->kproc()==0){
@@ -365,9 +365,9 @@ void ScaLBL_ColorModel::Run(){
// Compute the Phase indicator field // Compute the Phase indicator field
// Read for Aq, Bq happens in this routine (requires communication) // Read for Aq, Bq happens in this routine (requires communication)
ScaLBL_Comm->BiSendD3Q7AA(Aq,Bq); //READ FROM NORMAL ScaLBL_Comm->BiSendD3Q7AA(Aq,Bq); //READ FROM NORMAL
ScaLBL_D3Q7_AAodd_PhaseField(NeighborList, dvcMap, Aq, Bq, Den, Phi, ScaLBL_Comm->first_interior, ScaLBL_Comm->last_interior, Np); ScaLBL_D3Q7_AAodd_PhaseField(NeighborList, dvcMap, Aq, Bq, Den, Phi, ScaLBL_Comm->FirstInterior(), ScaLBL_Comm->LastInterior(), Np);
ScaLBL_Comm->BiRecvD3Q7AA(Aq,Bq); //WRITE INTO OPPOSITE ScaLBL_Comm->BiRecvD3Q7AA(Aq,Bq); //WRITE INTO OPPOSITE
ScaLBL_D3Q7_AAodd_PhaseField(NeighborList, dvcMap, Aq, Bq, Den, Phi, 0, ScaLBL_Comm->next, Np); ScaLBL_D3Q7_AAodd_PhaseField(NeighborList, dvcMap, Aq, Bq, Den, Phi, 0, ScaLBL_Comm->LastExterior(), Np);
// Perform the collision operation // Perform the collision operation
ScaLBL_Comm->SendD3Q19AA(fq); //READ FROM NORMAL ScaLBL_Comm->SendD3Q19AA(fq); //READ FROM NORMAL
@@ -375,7 +375,7 @@ void ScaLBL_ColorModel::Run(){
ScaLBL_Comm_Regular->SendHalo(Phi); ScaLBL_Comm_Regular->SendHalo(Phi);
ScaLBL_D3Q19_AAodd_Color(NeighborList, dvcMap, fq, Aq, Bq, Den, Phi, Velocity, rhoA, rhoB, tauA, tauB, ScaLBL_D3Q19_AAodd_Color(NeighborList, dvcMap, fq, Aq, Bq, Den, Phi, Velocity, rhoA, rhoB, tauA, tauB,
alpha, beta, Fx, Fy, Fz, Nx, Nx*Ny, ScaLBL_Comm->first_interior, ScaLBL_Comm->last_interior, Np); alpha, beta, Fx, Fy, Fz, Nx, Nx*Ny, ScaLBL_Comm->FirstInterior(), ScaLBL_Comm->LastInterior(), Np);
ScaLBL_Comm_Regular->RecvHalo(Phi); ScaLBL_Comm_Regular->RecvHalo(Phi);
ScaLBL_Comm->RecvD3Q19AA(fq); //WRITE INTO OPPOSITE ScaLBL_Comm->RecvD3Q19AA(fq); //WRITE INTO OPPOSITE
// Set BCs // Set BCs
@@ -392,23 +392,23 @@ void ScaLBL_ColorModel::Run(){
ScaLBL_Comm->D3Q19_Pressure_BC_Z(NeighborList, fq, dout, timestep); ScaLBL_Comm->D3Q19_Pressure_BC_Z(NeighborList, fq, dout, timestep);
} }
ScaLBL_D3Q19_AAodd_Color(NeighborList, dvcMap, fq, Aq, Bq, Den, Phi, Velocity, rhoA, rhoB, tauA, tauB, ScaLBL_D3Q19_AAodd_Color(NeighborList, dvcMap, fq, Aq, Bq, Den, Phi, Velocity, rhoA, rhoB, tauA, tauB,
alpha, beta, Fx, Fy, Fz, Nx, Nx*Ny, 0, ScaLBL_Comm->next, Np); alpha, beta, Fx, Fy, Fz, Nx, Nx*Ny, 0, ScaLBL_Comm->LastExterior(), Np);
ScaLBL_DeviceBarrier(); MPI_Barrier(comm); ScaLBL_DeviceBarrier(); MPI_Barrier(comm);
// *************EVEN TIMESTEP************* // *************EVEN TIMESTEP*************
timestep++; timestep++;
// Compute the Phase indicator field // Compute the Phase indicator field
ScaLBL_Comm->BiSendD3Q7AA(Aq,Bq); //READ FROM NORMAL ScaLBL_Comm->BiSendD3Q7AA(Aq,Bq); //READ FROM NORMAL
ScaLBL_D3Q7_AAeven_PhaseField(dvcMap, Aq, Bq, Den, Phi, ScaLBL_Comm->first_interior, ScaLBL_Comm->last_interior, Np); ScaLBL_D3Q7_AAeven_PhaseField(dvcMap, Aq, Bq, Den, Phi, ScaLBL_Comm->FirstInterior(), ScaLBL_Comm->LastInterior(), Np);
ScaLBL_Comm->BiRecvD3Q7AA(Aq,Bq); //WRITE INTO OPPOSITE ScaLBL_Comm->BiRecvD3Q7AA(Aq,Bq); //WRITE INTO OPPOSITE
ScaLBL_D3Q7_AAeven_PhaseField(dvcMap, Aq, Bq, Den, Phi, 0, ScaLBL_Comm->next, Np); ScaLBL_D3Q7_AAeven_PhaseField(dvcMap, Aq, Bq, Den, Phi, 0, ScaLBL_Comm->LastExterior(), Np);
// Perform the collision operation // Perform the collision operation
ScaLBL_Comm->SendD3Q19AA(fq); //READ FORM NORMAL ScaLBL_Comm->SendD3Q19AA(fq); //READ FORM NORMAL
// Halo exchange for phase field // Halo exchange for phase field
ScaLBL_Comm_Regular->SendHalo(Phi); ScaLBL_Comm_Regular->SendHalo(Phi);
ScaLBL_D3Q19_AAeven_Color(dvcMap, fq, Aq, Bq, Den, Phi, Velocity, rhoA, rhoB, tauA, tauB, ScaLBL_D3Q19_AAeven_Color(dvcMap, fq, Aq, Bq, Den, Phi, Velocity, rhoA, rhoB, tauA, tauB,
alpha, beta, Fx, Fy, Fz, Nx, Nx*Ny, ScaLBL_Comm->first_interior, ScaLBL_Comm->last_interior, Np); alpha, beta, Fx, Fy, Fz, Nx, Nx*Ny, ScaLBL_Comm->FirstInterior(), ScaLBL_Comm->LastInterior(), Np);
ScaLBL_Comm_Regular->RecvHalo(Phi); ScaLBL_Comm_Regular->RecvHalo(Phi);
ScaLBL_Comm->RecvD3Q19AA(fq); //WRITE INTO OPPOSITE ScaLBL_Comm->RecvD3Q19AA(fq); //WRITE INTO OPPOSITE
// Set boundary conditions // Set boundary conditions
@@ -425,7 +425,7 @@ void ScaLBL_ColorModel::Run(){
ScaLBL_Comm->D3Q19_Pressure_BC_Z(NeighborList, fq, dout, timestep); ScaLBL_Comm->D3Q19_Pressure_BC_Z(NeighborList, fq, dout, timestep);
} }
ScaLBL_D3Q19_AAeven_Color(dvcMap, fq, Aq, Bq, Den, Phi, Velocity, rhoA, rhoB, tauA, tauB, ScaLBL_D3Q19_AAeven_Color(dvcMap, fq, Aq, Bq, Den, Phi, Velocity, rhoA, rhoB, tauA, tauB,
alpha, beta, Fx, Fy, Fz, Nx, Nx*Ny, 0, ScaLBL_Comm->next, Np); alpha, beta, Fx, Fy, Fz, Nx, Nx*Ny, 0, ScaLBL_Comm->LastExterior(), Np);
ScaLBL_DeviceBarrier(); MPI_Barrier(comm); ScaLBL_DeviceBarrier(); MPI_Barrier(comm);
//************************************************************************ //************************************************************************

55
models/MRTModel.cpp
View File

@@ -124,6 +124,8 @@ void ScaLBL_MRTModel::Create(){
// copy the neighbor list // copy the neighbor list
ScaLBL_CopyToDevice(NeighborList, neighborList, neighborSize); ScaLBL_CopyToDevice(NeighborList, neighborList, neighborSize);
MPI_Barrier(comm); MPI_Barrier(comm);
Morphology = new Minkowski(Mask);
} }
void ScaLBL_MRTModel::Initialize(){ void ScaLBL_MRTModel::Initialize(){
@@ -162,21 +164,56 @@ void ScaLBL_MRTModel::Run(){
} }
//************************************************************************/ //************************************************************************/
stoptime = MPI_Wtime(); stoptime = MPI_Wtime();
// cout << "CPU time: " << (stoptime - starttime) << " seconds" << endl; if (rank==0) printf("-------------------------------------------------------------------\n");
cputime = stoptime - starttime; // Compute the walltime per timestep
// cout << "Lattice update rate: "<< double(Nx*Ny*Nz*timestep)/cputime/1000000 << " MLUPS" << endl; cputime = (stoptime - starttime)/timestep;
double MLUPS = double(Np*timestep)/cputime/1000000; // Performance obtained from each node
// if (rank==0) printf("********************************************************\n"); double MLUPS = double(Np)/cputime/1000000;
// if (rank==0) printf("CPU time = %f \n", cputime);
// if (rank==0) printf("Lattice update rate (per process)= %f MLUPS \n", MLUPS); if (rank==0) printf("********************************************************\n");
if (rank==0) printf("CPU time = %f \n", cputime);
if (rank==0) printf("Lattice update rate (per core)= %f MLUPS \n", MLUPS);
MLUPS *= nprocs; MLUPS *= nprocs;
if (rank==0) printf("Lattice update rate (total)= %f MLUPS \n", MLUPS);
if (rank==0) printf("********************************************************\n");
} }
void ScaLBL_MRTModel::VelocityField(double *Vz){ void ScaLBL_MRTModel::VelocityField(double *VELOCITY){
int SIZE=Np*sizeof(double); int SIZE=Np*sizeof(double);
ScaLBL_D3Q19_Momentum(fq,Velocity, Np); ScaLBL_D3Q19_Momentum(fq,Velocity, Np);
ScaLBL_DeviceBarrier(); MPI_Barrier(comm); ScaLBL_DeviceBarrier(); MPI_Barrier(comm);
ScaLBL_CopyToHost(&Vz[0],&Velocity[2*Np],SIZE); ScaLBL_CopyToHost(&VELOCITY[0],&Velocity[0],SIZE);
memcpy(Morphology.SDn.data(), Distance.data(), N*sizeof(double));
Morphology.Initialize();
Morphology.UpdateMeshValues();
Morphology.ComputeLocal();
Morphology.Reduce();
unsigned long int count_loc=0;
unsigned long int count;
double vax,vay,vaz;
double vax_loc,vay_loc,vaz_loc;
vax_loc = vay_loc = vaz_loc = 0.f;
for (int n=ScaLBL_Comm->0; n<ScaLBL_Comm->LastExterior(); n++){
vax_loc += VELOCITY[n];
vay_loc += VELOCITY[Np+n];
vaz_loc += VELOCITY[2*Np+n];
}
for (int n=ScaLBL_Comm->FirstInterior(); n<ScaLBL_Comm->LastInterior(); n++){
vax_loc += VELOCITY[n];
vay_loc += VELOCITY[Np+n];
vaz_loc += VELOCITY[2*Np+n];
}
MPI_Allreduce(vax_loc,vax,1,MPI_DOUBLE,MPI_SUM,Mask->Comm);
MPI_Allreduce(vay_loc,vay,1,MPI_DOUBLE,MPI_SUM,Mask->Comm);
MPI_Allreduce(vaz_loc,vaz,1,MPI_DOUBLE,MPI_SUM,Mask->Comm);
if (rank==0) printf("Vs As Js Xs vx vy vz\n");
if (rank==0) printf("%.8g %.8g %.8g %.8g %.8g %.8g %.8g\n",Morphology.V(),Morphology.A(),Morphology.J(),Morphology.X(),vax,vay,vaz);
} }

3
models/MRTModel.h
View File

@@ -12,6 +12,7 @@
#include "common/ScaLBL.h" #include "common/ScaLBL.h"
#include "common/Communication.h" #include "common/Communication.h"
#include "common/MPI_Helpers.h" #include "common/MPI_Helpers.h"
#include "analysis/Minkowski.h"
#include "ProfilerApp.h" #include "ProfilerApp.h"
class ScaLBL_MRTModel{ class ScaLBL_MRTModel{
@@ -55,6 +56,8 @@ public:
double *Velocity; double *Velocity;
double *Pressure; double *Pressure;
Minkowski Morphology;
private: private:
MPI_Comm comm; MPI_Comm comm;