This commit is contained in:
James E McClure 2015-07-04 16:34:39 -04:00
parent 9f2840ccdc
commit 5299fd64ea
2 changed files with 143 additions and 134 deletions

View File

@ -1110,87 +1110,98 @@ __global__ void dvc_MassColorCollideD3Q7(char *ID, double *A_even, double *A_od
id = ID[n];
if ( id != 0){
//.....Load the Color gradient.........
nx = ColorGrad[n];
ny = ColorGrad[N+n];
nz = ColorGrad[2*N+n];
C = sqrt(nx*nx+ny*ny+nz*nz);
if (C == 0.0) C=1.0;
nx = nx/C;
ny = ny/C;
nz = nz/C;
//....Load the flow velocity...........
ux = Velocity[n];
uy = Velocity[N+n];
uz = Velocity[2*N+n];
//........................................................................
// READ THE DISTRIBUTIONS
// (read from opposite array due to previous swap operation)
//........................................................................
f2 = A_odd[n];
f4 = A_odd[N+n];
f6 = A_odd[2*N+n];
f0 = A_even[n];
f1 = A_even[N+n];
f3 = A_even[2*N+n];
f5 = A_even[3*N+n];
na = f0+f1+f2+f3+f4+f5+f6;
//........................................................................
f2 = B_odd[n];
f4 = B_odd[N+n];
f6 = B_odd[2*N+n];
f0 = B_even[n];
f1 = B_even[N+n];
f3 = B_even[2*N+n];
f5 = B_even[3*N+n];
nb = f0+f1+f2+f3+f4+f5+f6;
nab = 1.0/(na+nb);
//........................................................................
delta = beta*na*nb*nab*0.1111111111111111*nx;
if (na*nb*nab<0.0)) delta=0.0;
a1 = na*(0.1111111111111111*(1+4.5*ux))+delta;
b1 = nb*(0.1111111111111111*(1+4.5*ux))-delta;
a2 = na*(0.1111111111111111*(1-4.5*ux))-delta;
b2 = nb*(0.1111111111111111*(1-4.5*ux))+delta;
//.....Load the Color gradient.........
nx = ColorGrad[n];
ny = ColorGrad[N+n];
nz = ColorGrad[2*N+n];
C = sqrt(nx*nx+ny*ny+nz*nz);
if (C==0.0) C=1.0;
nx = nx/C;
ny = ny/C;
nz = nz/C;
//....Load the flow velocity...........
ux = Velocity[n];
uy = Velocity[N+n];
uz = Velocity[2*N+n];
//........................................................................
// READ THE DISTRIBUTIONS
// (read from opposite array due to previous swap operation)
//........................................................................
f2 = A_odd[n];
f4 = A_odd[N+n];
f6 = A_odd[2*N+n];
f0 = A_even[n];
f1 = A_even[N+n];
f3 = A_even[2*N+n];
f5 = A_even[3*N+n];
na = f0+f1+f2+f3+f4+f5+f6;
//........................................................................
f2 = B_odd[n];
f4 = B_odd[N+n];
f6 = B_odd[2*N+n];
f0 = B_even[n];
f1 = B_even[N+n];
f3 = B_even[2*N+n];
f5 = B_even[3*N+n];
nb = f0+f1+f2+f3+f4+f5+f6;
nab = 1.0/(na+nb);
//........................................................................
//....Instantiate the density distributions
// Generate Equilibrium Distributions and stream
// Stationary value - distribution 0
A_even[n] = 0.3333333333333333*na;
B_even[n] = 0.3333333333333333*nb;
// Non-Stationary equilibrium distributions
//feq[0] = 0.1111111111111111*(1+4.5*ux);
//feq[1] = 0.1111111111111111*(1-4.5*ux);
//feq[2] = 0.1111111111111111*(1+4.5*uy);
//feq[3] = 0.1111111111111111*(1-4.5*uy);
//feq[4] = 0.1111111111111111*(1+4.5*uz);
//feq[5] = 0.1111111111111111*(1-4.5*uz);
A_odd[n] = a1;
A_even[N+n] = a2;
B_odd[n] = b1;
B_even[N+n] = b2;
//........................................................................
// q = 2
// Cq = {0,1,0}
delta = beta*na*nb*nab*0.1111111111111111*ny;
if (na*nb*nab<0.0)) delta=0.0;
//...............................................
// q = 0,2,4
// Cq = {1,0,0}, {0,1,0}, {0,0,1}
delta = beta*na*nb*nab*0.1111111111111111*nx;
if (!(na*nb*nab>0)) delta=0;
a1 = na*(0.1111111111111111*(1+4.5*ux))+delta;
b1 = nb*(0.1111111111111111*(1+4.5*ux))-delta;
a2 = na*(0.1111111111111111*(1-4.5*ux))-delta;
b2 = nb*(0.1111111111111111*(1-4.5*ux))+delta;
a1 = na*(0.1111111111111111*(1+4.5*uy))+delta;
b1 = nb*(0.1111111111111111*(1+4.5*uy))-delta;
a2 = na*(0.1111111111111111*(1-4.5*uy))-delta;
b2 = nb*(0.1111111111111111*(1-4.5*uy))+delta;
A_odd[n] = a1;
A_even[N+n] = a2;
B_odd[n] = b1;
B_even[N+n] = b2;
//...............................................
// q = 2
// Cq = {0,1,0}
delta = beta*na*nb*nab*0.1111111111111111*ny;
if (!(na*nb*nab>0)) delta=0;
a1 = na*(0.1111111111111111*(1+4.5*uy))+delta;
b1 = nb*(0.1111111111111111*(1+4.5*uy))-delta;
a2 = na*(0.1111111111111111*(1-4.5*uy))-delta;
b2 = nb*(0.1111111111111111*(1-4.5*uy))+delta;
A_odd[N+n] = a1;
A_even[2*N+n] = a2;
B_odd[N+n] = b1;
B_even[2*N+n] = b2;
//........................................................................
// q = 4
// Cq = {0,0,1}
delta = beta*na*nb*nab*0.1111111111111111*nz;
if (na*nb*nab<0.0)) delta=0.0;
A_odd[N+n] = a1;
A_even[2*N+n] = a2;
B_odd[N+n] = b1;
B_even[2*N+n] = b2;
//...............................................
// q = 4
// Cq = {0,0,1}
delta = beta*na*nb*nab*0.1111111111111111*nz;
if (!(na*nb*nab>0)) delta=0;
a1 = na*(0.1111111111111111*(1+4.5*uz))+delta;
b1 = nb*(0.1111111111111111*(1+4.5*uz))-delta;
a2 = na*(0.1111111111111111*(1-4.5*uz))-delta;
b2 = nb*(0.1111111111111111*(1-4.5*uz))+delta;
A_odd[2*N+n] = a1;
A_even[3*N+n] = a2;
B_odd[2*N+n] = b1;
B_even[3*N+n] = b2;
a1 = na*(0.1111111111111111*(1+4.5*uz))+delta;
b1 = nb*(0.1111111111111111*(1+4.5*uz))-delta;
a2 = na*(0.1111111111111111*(1-4.5*uz))-delta;
b2 = nb*(0.1111111111111111*(1-4.5*uz))+delta;
A_odd[2*N+n] = a1;
A_even[3*N+n] = a2;
B_odd[2*N+n] = b1;
B_even[3*N+n] = b2;
//........................................................................
}
}
}

View File

@ -109,57 +109,55 @@ __global__ void dvc_SwapD3Q7(char *ID, double *disteven, double *distodd, int N
if (n<N ){
id = ID[n];
if (id > 0){
//.......Back out the 3-D indices for node n..............
k = n/(Nx*Ny);
j = (n-Nx*Ny*k)/Nx;
i = n-Nx*Ny*k-Nz*j;
//........................................................................
// Retrieve even distributions from the local node (swap convention)
// f0 = disteven[n]; // Does not particupate in streaming
f1 = distodd[n];
f3 = distodd[N+n];
f5 = distodd[2*N+n];
//........................................................................
//.......Back out the 3-D indices for node n..............
k = n/(Nx*Ny);
j = (n-Nx*Ny*k)/Nx;
i = n-Nx*Ny*k-Nz*j;
//........................................................................
// Retrieve even distributions from the local node (swap convention)
// f0 = disteven[n]; // Does not particupate in streaming
f1 = distodd[n];
f3 = distodd[N+n];
f5 = distodd[2*N+n];
//........................................................................
//........................................................................
// Retrieve odd distributions from neighboring nodes (swap convention)
//........................................................................
nn = n+1; // neighbor index (pull convention)
if (!(i+1<Nx)) nn -= Nx; // periodic BC along the x-boundary
//if (i+1<Nx){
f2 = disteven[N+nn]; // pull neighbor for distribution 2
if (!(f2 < 0.0)){
distodd[n] = f2;
disteven[N+nn] = f1;
}
//}
//........................................................................
nn = n+Nx; // neighbor index (pull convention)
if (!(j+1<Ny)) nn -= Nx*Ny; // Perioidic BC along the y-boundary
//if (j+1<Ny){
f4 = disteven[2*N+nn]; // pull neighbor for distribution 4
if (!(f4 < 0.0)){
distodd[N+n] = f4;
disteven[2*N+nn] = f3;
// }
}
//........................................................................
nn = n+Nx*Ny; // neighbor index (pull convention)
if (!(k+1<Nz)) nn -= Nx*Ny*Nz; // Perioidic BC along the z-boundary
//if (k+1<Nz){
f6 = disteven[3*N+nn]; // pull neighbor for distribution 6
if (!(f6 < 0.0)){
distodd[2*N+n] = f6;
disteven[3*N+nn] = f5;
// }
}
//........................................................................
// Retrieve odd distributions from neighboring nodes (swap convention)
//........................................................................
nn = n+1; // neighbor index (pull convention)
if (!(i+1<Nx)) nn -= Nx; // periodic BC along the x-boundary
//if (i+1<Nx){
f2 = disteven[N+nn]; // pull neighbor for distribution 2
if (!(f2 < 0.0)){
distodd[n] = f2;
disteven[N+nn] = f1;
}
//}
//........................................................................
nn = n+Nx; // neighbor index (pull convention)
if (!(j+1<Ny)) nn -= Nx*Ny; // Perioidic BC along the y-boundary
//if (j+1<Ny){
f4 = disteven[2*N+nn]; // pull neighbor for distribution 4
if (!(f4 < 0.0)){
distodd[N+n] = f4;
disteven[2*N+nn] = f3;
}
//........................................................................
nn = n+Nx*Ny; // neighbor index (pull convention)
if (!(k+1<Nz)) nn -= Nx*Ny*Nz; // Perioidic BC along the z-boundary
//if (k+1<Nz){
f6 = disteven[3*N+nn]; // pull neighbor for distribution 6
if (!(f6 < 0.0)){
distodd[2*N+n] = f6;
disteven[3*N+nn] = f5;
}
}
}
}
}
//*************************************************************************
__global__ void dvc_ComputeDensityD3Q7(char *ID, double *disteven, double *distodd, double *Den,
__global__ void dvc_ComputeDensityD3Q7(char *ID, double *disteven, double *distodd, double *Den,
int Nx, int Ny, int Nz)
{
char id;
@ -172,19 +170,19 @@ __global__ void dvc_ComputeDensityD3Q7(char *ID, double *disteven, double *dist
//........Get 1-D index for this thread....................
n = S*blockIdx.x*blockDim.x + s*blockDim.x + threadIdx.x;
if (n<N){
id = ID[n];
if ( id > 0 ){
// Read the distributions
f0 = disteven[n];
f2 = disteven[N+n];
f4 = disteven[2*N+n];
f6 = disteven[3*N+n];
f1 = distodd[n];
f3 = distodd[N+n];
f5 = distodd[2*N+n];
// Compute the density
Den[n] = f0+f1+f2+f3+f4+f5+f6;
}
id = ID[n];
if (id > 0 ){
// Read the distributions
f0 = disteven[n];
f2 = disteven[N+n];
f4 = disteven[2*N+n];
f6 = disteven[3*N+n];
f1 = distodd[n];
f3 = distodd[N+n];
f5 = distodd[2*N+n];
// Compute the density
Den[n] = f0+f1+f2+f3+f4+f5+f6;
}
}
}
}