tenabled debug output of corrector term - 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 cc9ba7079e9c0d8dc97135b0f977c0070e2b89bc
 (DIR) parent a4c60b5f703e14c9eceaebfa453254a7b6907a8c
 (HTM) Author: Anders Damsgaard <anders.damsgaard@geo.au.dk>
       Date:   Tue, 21 Oct 2014 09:25:48 +0200
       
       enabled debug output of corrector term
       
       Diffstat:
         M src/debug.h                         |      10 +++++++---
         M src/navierstokes.cuh                |       4 ++++
         M src/sphere.cpp                      |       2 +-
         M tests/cfd_tests_neumann.py          |       7 +++++--
       
       4 files changed, 17 insertions(+), 6 deletions(-)
       ---
 (DIR) diff --git a/src/debug.h b/src/debug.h
       t@@ -31,7 +31,8 @@ const int write_conv_log = 1;
        
        // The interval between iteration number reporting in 'output/<sid>-conv.log'
        //const int conv_log_interval = 10;
       -const int conv_log_interval = 4;
       +//const int conv_log_interval = 4;
       +const int conv_log_interval = 1;
        
        // Enable drag force and particle fluid coupling
        #define CFD_DEM_COUPLING
       t@@ -40,10 +41,13 @@ const int conv_log_interval = 4;
        #define CHECK_NS_FINITE
        
        // Enable reporting of velocity prediction components to stdout
       -//#define REPORT_V_P_COMPONENTS
       +#define REPORT_V_P_COMPONENTS
       +
       +// Enable reporting of velocity correction components to stdout
       +#define REPORT_V_C_COMPONENTS
        
        // Enable reporting of forcing function terms to stdout
       -//#define REPORT_FORCING_TERMS
       +#define REPORT_FORCING_TERMS
        
        // Choose solver model (see Zhou et al. 2010 "Discrete particle simulation of
        // particle-fluid flow: model formulations and their applicability", table. 1.
 (DIR) diff --git a/src/navierstokes.cuh b/src/navierstokes.cuh
       t@@ -2996,6 +2996,10 @@ __global__ void updateNSvelocity(
                Float3 v = v_p - ndem*devC_dt*c_grad_p/devC_params.rho_f*grad_epsilon;
        #endif
        
       +#define REPORT_V_C_COMPONENTS
       +        printf("[%d,%d,%d] v_c = %f\t%f\t%f\n",
       +                v.x-v_p.x, v.y-v_p.y, v.z-v_p.z);
       +#endif
                // Print values for debugging
                /* if (z == 0) {
                   Float e_up = dev_ns_epsilon[idx(x,y,z+1)];
 (DIR) diff --git a/src/sphere.cpp b/src/sphere.cpp
       t@@ -149,7 +149,7 @@ void DEM::checkValues(void)
        
            // Check that we have a positive number of particles
            if (np < 1) {
       -        cerr << "Warning: No particles are being simulated (np = " << np
       +        cerr << "Info: No particles are being simulated (np = " << np
                    << ")" << endl;
            }
        
 (DIR) diff --git a/tests/cfd_tests_neumann.py b/tests/cfd_tests_neumann.py
       t@@ -10,6 +10,7 @@ print('### CFD tests - Dirichlet/Neumann BCs ###')
        
        print('''# Neumann bottom, Dirichlet top BC.
        # No gravity, no pressure gradients => no flow''')
       +'''
        orig = sphere.sim("neumann", fluid = True)
        cleanup(orig)
        orig.defaultParams(mu_s = 0.4, mu_d = 0.4)
       t@@ -37,6 +38,7 @@ else:
            print(numpy.mean(py.v_f))
            print(numpy.max(py.v_f))
            raise Exception("Failed")
       +'''
        
        print('''# Neumann bottom, Dirichlet top BC.
        # Gravity, pressure gradients => transient flow''')
       t@@ -44,13 +46,14 @@ orig = sphere.sim("neumann", fluid = True)
        orig.defaultParams(mu_s = 0.4, mu_d = 0.4)
        orig.defineWorldBoundaries([0.4, 0.4, 1], dx = 0.1)
        orig.initFluid(mu = 8.9e-4)
       -orig.initTemporal(total = 0.05, file_dt = 0.005, dt = 1.0e-4)
       +#orig.initTemporal(total = 0.05, file_dt = 0.005, dt = 1.0e-4)
       +orig.initTemporal(total = 1.0e-3, file_dt = 1.0e-4, dt = 1.0e-4)
        py = sphere.sim(sid = orig.sid, fluid = True)
        orig.g[2] = -10.0
        orig.bc_bot[0] = 1      # No-flow BC at bottom (Neumann)
        #orig.run(dry=True)
        orig.run(verbose=False)
       -#orig.writeVTKall()
       +orig.writeVTKall()
        py.readlast(verbose = False)
        #ideal_grad_p_z = numpy.linspace(
        #        orig.p_f[0,0,0] + orig.L[2]*orig.rho_f*numpy.abs(orig.g[2]),