tChanged name of force sum to F instead of N - 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 c3c603cccf305f9764fc3c3179bfdb9976e531b3
(DIR) parent 56390347e098fcaf5b75899e0501283f20352d6f
(HTM) Author: Anders Damsgaard <adc@geo.au.dk>
Date: Mon, 3 Sep 2012 14:21:26 +0200
Changed name of force sum to F instead of N
Diffstat:
M src/contactmodels.cuh | 14 +++++++-------
M src/contactsearch.cuh | 26 +++++++++++++-------------
2 files changed, 20 insertions(+), 20 deletions(-)
---
(DIR) diff --git a/src/contactmodels.cuh b/src/contactmodels.cuh
t@@ -7,7 +7,7 @@
// Linear viscoelastic contact model for particle-wall interactions
// with tangential friction and rolling resistance
-__device__ Float contactLinear_wall(Float3* N, Float3* T, Float* es_dot, Float* p,
+__device__ Float contactLinear_wall(Float3* F, Float3* T, Float* es_dot, Float* p,
unsigned int idx_a, Float radius_a,
Float4* dev_vel_sorted, Float4* dev_angvel_sorted,
Float3 n, Float delta, Float wvel)
t@@ -87,7 +87,7 @@ __device__ Float contactLinear_wall(Float3* N, Float3* T, Float* es_dot, Float*
}*/
// Total force from wall
- *N += f_n + f_t;
+ *F += f_n + f_t;
// Total torque from wall
*T += -radius_a * cross(n, f_t) + T_res;
t@@ -96,14 +96,14 @@ __device__ Float contactLinear_wall(Float3* N, Float3* T, Float* es_dot, Float*
*p += f_n_length / (4.0f * PI * radius_a*radius_a);
// Return force excerted onto the wall
- //return -dot(*N, n);
+ //return -dot(*F, n);
return dot(f_n, n);
}
// Linear vicoelastic contact model for particle-particle interactions
// with tangential friction and rolling resistance
-__device__ void contactLinearViscous(Float3* N, Float3* T, Float* es_dot, Float* p,
+__device__ void contactLinearViscous(Float3* F, Float3* T, Float* es_dot, Float* p,
unsigned int idx_a, unsigned int idx_b,
Float4* dev_vel_sorted,
Float4* dev_angvel_sorted,
t@@ -234,7 +234,7 @@ __device__ void contactLinearViscous(Float3* N, Float3* T, Float* es_dot, Float*
*/
// Add force components from this collision to total force for particle
- *N += f_n + f_t + f_c;
+ *F += f_n + f_t + f_c;
*T += -R_bar * cross(n_ab, f_t) + T_res;
// Pressure excerted onto the particle from this contact
t@@ -244,7 +244,7 @@ __device__ void contactLinearViscous(Float3* N, Float3* T, Float* es_dot, Float*
// Linear elastic contact model for particle-particle interactions
-__device__ void contactLinear(Float3* N, Float3* T,
+__device__ void contactLinear(Float3* F, Float3* T,
Float* es_dot, Float* p,
unsigned int idx_a_orig,
unsigned int idx_b_orig,
t@@ -396,7 +396,7 @@ __device__ void contactLinear(Float3* N, Float3* T,
}
// Add force components from this collision to total force for particle
- *N += f_n + f_t + f_c;
+ *F += f_n + f_t + f_c;
*T += -R_bar * cross(n_ab, f_t) + T_res;
// Pressure excerted onto the particle from this contact
(DIR) diff --git a/src/contactsearch.cuh b/src/contactsearch.cuh
t@@ -408,7 +408,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
Float p = 0.0f;
// Allocate memory for temporal force/torque vector values
- Float3 N = MAKE_FLOAT3(0.0f, 0.0f, 0.0f);
+ Float3 F = MAKE_FLOAT3(0.0f, 0.0f, 0.0f);
Float3 T = MAKE_FLOAT3(0.0f, 0.0f, 0.0f);
// Apply linear elastic, frictional contact model to registered contacts
t@@ -447,7 +447,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
if (delta_n < 0.0f) {
//cuPrintf("\nProcessing contact, idx_a_orig = %u, idx_b_orig = %u, contact = %d, delta_n = %f\n",
// idx_a_orig, idx_b_orig, i, delta_n);
- /*contactLinearViscous(&N, &T, &es_dot, &p,
+ /*contactLinearViscous(&F, &T, &es_dot, &p,
idx_a_orig, idx_b_orig,
dev_vel,
dev_angvel,
t@@ -455,7 +455,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
x_ab, x_ab_length,
delta_n, devC_kappa);
dev_delta_t[mempos] = MAKE_FLOAT4(0.0f, 0.0f, 0.0f, 0.0f);*/
- contactLinear(&N, &T, &es_dot, &p,
+ contactLinear(&F, &T, &es_dot, &p,
idx_a_orig,
idx_b_orig,
vel_a,
t@@ -501,7 +501,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
for (int x_dim=-1; x_dim<2; ++x_dim) { // x-axis
targetPos = gridPos + make_int3(x_dim, y_dim, z_dim);
overlapsInCell(targetPos, idx_a, x_a, radius_a,
- &N, &T, &es_dot, &p,
+ &F, &T, &es_dot, &p,
dev_x_sorted, dev_radius_sorted,
dev_vel_sorted, dev_angvel_sorted,
dev_cellStart, dev_cellEnd,
t@@ -521,7 +521,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
delta_w = w_up_nx.w - (x_a.z + radius_a);
w_n = MAKE_FLOAT3(0.0f, 0.0f, -1.0f);
if (delta_w < 0.0f) {
- w_force = contactLinear_wall(&N, &T, &es_dot, &p, idx_a, radius_a,
+ w_force = contactLinear_wall(&F, &T, &es_dot, &p, idx_a, radius_a,
dev_vel_sorted, dev_angvel_sorted,
w_n, delta_w, w_up_mvfd.y);
}
t@@ -530,7 +530,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
delta_w = x_a.z - radius_a - origo.z;
w_n = MAKE_FLOAT3(0.0f, 0.0f, 1.0f);
if (delta_w < 0.0f) {
- (void)contactLinear_wall(&N, &T, &es_dot, &p, idx_a, radius_a,
+ (void)contactLinear_wall(&F, &T, &es_dot, &p, idx_a, radius_a,
dev_vel_sorted, dev_angvel_sorted,
w_n, delta_w, 0.0f);
}
t@@ -542,7 +542,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
delta_w = x_a.x - radius_a - origo.x;
w_n = MAKE_FLOAT3(1.0f, 0.0f, 0.0f);
if (delta_w < 0.0f) {
- (void)contactLinear_wall(&N, &T, &es_dot, &p, idx_a, radius_a,
+ (void)contactLinear_wall(&F, &T, &es_dot, &p, idx_a, radius_a,
dev_vel_sorted, dev_angvel_sorted,
w_n, delta_w, 0.0f);
}
t@@ -551,7 +551,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
delta_w = L.x - (x_a.x + radius_a);
w_n = MAKE_FLOAT3(-1.0f, 0.0f, 0.0f);
if (delta_w < 0.0f) {
- (void)contactLinear_wall(&N, &T, &es_dot, &p, idx_a, radius_a,
+ (void)contactLinear_wall(&F, &T, &es_dot, &p, idx_a, radius_a,
dev_vel_sorted, dev_angvel_sorted,
w_n, delta_w, 0.0f);
}
t@@ -560,7 +560,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
delta_w = x_a.y - radius_a - origo.y;
w_n = MAKE_FLOAT3(0.0f, 1.0f, 0.0f);
if (delta_w < 0.0f) {
- (void)contactLinear_wall(&N, &T, &es_dot, &p, idx_a, radius_a,
+ (void)contactLinear_wall(&F, &T, &es_dot, &p, idx_a, radius_a,
dev_vel_sorted, dev_angvel_sorted,
w_n, delta_w, 0.0f);
}
t@@ -569,7 +569,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
delta_w = L.y - (x_a.y + radius_a);
w_n = MAKE_FLOAT3(0.0f, -1.0f, 0.0f);
if (delta_w < 0.0f) {
- (void)contactLinear_wall(&N, &T, &es_dot, &p, idx_a, radius_a,
+ (void)contactLinear_wall(&F, &T, &es_dot, &p, idx_a, radius_a,
dev_vel_sorted, dev_angvel_sorted,
w_n, delta_w, 0.0f);
}
t@@ -580,7 +580,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
delta_w = x_a.y - radius_a - origo.y;
w_n = MAKE_FLOAT3(0.0f, 1.0f, 0.0f);
if (delta_w < 0.0f) {
- (void)contactLinear_wall(&N, &T, &es_dot, &p, idx_a, radius_a,
+ (void)contactLinear_wall(&F, &T, &es_dot, &p, idx_a, radius_a,
dev_vel_sorted, dev_angvel_sorted,
w_n, delta_w, 0.0f);
}
t@@ -589,7 +589,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
delta_w = L.y - (x_a.y + radius_a);
w_n = MAKE_FLOAT3(0.0f, -1.0f, 0.0f);
if (delta_w < 0.0f) {
- (void)contactLinear_wall(&N, &T, &es_dot, &p, idx_a, radius_a,
+ (void)contactLinear_wall(&F, &T, &es_dot, &p, idx_a, radius_a,
dev_vel_sorted, dev_angvel_sorted,
w_n, delta_w, 0.0f);
}
t@@ -601,7 +601,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
// Write force to unsorted position
unsigned int orig_idx = dev_gridParticleIndex[idx_a];
- dev_force[orig_idx] = MAKE_FLOAT4(N.x, N.y, N.z, 0.0f);
+ dev_force[orig_idx] = MAKE_FLOAT4(F.x, F.y, F.z, 0.0f);
dev_torque[orig_idx] = MAKE_FLOAT4(T.x, T.y, T.z, 0.0f);
dev_es_dot[orig_idx] = es_dot;
dev_es[orig_idx] += es_dot * devC_dt;