tHertzian nonlinear contact model added. Untested. - 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 f9a48cab2eecf0dcce6e2e6ddf1ca837a14a41cc
(DIR) parent 1569cbad46a0701483843304614911364464542a
(HTM) Author: Anders Damsgaard <adc@geo.au.dk>
Date: Tue, 16 Oct 2012 16:26:30 +0200
Hertzian nonlinear contact model added. Untested.
Diffstat:
M src/contactsearch.cuh | 45 ++++++++++++++++++-------------
M src/device.cu | 4 ++--
2 files changed, 28 insertions(+), 21 deletions(-)
---
(DIR) diff --git a/src/contactsearch.cuh b/src/contactsearch.cuh
t@@ -426,7 +426,7 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
Float3 T = MAKE_FLOAT3(0.0f, 0.0f, 0.0f);
// Apply linear elastic, frictional contact model to registered contacts
- if (devC_shearmodel == 2) {
+ if (devC_shearmodel == 2 || devC_shearmodel == 3) {
unsigned int idx_b_orig, mempos;
Float delta_n, x_ab_length, radius_b;
Float3 x_ab;
t@@ -458,24 +458,31 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
// Process collision if the particles are overlapping
if (delta_n < 0.0f) {
- /*contactLinearViscous(&F, &T, &es_dot, &ev_dot, &p,
- idx_a_orig, idx_b_orig,
- dev_vel,
- dev_angvel,
- radius_a, radius_b,
- x_ab, x_ab_length,
- delta_n, devC_kappa);*/
- contactLinear(&F, &T, &es_dot, &ev_dot, &p,
- idx_a_orig,
- idx_b_orig,
- vel_a,
- dev_vel,
- angvel_a,
- dev_angvel,
- radius_a, radius_b,
- x_ab, x_ab_length,
- delta_n, dev_delta_t,
- mempos);
+ if (devC_shearmodel == 2) {
+ contactLinear(&F, &T, &es_dot, &ev_dot, &p,
+ idx_a_orig,
+ idx_b_orig,
+ vel_a,
+ dev_vel,
+ angvel_a,
+ dev_angvel,
+ radius_a, radius_b,
+ x_ab, x_ab_length,
+ delta_n, dev_delta_t,
+ mempos);
+ } else if (devC_shearmodel == 3) {
+ contactHertz(&F, &T, &es_dot, &ev_dot, &p,
+ idx_a_orig,
+ idx_b_orig,
+ vel_a,
+ dev_vel,
+ angvel_a,
+ dev_angvel,
+ radius_a, radius_b,
+ x_ab, x_ab_length,
+ delta_n, dev_delta_t,
+ mempos);
+ }
} else {
__syncthreads();
// Remove this contact (there is no particle with index=np)
(DIR) diff --git a/src/device.cu b/src/device.cu
t@@ -486,8 +486,8 @@ __host__ void gpuMain(Float4* host_x,
// The contact search in topology() is only necessary for determining
// the accumulated shear distance needed in the linear elastic
- // shear force model
- if (params->shearmodel == 2) {
+ // and nonlinear contact force model
+ if (params->shearmodel == 2 || params->shearmodel == 3) {
// For each particle: Search contacts in neighbor cells
if (PROFILING == 1)
startTimer(&kernel_tic);