[DO NOT MERGE!] cheating: replace negative rv with small positive values

src/detail/negtozero.hpp

@@ -0,0 +1,34 @@
+#pragma once
+
+// functor for testing not-positive values
+struct isneg { 
+  template<class real_t>
+  BOOST_GPU_ENABLED
+  bool operator()(const real_t a) const {
+    return a <= real_t(0);
+  }
+};
+
+#define _LIBCLOUDPHXX_SMALL_POSITIVE_VALUE 1e-10
+
+// make not-positive values slightly positive...
+#ifdef NDEBUG
+#define negtozero(arr, name) {\
+  thrust::replace_if(\
+    arr.begin(), arr.end(), \
+    isneg(), \
+    real_t(_LIBCLOUDPHXX_SMALL_POSITIVE_VALUE) \
+  );}
+#else
+#define negtozero(arr, name) {\
+  auto min_position = thrust::min_element(arr.begin(), arr.end());\
+  if(*min_position <= 0.)\
+  {\
+    printf("A non-positive number (%g) detected in: " name " at position %g\n", *min_position, min_position);\
+    printf("Cheating in libcloudphxx: replacing non-positive values with small positive ones\n");\
+    thrust::replace_if(\
+      arr.begin(), arr.end(), \
+      isneg(), \
+      real_t(_LIBCLOUDPHXX_SMALL_POSITIVE_VALUE) \
+    );}}
+#endif

src/impl/particles_impl_sstp.ipp

@@ -88,6 +88,7 @@ namespace libcloudphxx
             arg::_1 + arg::_2 / sstp              // op: rv = rv + drv_adv / sstp
           );
         }
+        negtozero((*scl[ix]), "Eulerian field after per-cell sstp_step");
       }
     }
 
@@ -138,6 +139,7 @@ namespace libcloudphxx
       	    arg::_1 + arg::_2                     // op: rv = rv + drv_adv
       	  );
       	}
+        negtozero((*tmp[ix]), "Eulerian field after per-particle sstp_step");
       }
     }
   };  

src/impl/particles_impl_update_th_rv.ipp

@@ -60,7 +60,7 @@ namespace libcloudphxx
         rv.begin(),            // output
         thrust::plus<real_t>() 
       );
-      assert(*thrust::min_element(rv.begin(), rv.end()) >= 0);
+      negtozero(rv, "rv after post-cond updating rv");
       nancheck(rv, "update_th_rv: rv after update");
 
       // updating th

src/particles.tpp

@@ -18,6 +18,7 @@
 #include "detail/wang_collision_enhancement.hpp"
 #include "detail/kernel_onishi_nograv.hpp"
 #include "detail/checknan.hpp"
+#include "detail/negtozero.hpp"
 #include "detail/formatter.cpp"
 #include "detail/tpl_calc_wrapper.hpp"
 #include "detail/kernels.hpp"