tuse local neighbor cell center coordinate - 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 59a0f423db8de79044bde53634871bcb97a266cc
 (DIR) parent f1d944982c3295e7505fcc008917f6863aa11ffc
 (HTM) Author: Anders Damsgaard <anders.damsgaard@geo.au.dk>
       Date:   Fri, 20 Mar 2015 15:34:05 +0100
       
       use local neighbor cell center coordinate
       
       Diffstat:
         M src/darcy.cuh                       |      23 ++++++++++++++---------
         M tests/fluid_particle_interaction_d… |      11 +++++++----
       
       2 files changed, 21 insertions(+), 13 deletions(-)
       ---
 (DIR) diff --git a/src/darcy.cuh b/src/darcy.cuh
       t@@ -378,6 +378,8 @@ __global__ void findDarcyPorositiesLinear(
                    // The cell 3d index
                    const int3 gridPos = make_int3((int)x,(int)y,(int)z);
        
       +            Float3 X_n;
       +
                    // The neighbor cell 3d index
                    int3 targetCell;
        
       t@@ -412,6 +414,9 @@ __global__ void findDarcyPorositiesLinear(
                                    // Make sure cell is not empty
                                    if (startIdx != 0xffffffff) {
        
       +                                X_n = X
       +                                    + MAKE_FLOAT3(x_dim*dx, y_dim*dy, z_dim*dz);
       +
                                        // Highest particle index in cell
                                        __syncthreads();
                                        endIdx = dev_cellEnd[cellID];
       t@@ -442,34 +447,34 @@ __global__ void findDarcyPorositiesLinear(
                                            // nodes of component-wise velocity
                                            x3 += distmod;
                                            s = weight(x3,
       -                                            X + MAKE_FLOAT3(-0.5*dx, 0.0, 0.0),
       +                                            X_n + MAKE_FLOAT3(-0.5*dx, 0., 0.),
                                                    dx, dy, dz);
                                            v_p_xn += s*vol_p*v3.x/(s*vol_p + 1.0e-16);
        
                                            s = weight(x3,
       -                                            X + MAKE_FLOAT3( 0.5*dx, 0.0, 0.0),
       +                                            X_n + MAKE_FLOAT3( 0.5*dx, 0., 0.),
                                                    dx, dy, dz);
       -                                    v_p_xp += s*vol_p*v3.x/(s*vol_p - 1.0e-16);
       +                                    v_p_xp += s*vol_p*v3.x/(s*vol_p + 1.0e-16);
        
                                            s = weight(x3,
       -                                            X + MAKE_FLOAT3( 0.0, -0.5*dy, 0.0),
       +                                            X_n + MAKE_FLOAT3( 0., -0.5*dy, 0.),
                                                    dx, dy, dz);
                                            v_p_yn += s*vol_p*v3.y/(s*vol_p + 1.0e-16);
        
                                            s = weight(x3,
       -                                            X + MAKE_FLOAT3( 0.0, 0.5*dy, 0.0),
       +                                            X_n + MAKE_FLOAT3( 0., 0.5*dy, 0.),
                                                    dx, dy, dz);
       -                                    v_p_yp += s*vol_p*v3.y/(s*vol_p - 1.0e-16);
       +                                    v_p_yp += s*vol_p*v3.y/(s*vol_p + 1.0e-16);
        
                                            s = weight(x3,
       -                                            X + MAKE_FLOAT3( 0.0, 0.0, -0.5*dz),
       +                                            X_n + MAKE_FLOAT3( 0., 0., -0.5*dz),
                                                    dx, dy, dz);
                                            v_p_zn += s*vol_p*v3.z/(s*vol_p + 1.0e-16);
        
                                            s = weight(x3,
       -                                            X + MAKE_FLOAT3( 0.0, 0.0, 0.5*dz),
       +                                            X_n + MAKE_FLOAT3( 0., 0., 0.5*dz),
                                                    dx, dy, dz);
       -                                    v_p_zp += s*vol_p*v3.z/(s*vol_p - 1.0e-16);
       +                                    v_p_zp += s*vol_p*v3.z/(s*vol_p + 1.0e-16);
                                        }
                                    }
                                }
 (DIR) diff --git a/tests/fluid_particle_interaction_darcy.py b/tests/fluid_particle_interaction_darcy.py
       t@@ -12,11 +12,12 @@ sim.initFluid(cfd_solver = 1)
        # No gravity, pressure gradient enforced by Dirichlet boundaries.
        # The particle should be sucked towards the low pressure
        print('# Test 1: Test pressure gradient force')
       -sim.p_f[:,:,-1] = 1.0
       +#sim.p_f[:,:,-1] = 1.0
        sim.addParticle([0.5, 0.5, 0.5], 0.01)
       +sim.vel[0,1] = 0.01
        sim.initTemporal(total=0.001, file_dt=0.0001)
       -#sim.time_file_dt[0] = sim.time_dt[0]
       -#sim.time_total[0] = sim.time_dt[0]
       +sim.time_file_dt[0] = sim.time_dt[0]
       +sim.time_total[0] = sim.time_dt[0]
        
        #sim.g[2] = -10.
        sim.run(verbose=False)
       t@@ -28,9 +29,10 @@ sim.run(verbose=False)
        sim.readlast(verbose=False)
        test(sim.vel[0,2] < 0.0, 'Particle velocity:')
        
       -sim.cleanup()
       +#sim.cleanup()
        
        
       +'''
        # Gravity, pressure gradient enforced by Dirichlet boundaries.
        # The particle should be sucked towards the low pressure
        print('# Test 2: Test pressure gradient force from buoyancy')
       t@@ -55,3 +57,4 @@ sim.readlast(verbose=False)
        test(sim.vel[0,2] < 0.0, 'Particle velocity:')
        
        sim.cleanup()
       +'''