trenamed interaction force declaration - 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 bcc75acd3cd0e3e2b9e32510e870232204eb3cc5
(DIR) parent e2bfb9f3fdc7dc794896d1f001a937da9d68cc56
(HTM) Author: Anders Damsgaard <anders.damsgaard@geo.au.dk>
Date: Fri, 6 Jun 2014 14:06:08 +0200
renamed interaction force declaration
Diffstat:
M src/device.cu | 11 ++++++++---
M src/navierstokes.cuh | 13 +++++++++----
M src/sphere.h | 4 ++--
3 files changed, 19 insertions(+), 9 deletions(-)
---
(DIR) diff --git a/src/device.cu b/src/device.cu
t@@ -947,13 +947,18 @@ __host__ void DEM::startTime()
cudaThreadSynchronize();
checkForCudaErrorsIter("Post findInteractionForce", iter);
+ setNSghostNodes<Float><<<dimGridFluid, dimBlockFluid>>>(
+ dev_ns_p, ns.bc_bot, ns.bc_top);
+ cudaThreadSynchronize();
+ checkForCudaErrorsIter("Post setNSghostNodes(dev_ns_p)", iter);
+
// Apply fluid-particle interaction force to the fluid
applyInteractionForceToFluid<<<dimGridFluid, dimBlockFluid>>>(
dev_gridParticleIndex,
dev_cellStart,
dev_cellEnd,
dev_ns_f_pf,
- dev_ns_fi);
+ dev_ns_F_pf);
cudaThreadSynchronize();
checkForCudaErrorsIter("Post applyInteractionForceToFluid",
iter);
t@@ -1047,7 +1052,7 @@ __host__ void DEM::startTime()
dev_ns_v, ns.bc_bot, ns.bc_top);
cudaThreadSynchronize();
checkForCudaErrorsIter("Post setNSghostNodes(v)", iter);
-
+
// Find the divergence of phi*vi*v, needed for predicting the
// fluid velocities
if (PROFILING == 1)
t@@ -1080,7 +1085,7 @@ __host__ void DEM::startTime()
ns.bc_bot,
ns.bc_top,
ns.beta,
- dev_ns_fi,
+ dev_ns_F_pf,
ns.ndem,
dev_ns_v_p_x,
dev_ns_v_p_y,
(DIR) diff --git a/src/navierstokes.cuh b/src/navierstokes.cuh
t@@ -75,8 +75,10 @@ void DEM::initNSmemDev(void)
cudaMalloc((void**)&dev_ns_v_p_z, memSizeFface); // pred. vel. in stag. grid
cudaMalloc((void**)&dev_ns_vp_avg, memSizeF*3); // avg. particle velocity
cudaMalloc((void**)&dev_ns_d_avg, memSizeF); // avg. particle diameter
- //cudaMalloc((void**)&dev_ns_F_pf, memSizeF*3); // interaction force
- cudaMalloc((void**)&dev_ns_F_pf, memSizeFface); // interaction force
+ cudaMalloc((void**)&dev_ns_F_pf, memSizeF*3); // interaction force
+ //cudaMalloc((void**)&dev_ns_F_pf_x, memSizeFface); // interaction force
+ //cudaMalloc((void**)&dev_ns_F_pf_y, memSizeFface); // interaction force
+ //cudaMalloc((void**)&dev_ns_F_pf_z, memSizeFface); // interaction force
cudaMalloc((void**)&dev_ns_phi, memSizeF); // cell porosity
cudaMalloc((void**)&dev_ns_dphi, memSizeF); // cell porosity change
//cudaMalloc((void**)&dev_ns_div_phi_v_v, memSizeF*3); // div(phi v v)
t@@ -87,7 +89,7 @@ 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_tau, memSizeF*6); // stress tensor
//cudaMalloc((void**)&dev_ns_div_tau, memSizeF*3); // div(tau), cell center
cudaMalloc((void**)&dev_ns_div_tau_x, memSizeFface); // div(tau), cell face
cudaMalloc((void**)&dev_ns_div_tau_y, memSizeFface); // div(tau), cell face
t@@ -114,6 +116,9 @@ void DEM::freeNSmemDev()
cudaFree(dev_ns_vp_avg);
cudaFree(dev_ns_d_avg);
cudaFree(dev_ns_F_pf);
+ //cudaFree(dev_ns_F_pf_x);
+ //cudaFree(dev_ns_F_pf_y);
+ //cudaFree(dev_ns_F_pf_z);
cudaFree(dev_ns_phi);
cudaFree(dev_ns_dphi);
//cudaFree(dev_ns_div_phi_v_v);
t@@ -3034,7 +3039,7 @@ __global__ void applyInteractionForceToFluid(
unsigned int* dev_cellStart, // in
unsigned int* dev_cellEnd, // in
Float3* dev_ns_f_pf, // in
- Float3* dev_ns_F_pf) // out
+ Float3* dev_ns_F_pf) // out
{
// 3D thread index
const unsigned int x = blockDim.x * blockIdx.x + threadIdx.x;
(DIR) diff --git a/src/sphere.h b/src/sphere.h
t@@ -181,7 +181,7 @@ class DEM {
Float* dev_ns_v_p_z; // Predicted cell fluid velocity in st. gr.
Float3* dev_ns_vp_avg; // Average particle velocity in cell
Float* dev_ns_d_avg; // Average particle diameter in cell
- Float3* dev_ns_fi; // Particle-fluid interaction force
+ Float3* dev_ns_F_pf; // Interaction force on fluid
Float* dev_ns_phi; // Cell porosity
Float* dev_ns_dphi; // Cell porosity change
//Float3* dev_ns_div_phi_v_v; // Divegence used in velocity prediction
t@@ -200,7 +200,7 @@ class DEM {
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
- Float3* dev_ns_f_pf; // Particle-fluid interaction force
+ Float3* dev_ns_f_pf; // Interaction force on particles
//// Navier Stokes functions