dfh in progress

This commit is contained in:
James E McClure
2018-04-30 10:27:38 -04:00
parent 6c54cbe4c8
commit 1173082022
4 changed files with 256 additions and 17 deletions

View File

@@ -5,6 +5,46 @@
#define NBLOCKS 1024
#define NTHREADS 256
__global__ void dvc_ScaLBL_DFH_Init(double *Phi, double *Den, double *Aq, double *Bq, int start, int finish, int Np){
int idx,n;
double phi,nA,nB;
int S = Np/NBLOCKS/NTHREADS + 1;
for (int s=0; s<S; s++){
//........Get 1-D index for this thread....................
idx = S*blockIdx.x*blockDim.x + s*blockDim.x + threadIdx.x+start;
if (idx<finish) {
phi = Phi[idx];
if (phi > 0.f){
nA = 1.0; nB = 0.f;
}
else{
nB = 1.0; nA = 0.f;
}
Den[idx] = nA;
Den[Np+idx] = nB;
Aq[idx]=0.3333333333333333*nA;
Aq[Np+idx]=0.1111111111111111*nA;
Aq[2*Np+idx]=0.1111111111111111*nA;
Aq[3*Np+idx]=0.1111111111111111*nA;
Aq[4*Np+idx]=0.1111111111111111*nA;
Aq[5*Np+idx]=0.1111111111111111*nA;
Aq[6*Np+idx]=0.1111111111111111*nA;
Bq[idx]=0.3333333333333333*nB;
Bq[Np+idx]=0.1111111111111111*nB;
Bq[2*Np+idx]=0.1111111111111111*nB;
Bq[3*Np+idx]=0.1111111111111111*nB;
Bq[4*Np+idx]=0.1111111111111111*nB;
Bq[5*Np+idx]=0.1111111111111111*nB;
Bq[6*Np+idx]=0.1111111111111111*nB;
}
}
}
// LBM based on density functional hydrodynamics
__global__ void dvc_ScaLBL_D3Q19_AAeven_DFH(int *neighborList, double *dist, double *Aq, double *Bq, double *Den, double *Phi,
double *SolidPotential, double rhoA, double rhoB, double tauA, double tauB, double alpha, double beta,
@@ -1359,6 +1399,75 @@ __global__ void dvc_ScaLBL_D3Q7_AAeven_DFH(double *Aq, double *Bq, double *Den,
}
}
__global__ void dvc_ScaLBL_D3Q19_Gradient_DFH(int *neighborList, double *Phi, double *ColorGrad, int start, int finish, int Np){
int n,nn;
// distributions
double m1,m2,m3,m4,m5,m6,m7,m8,m9;
double m10,m11,m12,m13,m14,m15,m16,m17,m18;
double nx,ny,nz;
int S = Np/NBLOCKS/NTHREADS + 1;
for (int s=0; s<S; s++){
//........Get 1-D index for this thread....................
n = S*blockIdx.x*blockDim.x + s*blockDim.x + threadIdx.x + start;
if (n<finish) {
nn = neighborList[n+Np]%Np;
m1 = Phi[nn];
nn = neighborList[n]%Np;
m2 = Phi[nn];
nn = neighborList[n+3*Np]%Np;
m3 = Phi[nn];
nn = neighborList[n+2*Np]%Np;
m4 = Phi[nn];
nn = neighborList[n+5*Np]%Np;
m5 = Phi[nn];
nn = neighborList[n+4*Np]%Np;
m6 = Phi[nn];
nn = neighborList[n+7*Np]%Np;
m7 = Phi[nn];
nn = neighborList[n+6*Np]%Np;
m8 = Phi[nn];
nn = neighborList[n+9*Np]%Np;
m9 = Phi[nn];
nn = neighborList[n+8*Np]%Np;
m10 = Phi[nn];
nn = neighborList[n+11*Np]%Np;
m11 = Phi[nn];
nn = neighborList[n+10*Np]%Np;
m12 = Phi[nn];
nn = neighborList[n+13*Np]%Np;
m13 = Phi[nn];
nn = neighborList[n+12*Np]%Np;
m14 = Phi[nn];
nn = neighborList[n+15*Np]%Np;
m15 = Phi[nn];
nn = neighborList[n+14*Np]%Np;
m16 = Phi[nn];
nn = neighborList[n+17*Np]%Np;
m17 = Phi[nn];
nn = neighborList[n+16*Np]%Np;
m18 = Phi[nn];
nx = -(m1-m2+0.5*(m7-m8+m9-m10+m11-m12+m13-m14));
ny = -(m3-m4+0.5*(m7-m8-m9+m10+m15-m16+m17-m18));
nz = -(m5-m6+0.5*(m11-m12-m13+m14+m15-m16-m17+m18));
//...............................................
//...Store the Color Gradient....................
ColorGrad[n] = nx;
ColorGrad[Np+n] = ny;
ColorGrad[2*Np+n] = nz;
//...............................................
}
}
}
extern "C" void ScaLBL_DFH_Init(double *Phi, double *Den, double *Aq, double *Bq, int start, int finish, int Np){
dvc_ScaLBL_DFH_Init<<<NBLOCKS,NTHREADS >>>(Phi, Den, Aq, Bq, start, finish, Np);
cudaError_t err = cudaGetLastError();
if (cudaSuccess != err){
printf("CUDA error in ScaLBL_DFH_Init: %s \n",cudaGetErrorString(err));
}
}
extern "C" void ScaLBL_D3Q19_AAeven_DFH(int *neighborList, double *dist, double *Aq, double *Bq, double *Den, double *Phi,
double *SolidPotential, double rhoA, double rhoB, double tauA, double tauB, double alpha, double beta,
double Fx, double Fy, double Fz, int start, int finish, int Np){
@@ -1419,3 +1528,14 @@ extern "C" void ScaLBL_D3Q7_AAeven_DFH(double *Aq, double *Bq, double *Den, doub
}
extern "C" void ScaLBL_D3Q19_Gradient_DFH(int *NeighborList, double *Phi, double *ColorGrad, int start, int finish, int Np){
int strideY=Nx;
int strideZ=Nx*Ny;
dvc_ScaLBL_D3Q19_Gradient_DFH<<<NBLOCKS,NTHREADS >>>(NeighborList, Phi, ColorGrad, 0, Np, Np);
cudaError_t err = cudaGetLastError();
if (cudaSuccess != err){
printf("CUDA error in ScaLBL_D3Q19_Gradient_DFH: %s \n",cudaGetErrorString(err));
}
}