tadded finite checks to interaction forces - 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 396db48750933e4dbbbecf607fb7c9c92173887b
 (DIR) parent d4d46bfe646485d0011d5dd7d342995be9a6b076
 (HTM) Author: Anders Damsgaard <anders.damsgaard@geo.au.dk>
       Date:   Fri, 23 May 2014 10:50:55 +0200
       
       added finite checks to interaction forces
       
       Diffstat:
         M src/debug.h                         |       4 ++--
         M src/navierstokes.cuh                |      10 ++++++++++
       
       2 files changed, 12 insertions(+), 2 deletions(-)
       ---
 (DIR) diff --git a/src/debug.h b/src/debug.h
       t@@ -52,7 +52,7 @@ const int conv_log_interval = 10;
        // simulation of particulate systems: Theoretical developments".
        // SET_2 corresponds approximately to Model A in Zhu et al. 2007.
        // Choose exactly one.
       -//#define SET_1
       -#define SET_2
       +#define SET_1
       +//#define SET_2
        
        #endif
 (DIR) diff --git a/src/navierstokes.cuh b/src/navierstokes.cuh
       t@@ -2851,12 +2851,19 @@ __global__ void findInteractionForce(
                // Interaction force on particle (force) and fluid (f_pf)
                __syncthreads();
                const Float3 f_pf = f_d + f_p + f_v;
       +
       +#ifdef CHECK_NS_FINITE
       +        checkFiniteFloat3("f_pf", i_x, i_y, i_z, f_pf);
       +#endif
       +
        #ifdef SET_1
                dev_ns_f_pf[i] = f_pf;
        #endif
       +
        #ifdef SET_2
                dev_ns_f_pf[i] = f_d;
        #endif
       +
                dev_force[i] += MAKE_FLOAT4(f_pf.x, f_pf.y, f_pf.z, 0.0);
            }
        }
       t@@ -2913,6 +2920,9 @@ __global__ void applyInteractionForceToFluid(
                }
        
                __syncthreads();
       +#ifdef CHECK_NS_FINITE
       +        checkFiniteFloat3("fi", x, y, z, fi/(dx*dy*dz));
       +#endif
                dev_ns_fi[cellidx] = fi/(dx*dy*dz);
            }
        }