tfirst corrections for velocity ghost node writing - 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 9216f53f507b61c0f492f9903acfe0ab70058d37
 (DIR) parent 4b18f7a7031d12e356f96a9784dbf9118802d2f8
 (HTM) Author: Anders Damsgaard <anders.damsgaard@geo.au.dk>
       Date:   Thu,  5 Jun 2014 12:15:53 +0200
       
       first corrections for velocity ghost node writing
       
       Diffstat:
         M src/navierstokes.cpp                |       3 ++-
         M src/navierstokes.cuh                |      68 ++++++++++++++++++++++++-------
       
       2 files changed, 55 insertions(+), 16 deletions(-)
       ---
 (DIR) diff --git a/src/navierstokes.cpp b/src/navierstokes.cpp
       t@@ -97,7 +97,8 @@ unsigned int DEM::vidx(
                const int y,
                const int z)
        {
       -    return x + (ns.nx+1)*y + (ns.nx+1)*(ns.ny+1)*z;
       +    //return x + (ns.nx+1)*y + (ns.nx+1)*(ns.ny+1)*z; // without ghost nodes
       +    return x + (ns.nx+3)*y + (ns.nx+3)*(ns.ny+3)*z; // with ghost nodes
        }
        
        // Determine if the FTCS (forward time, central space) solution of the Navier
 (DIR) diff --git a/src/navierstokes.cuh b/src/navierstokes.cuh
       t@@ -658,38 +658,76 @@ __global__ void setNSghostNodesFace(
            const unsigned int nz = devC_grid.num[2];
        
            // check that we are not outside the fluid grid
       -    if (x <= nx && y <= ny && z <= nz) {
       +    //if (x <= nx && y <= ny && z <= nz) {
       +    if (x < nx && y < ny && z < nz) {
        
                const T val_x = dev_scalarfield_x[vidx(x,y,z)];
                const T val_y = dev_scalarfield_y[vidx(x,y,z)];
                const T val_z = dev_scalarfield_z[vidx(x,y,z)];
        
                // x
       -        if (x == 0)
       +        if (x == 0) {
       +            dev_scalarfield_x[vidx(nx,y,z)] = val_x;
       +            dev_scalarfield_y[vidx(nx,y,z)] = val_y;
       +            dev_scalarfield_z[vidx(nx,y,z)] = val_z;
       +        }
       +        if (x == 1) {
                    dev_scalarfield_x[vidx(nx+1,y,z)] = val_x;
       -        if (x == nx)
       +        }
       +        if (x == nx-1) {
                    dev_scalarfield_x[vidx(-1,y,z)] = val_x;
       +            dev_scalarfield_y[vidx(-1,y,z)] = val_y;
       +            dev_scalarfield_z[vidx(-1,y,z)] = val_z;
       +        }
        
                // y
       -        if (y == 0)
       +        if (y == 0) {
       +            dev_scalarfield_x[vidx(x,ny,z)] = val_x;
       +            dev_scalarfield_y[vidx(x,ny,z)] = val_y;
       +            dev_scalarfield_z[vidx(x,ny,z)] = val_z;
       +        }
       +        if (y == 1) {
                    dev_scalarfield_y[vidx(x,ny+1,z)] = val_y;
       -        if (y == ny)
       +        }
       +        if (y == ny-1) {
       +            dev_scalarfield_x[vidx(x,-1,z)] = val_x;
                    dev_scalarfield_y[vidx(x,-1,z)] = val_y;
       +            dev_scalarfield_z[vidx(x,-1,z)] = val_z;
       +        }
        
                // z
       -        if (z == 0 && bc_bot == 0)
       -            dev_scalarfield_z[vidx(x,y,-1)] = val_z;     // Dirichlet
       -        if (z == 1 && bc_bot == 1)
       -            dev_scalarfield_z[vidx(x,y,-1)] = val_z;     // Neumann
       -        if (z == 0 && bc_bot == 2)
       +        if (z == 0 && bc_bot == 0) {
       +            dev_scalarfield_x[vidx(x,y,-1)] = val_y;     // Dirichlet -z
       +            dev_scalarfield_y[vidx(x,y,-1)] = val_x;     // Dirichlet -z
       +            dev_scalarfield_z[vidx(x,y,-1)] = val_z;     // Dirichlet -z
       +        }
       +        if (z == 1 && bc_bot == 1) {
       +            dev_scalarfield_x[vidx(x,y,-1)] = val_x;     // Neumann -z
       +            dev_scalarfield_y[vidx(x,y,-1)] = val_y;     // Neumann -z
       +            dev_scalarfield_z[vidx(x,y,-1)] = val_z;     // Neumann -z
       +        }
       +        if (z == 0 && bc_bot == 2) {
       +            dev_scalarfield_x[vidx(x,y,nz)] = val_x;     // Periodic -z
       +            dev_scalarfield_y[vidx(x,y,nz)] = val_y;     // Periodic -z
       +            dev_scalarfield_z[vidx(x,y,nz)] = val_z;     // Periodic -z
       +        }
       +        if (z == 1 && bc_bot == 2) {
                    dev_scalarfield_z[vidx(x,y,nz+1)] = val_z;   // Periodic -z
       +        }
        
       -        if (z == nz && bc_top == 0)
       -            dev_scalarfield_z[vidx(x,y,nz)] = val_z;     // Dirichlet
       -        if (z == nz-1 && bc_top == 1)
       -            dev_scalarfield_z[vidx(x,y,nz)] = val_z;     // Neumann
       -        if (z == nz && bc_top == 2)
       +        if (z == nz-1 && bc_top == 0) {
       +            dev_scalarfield_z[vidx(x,y,nz)] = val_z;     // Dirichlet +z
       +        }
       +        if (z == nz-2 && bc_top == 1) {
       +            dev_scalarfield_x[vidx(x,y,nz)] = val_x;     // Neumann +z
       +            dev_scalarfield_y[vidx(x,y,nz)] = val_y;     // Neumann +z
       +            dev_scalarfield_z[vidx(x,y,nz)] = val_z;     // Neumann +z
       +        }
       +        if (z == nz-1 && bc_top == 2) {
       +            dev_scalarfield_x[vidx(x,y,-1)] = val_x;     // Periodic +z
       +            dev_scalarfield_y[vidx(x,y,-1)] = val_y;     // Periodic +z
                    dev_scalarfield_z[vidx(x,y,-1)] = val_z;     // Periodic +z
       +        }
            }
        }