tRemoved bug that caused unspec. launch failure when nw=0 - 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 838ab5a7172daca802dfb5f095f7969b0dc546b6
 (DIR) parent 6570fb41207713393b1bb49301116c3e2c20dfd7
 (HTM) Author: Anders Damsgaard <adc@geo.au.dk>
       Date:   Tue, 11 Dec 2012 09:38:26 +0100
       
       Removed bug that caused unspec. launch failure when nw=0
       
       Diffstat:
         M src/contactsearch.cuh               |      29 +++++++++++++++++------------
         M src/sphere.cpp                      |      39 ++++++++++++++++++++++++++++---
       
       2 files changed, 53 insertions(+), 15 deletions(-)
       ---
 (DIR) diff --git a/src/contactsearch.cuh b/src/contactsearch.cuh
       t@@ -451,20 +451,25 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
                        mempos = (unsigned int)(idx_a_orig * devC_nc + i);
                        __syncthreads();
                        idx_b_orig = dev_contacts[mempos];
       -                distmod    = dev_distmod[mempos];
       -                x_b        = dev_x[idx_b_orig];
       -                radius_b   = distmod.w;
       +                //radius_b   = distmod.w;
        
       -                // Inter-particle vector, corrected for periodic boundaries
       -                x_ab = MAKE_FLOAT3(x_a.x - x_b.x + distmod.x,
       -                        x_a.y - x_b.y + distmod.y,
       -                        x_a.z - x_b.z + distmod.z);
        
       -                x_ab_length = length(x_ab);
       -                delta_n = x_ab_length - (radius_a + radius_b);
       +                if (idx_b_orig != (unsigned int)devC_np) {
        
       +                    // Read inter-particle distance correction vector
       +                    distmod    = dev_distmod[mempos];
        
       -                if (idx_b_orig != (unsigned int)devC_np) {
       +                    // Read particle b position and radius
       +                    x_b = dev_x[idx_b_orig];
       +                    radius_b = x_b.w;
       +
       +                    // Inter-particle vector, corrected for periodic boundaries
       +                    x_ab = MAKE_FLOAT3(x_a.x - x_b.x + distmod.x,
       +                            x_a.y - x_b.y + distmod.y,
       +                            x_a.z - x_b.z + distmod.z);
       +
       +                    x_ab_length = length(x_ab);
       +                    delta_n = x_ab_length - (radius_a + radius_b);
        
                            // Process collision if the particles are overlapping
                            if (delta_n < 0.0f) {
       t@@ -553,7 +558,6 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
                    w_force = contactLinear_wall(&F, &T, &es_dot, &ev_dot, &p, idx_a, radius_a,
                            dev_vel_sorted, dev_angvel_sorted,
                            w_n, delta_w, w_up_mvfd.y);
       -            //cuPrintf("particle %d collides with upper wall, wforce = %f\n", idx_a, w_force);
                }
        
                // Lower wall (force on wall not stored)
       t@@ -638,7 +642,8 @@ __global__ void interact(unsigned int* dev_gridParticleIndex, // Input: Unsorted
                dev_es[orig_idx]     += es_dot * devC_dt;
                dev_ev[orig_idx]     += ev_dot * devC_dt;
                dev_p[orig_idx]       = p;
       -        dev_walls_force_pp[orig_idx] = w_force;
       +        if (devC_nw > 0)
       +            dev_walls_force_pp[orig_idx] = w_force;
            }
        } // End of interact(...)
        
 (DIR) diff --git a/src/sphere.cpp b/src/sphere.cpp
       t@@ -146,13 +146,46 @@ void DEM::checkValues(void)
                exit(1);
            }
        
       -    // Check that radii are positive values
       +
       +    // Per-particle checks
       +    Float4 x;
            for (i = 0; i < np; ++i) {
       -        if (k.x[i].w <= 0.0) {
       +
       +        // Read value into register
       +        x = k.x[i];
       +
       +        // Check that radii are positive values
       +        if (x.w <= 0.0) {
                    cerr << "Error: Particle " << i << " has a radius of "
       -                << k.x[i].w << " m.";
       +                << k.x[i].w << " m." << std::endl;
                    exit(1);
                }
       +
       +        // Check that all particles are inside of the grid
       +        if (x.x < grid.origo[0] ||
       +                x.y < grid.origo[1] ||
       +                x.z < grid.origo[2] ||
       +                x.x > grid.L[0] ||
       +                x.y > grid.L[1] ||
       +                x.z > grid.L[2]) {
       +            cerr << "Error: Particle " << i << " is outside of "
       +                << "the computational grid\n"
       +                << "k.x[i] = ["
       +                << x.x << ", "
       +                << x.y << ", "
       +                << x.z << "]\n"
       +                << "grid.origo = ["
       +                << grid.origo[0] << ", "
       +                << grid.origo[1] << ", "
       +                << grid.origo[2] << "], "
       +                << "grid.L = ["
       +                << grid.L[0] << ", "
       +                << grid.L[1] << ", "
       +                << grid.L[2] << "]."
       +                << std::endl;
       +            exit(1);
       +
       +        }
            }
        
            // If present, check that the upper wall is above all particles