ttangential viscous energy dissipated no matter if it slips or not - 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 ce8b1910d25021bdb2fd0795668b3b216ab111e9
 (DIR) parent 9d713d88e6e730aba71a8d153120d28f6584ee40
 (HTM) Author: Anders Damsgaard <anders.damsgaard@geo.au.dk>
       Date:   Fri, 20 Jun 2014 13:28:34 +0200
       
       ttangential viscous energy dissipated no matter if it slips or not
       
       Diffstat:
         M src/contactmodels.cuh               |       6 +++---
         M tests/contactmodel.py               |      27 +++++++++++++++------------
       
       2 files changed, 18 insertions(+), 15 deletions(-)
       ---
 (DIR) diff --git a/src/contactmodels.cuh b/src/contactmodels.cuh
       t@@ -414,6 +414,9 @@ __device__ void contactLinear(Float3* F, Float3* T,
                // Add tangential displacement to total tangential displacement
                delta_t = delta_t0 + vel_t * devC_dt;
        
       +        // Store energy dissipated in tangential viscous component
       +        *ev_dot += 0.5 * devC_params.gamma_t * vel_t_length * vel_t_length;
       +
                // If failure criterion is not met, contact is viscous-linear elastic.
                // If failure criterion is met, contact force is limited, 
                // resulting in a slip and energy dissipation
       t@@ -448,9 +451,6 @@ __device__ void contactLinear(Float3* F, Float3* T,
                        * (devC_params.mu_d * length(f_n-f_c) * t
                                + devC_params.gamma_t * vel_t);
        
       -            // Store energy dissipated in tangential viscous component
       -            *ev_dot += 0.5 * devC_params.gamma_t * vel_t_length * vel_t_length;
       -
                    // Shear friction heat production rate: 
                    // The energy lost from the tangential spring is dissipated as heat
                    //*es_dot += -dot(vel_t_ab, f_t);
 (DIR) diff --git a/tests/contactmodel.py b/tests/contactmodel.py
       t@@ -154,11 +154,11 @@ orig.initTemporal(total = 0.1, file_dt = 0.01)
        orig.run(verbose=False)
        after.readlast(verbose=False)
        pytestutils.test((after.angvel[:,1] < 0.0).all(),
       -                 "Oblique normal collision (1/5):")
       +                 "Oblique normal collision (1/8):")
        pytestutils.compareFloats(orig.totalKineticEnergy(),
                                  after.totalKineticEnergy()
                                  + after.totalRotationalEnergy(),
       -                          "Oblique normal collision (2/5):", tolerance=0.05)
       +                          "Oblique normal collision (2/8):", tolerance=0.05)
        
        # Normal impact, low angle, slip
        orig = sphere.sim(np=2, sid='contactmodeltest')
       t@@ -181,13 +181,13 @@ pytestutils.compareFloats(orig.totalKineticEnergy(),
                                  after.totalKineticEnergy()
                                  + after.totalRotationalEnergy()
                                  + after.totalFrictionalEnergy(),
       -                          "Oblique normal collision (3/5):", tolerance=0.05)
       +                          "Oblique normal collision (3/8):", tolerance=0.05)
        pytestutils.test((after.angvel[:,1] < 0.0).all(),
       -                 "Oblique normal collision (4/5):")
       +                 "Oblique normal collision (4/8):")
        pytestutils.test(after.totalFrictionalEnergy() > 0.0,
       -                 "Oblique normal collision (5/5):")
       +                 "Oblique normal collision (5/8):")
        
       -# Normal impact, low angle, slip
       +# Normal impact, low angle, slip, viscous damping tangentially
        orig = sphere.sim(np=2, sid='contactmodeltest')
        after = sphere.sim(np=2, sid='contactmodeltest')
        sphere.cleanup(orig)
       t@@ -199,19 +199,22 @@ orig.vel[0,2] = 1
        orig.vel[0,0] = 1
        orig.mu_s[0] = 0.3
        orig.mu_d[0] = 0.3
       -orig.gamma_t[0] = 1.0e6
       +orig.gamma_t[0] = 1.0e3
        orig.defineWorldBoundaries(L=[10,10,10])
        orig.initTemporal(total = 0.1, file_dt = 0.01)
        
        orig.run(verbose=False)
        after.readlast(verbose=False)
       -after.totalViscousEnergy()
       +print(after.totalViscousEnergy())
        pytestutils.compareFloats(orig.totalKineticEnergy(),
                                  after.totalKineticEnergy()
                                  + after.totalRotationalEnergy()
       -                          + after.totalFrictionalEnergy(),
       -                          "Oblique normal collision (6/5):", tolerance=0.05)
       +                          + after.totalFrictionalEnergy()
       +                          + after.totalViscousEnergy(),
       +                          "Oblique normal collision (6/8):", tolerance=0.05)
        pytestutils.test((after.angvel[:,1] < 0.0).all(),
       -                 "Oblique normal collision (7/5):")
       +                 "Oblique normal collision (7/8):")
       +pytestutils.test(after.totalFrictionalEnergy() > 0.0,
       +                 "Oblique normal collision (8/8):")
        pytestutils.test(after.totalFrictionalEnergy() > 0.0,
       -                 "Oblique normal collision (8/5):")
       +                 "Oblique normal collision (8/8):")