tfix minor implementation errors - 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 ea358821f5b6bd7e654e769b5b8e3b844d8a3618
(DIR) parent ab4f66ef7e057437d51be8d867b0a1eb4ffe6f21
(HTM) Author: Anders Damsgaard Christensen <adc@geo.au.dk>
Date: Tue, 16 Aug 2016 11:56:29 -0700
fix minor implementation errors
Diffstat:
M src/contactmodels.cuh | 19 ++++++++++++-------
1 file changed, 12 insertions(+), 7 deletions(-)
---
(DIR) diff --git a/src/contactmodels.cuh b/src/contactmodels.cuh
t@@ -57,10 +57,8 @@ __device__ Float contactLinear_wall(
// 2013) based on Young's modulus if params.E > 0.0.
// Use the calculated stiffness for normal and tangential components.
Float k_n = devC_params.k_n;
- Float k_t = devC_params.k_t;
if (devC_params.E > .001) {
- k_n = 3.141592654/2.0 * devC_params.E * (radius_a + radius_b)/2.;
- k_t = k_n;
+ k_n = 3.141592654/2.0 * devC_params.E * radius_a;
}
// Normal force component: Elastic - viscous damping
t@@ -182,10 +180,8 @@ __device__ void contactLinearViscous(
// 2013) based on Young's modulus if params.E > 0.0.
// Use the calculated stiffness for normal and tangential components.
Float k_n = devC_params.k_n;
- Float k_t = devC_params.k_t;
if (devC_params.E > .001) {
k_n = 3.141592654/2.0 * devC_params.E * (radius_a + radius_b)/2.;
- k_t = k_n;
}
// Normal force component: Elastic - viscous damping
t@@ -323,6 +319,16 @@ __device__ void contactLinear(
// New tangential displacement vector
Float3 delta_t;
+ // Use scale-invariant contact model (Ergenzinger et al 2010, Obermayr et al
+ // 2013) based on Young's modulus if params.E > 0.0.
+ // Use the calculated stiffness for normal and tangential components.
+ Float k_n = devC_params.k_n;
+ Float k_t = devC_params.k_t;
+ if (devC_params.E > .001) {
+ k_n = 3.141592654/2.0 * devC_params.E * (radius_a + radius_b)/2.;
+ k_t = k_n;
+ }
+
// Normal force component: Elastic - viscous damping
f_n = fmax(0.0, -k_n*delta + devC_params.gamma_n * vel_n) * n;
Float f_n_length = length(f_n);
t@@ -558,8 +564,7 @@ __device__ void contactHertz(
// Shear friction heat production rate:
// The energy lost from the tangential spring is dissipated as heat
//*es_dot += -dot(vel_t_ab, f_t);
- *es_dot += length(delta_t0 - delta_t) * k_t / devC_dt; // Seen in
- EsyS-Particle
+ *es_dot += length(delta_t0 - delta_t) * k_t / devC_dt; // Seen in EsyS-Particle
//*es_dot += fabs(dot(delta_t0 - delta_t, f_t)) / devC_dt;
} else { // Static case