tpossibly final staggered grid corrections - sphere - GPU-based 3D discrete element method algorithm with optional fluid coupling
(HTM) git clone git://src.adamsgaard.dk/sphere
(DIR) Log
(DIR) Files
(DIR) Refs
(DIR) LICENSE
---
(DIR) commit d009927b87ac9bdd75a225ef9af5a7433930b2bc
(DIR) parent 343ded598fc82cd76bf571514e691b4b59e38cbb
(HTM) Author: Anders Damsgaard <anders.damsgaard@geo.au.dk>
Date: Tue, 3 Jun 2014 16:18:15 +0200
possibly final staggered grid corrections
Diffstat:
M src/device.cu | 12 +++++++++++-
M src/navierstokes.cuh | 30 +++++++++++++-----------------
M src/sphere.h | 4 ++--
3 files changed, 26 insertions(+), 20 deletions(-)
---
(DIR) diff --git a/src/device.cu b/src/device.cu
t@@ -1079,7 +1079,16 @@ __host__ void DEM::startTime()
ns.bc_bot, ns.bc_top);
cudaThreadSynchronize();
checkForCudaErrorsIter(
- "Post setNSghostNodesFaceFloat3(dev_ns_v_p)", iter);
+ "Post setNSghostNodesFace(dev_ns_v_p)", iter);
+
+ interpolateFaceToCenter<<<dimGridFluid, dimBlockFluid>>>(
+ dev_ns_v_p_x,
+ dev_ns_v_p_y,
+ dev_ns_v_p_z,
+ dev_ns_v_p);
+ cudaThreadSynchronize();
+ checkForCudaErrorsIter("Post interpolateFaceToCenter", iter);
+
// In the first iteration of the sphere program, we'll need to
// manually estimate the values of epsilon. In the subsequent
t@@ -1195,6 +1204,7 @@ __host__ void DEM::startTime()
dev_ns_f,
dev_ns_phi,
dev_ns_dphi,
+ dev_ns_v_p,
dev_ns_v_x,
dev_ns_v_y,
dev_ns_v_z,
(DIR) diff --git a/src/navierstokes.cuh b/src/navierstokes.cuh
t@@ -86,8 +86,8 @@ void DEM::initNSmemDev(void)
cudaMalloc((void**)&dev_ns_f, memSizeF); // forcing function value
cudaMalloc((void**)&dev_ns_f1, memSizeF); // constant addition in forcing
cudaMalloc((void**)&dev_ns_f2, memSizeF*3); // constant slope in forcing
- cudaMalloc((void**)&dev_ns_tau, memSizeF*6); // stress tensor (symmetrical)
- cudaMalloc((void**)&dev_ns_div_tau, memSizeF3); // div(tau), cell center
+ //cudaMalloc((void**)&dev_ns_tau, memSizeF*6); // stress tensor (symmetrical)
+ //cudaMalloc((void**)&dev_ns_div_tau, memSizeF*3); // div(tau), cell center
cudaMalloc((void**)&dev_ns_div_tau_x, memSizeFvel); // div(tau), cell face
cudaMalloc((void**)&dev_ns_div_tau_y, memSizeFvel); // div(tau), cell face
cudaMalloc((void**)&dev_ns_div_tau_z, memSizeFvel); // div(tau), cell face
t@@ -123,7 +123,7 @@ void DEM::freeNSmemDev()
cudaFree(dev_ns_f);
cudaFree(dev_ns_f1);
cudaFree(dev_ns_f2);
- cudaFree(dev_ns_tau);
+ //cudaFree(dev_ns_tau);
cudaFree(dev_ns_div_phi_vi_v);
//cudaFree(dev_ns_div_phi_tau);
//cudaFree(dev_ns_div_tau);
t@@ -1408,18 +1408,18 @@ __device__ Float divergence(
{
// Read 6 cell-face values
__syncthreads();
- const Float3 xn = dev_vectorfield[idx(x,y,z)];
- const Float3 xp = dev_vectorfield[idx(x+1,y,z)];
- const Float3 yn = dev_vectorfield[idx(x,y,z)];
- const Float3 yp = dev_vectorfield[idx(x,y+1,z)];
- const Float3 zn = dev_vectorfield[idx(x,y,z)];
- const Float3 zp = dev_vectorfield[idx(x,y,z+1)];
+ const Float xn = dev_vectorfield_x[vidx(x,y,z)];
+ const Float xp = dev_vectorfield_x[vidx(x+1,y,z)];
+ const Float yn = dev_vectorfield_y[vidx(x,y,z)];
+ const Float yp = dev_vectorfield_y[vidx(x,y+1,z)];
+ const Float zn = dev_vectorfield_z[vidx(x,y,z)];
+ const Float zp = dev_vectorfield_z[vidx(x,y,z+1)];
// Calculate the central difference gradrients and the divergence
return
- (xp.x - xn.x)/dx +
- (yp.y - yn.y)/dy +
- (zp.z - zn.z)/dz;
+ (xp - xn)/dx +
+ (yp - yn)/dy +
+ (zp - zn)/dz;
}
// Find the divergence of a tensor field
t@@ -2109,11 +2109,6 @@ __global__ void findPredNSvelocities(
const Float dphi_yn = dev_ns_dphi[idx(x,y-1,z)];
const Float dphi_zn = dev_ns_dphi[idx(x,y,z-1)];
- const Float3 v_xn = dev_ns_v[idx(x-1,y,z)];
- const Float3 v_c = dev_ns_v[cellidx];
- const Float3 v_yn = dev_ns_v[idx(x,y-1,z)];
- const Float3 v_zn = dev_ns_v[idx(x,y,z-1)];
-
const Float3 div_phi_vi_v_xn = dev_ns_div_phi_vi_v[idx(x-1,y,z)];
const Float3 div_phi_vi_v_c = dev_ns_div_phi_vi_v[cellidx];
const Float3 div_phi_vi_v_yn = dev_ns_div_phi_vi_v[idx(x,y-1,z)];
t@@ -2253,6 +2248,7 @@ __global__ void findNSforcing(
Float* dev_ns_f,
Float* dev_ns_phi,
Float* dev_ns_dphi,
+ Float3* dev_ns_v_p,
Float* dev_ns_v_x,
Float* dev_ns_v_y,
Float* dev_ns_v_z,
(DIR) diff --git a/src/sphere.h b/src/sphere.h
t@@ -175,7 +175,7 @@ class DEM {
Float* dev_ns_v_x; // Cell fluid velocity in staggered grid
Float* dev_ns_v_y; // Cell fluid velocity in staggered grid
Float* dev_ns_v_z; // Cell fluid velocity in staggered grid
- //Float3* dev_ns_v_p; // Predicted cell fluid velocity
+ Float3* dev_ns_v_p; // Averaged predicted cell fluid velocity
Float* dev_ns_v_p_x; // Predicted cell fluid velocity in st. gr.
Float* dev_ns_v_p_y; // Predicted cell fluid velocity in st. gr.
Float* dev_ns_v_p_z; // Predicted cell fluid velocity in st. gr.
t@@ -196,7 +196,7 @@ class DEM {
//Float* dev_ns_tau; // Fluid stress tensor
Float3* dev_ns_div_phi_vi_v; // div(phi*vi*v)
//Float3* dev_ns_div_phi_tau; // div(phi*tau)
- Float3* dev_ns_div_tau; // div(tau)
+ //Float3* dev_ns_div_tau; // div(tau)
Float* dev_ns_div_tau_x; // div(tau) on x-face
Float* dev_ns_div_tau_y; // div(tau) on y-face
Float* dev_ns_div_tau_z; // div(tau) on z-face