Merge pull request #479 from AgnieszkaMakulska/ice

Ice in SD debug

Author: pdziekan

src/impl/ice/particles_impl_ice_dep.ipp

@@ -5,12 +5,6 @@
   * GPLv3+ (see the COPYING file or http://www.gnu.org/licenses/)
   */
 
-// #include <thrust/iterator/transform_iterator.h>
-// #include <libcloudph++/common/maxwell-mason.hpp>
-// #include <libcloudph++/common/kappa_koehler.hpp>
-// #include <libcloudph++/common/kelvin_term.hpp>
-// #include <libcloudph++/common/transition_regime.hpp>
-// #include <libcloudph++/common/ventil.hpp>
 
 namespace libcloudphxx
 {
@@ -31,8 +25,6 @@ namespace libcloudphxx
       hskpng_sort();
 
       // Vector to store 3rd moment
-      // auto d_ice_mass_g = tmp_device_real_cell.get_guard();
-      // thrust_device::vector<real_t> &d_ice_mass = d_ice_mass_g.get();
       if(step == 0)
         reset_guardp(ice_mass_gp, tmp_device_real_cell);
       thrust_device::vector<real_t> &ice_mass = ice_mass_gp->get();
@@ -137,8 +129,6 @@ namespace libcloudphxx
         ice_mass_gp.reset(); // destroy guard to tmp array that stored ice_mass
       }
 
-      // update th and rv according to change in third specific moment
-      // update_th_rv(d_ice_mass, impl::phase_change::deposition);
       }
     };
   }

src/particles_diag.ipp

@@ -250,13 +250,17 @@ namespace libcloudphxx
     template <typename real_t, backend_t device>
     void particles_t<real_t, device>::diag_ice_a_rng(const real_t &a_min, const real_t &a_max)
     {
+      if(pimpl->opts_init.ice_switch == false)
+        throw std::runtime_error("libcloudph++: ice is switched off in opts_init, but diag_ice was called");
       pimpl->moms_rng(a_min, a_max, pimpl->ice_a.begin(), false);
     }
 
     // selects particles with (ice_c >= c_min && ice_c < c_max)
     template <typename real_t, backend_t device>
     void particles_t<real_t, device>::diag_ice_c_rng(const real_t &c_min, const real_t &c_max)
     {
+      if(pimpl->opts_init.ice_switch == false)
+        throw std::runtime_error("libcloudph++: ice is switched off in opts_init, but diag_ice was called");
       pimpl->moms_rng(c_min, c_max, pimpl->ice_c.begin(), false);
     }
 
@@ -271,14 +275,16 @@ namespace libcloudphxx
     template <typename real_t, backend_t device>
     void particles_t<real_t, device>::diag_ice()
     {
+      if(pimpl->opts_init.ice_switch == false)
+        throw std::runtime_error("libcloudph++: ice is switched off in opts_init, but diag_ice was called");
       pimpl->moms_gt0(pimpl->ice_a.begin()); // ice_a greater than 0
     }
 
     // selects water particles 
     template <typename real_t, backend_t device>
     void particles_t<real_t, device>::diag_water()
     {
-      pimpl->moms_eq0(pimpl->ice_a.begin()); // ice_a equal to 0
+      pimpl->moms_gt0(pimpl->rw2.begin()); // rw2 greater than 0
     }
 
     // selects particles with (r_d >= r_min && r_d < r_max) from particles previously selected
@@ -305,13 +311,17 @@ namespace libcloudphxx
     template <typename real_t, backend_t device>
     void particles_t<real_t, device>::diag_ice_a_rng_cons(const real_t &a_min, const real_t &a_max)
     {
+      if(pimpl->opts_init.ice_switch == false)
+        throw std::runtime_error("libcloudph++: ice is switched off in opts_init, but diag_ice was called");
       pimpl->moms_rng(a_min, a_max, pimpl->ice_a.begin(), true);
     }
 
     // selects particles with (ice_c >= c_min && ice_c < c_max) from particles previously selected
     template <typename real_t, backend_t device>
     void particles_t<real_t, device>::diag_ice_c_rng_cons(const real_t &c_min, const real_t &c_max)
     {
+      if(pimpl->opts_init.ice_switch == false)
+        throw std::runtime_error("libcloudph++: ice is switched off in opts_init, but diag_ice was called");
       pimpl->moms_rng(c_min, c_max, pimpl->ice_c.begin(), true);
     }
 
@@ -326,14 +336,16 @@ namespace libcloudphxx
     template <typename real_t, backend_t device>
     void particles_t<real_t, device>::diag_ice_cons()
     {
+      if(pimpl->opts_init.ice_switch == false)
+        throw std::runtime_error("libcloudph++: ice is switched off in opts_init, but diag_ice was called");
       pimpl->moms_gt0(pimpl->ice_a.begin(), true); // ice_a greater than 0
     }
 
     // selects water particles from particles previously selected 
     template <typename real_t, backend_t device>
     void particles_t<real_t, device>::diag_water_cons()
     {
-      pimpl->moms_eq0(pimpl->ice_a.begin(), true); // ice_a equal to 0
+      pimpl->moms_gt0(pimpl->rw2.begin(), true); // rw2 greater than 0
     }
 
     // selects particles with RH >= Sc   (Sc - critical supersaturation)
@@ -412,20 +424,26 @@ namespace libcloudphxx
     template <typename real_t, backend_t device>
     void particles_t<real_t, device>::diag_ice_a_mom(const int &n)
     {
+      if(pimpl->opts_init.ice_switch == false)
+        throw std::runtime_error("libcloudph++: ice is switched off in opts_init, but diag_ice was called");
       pimpl->moms_calc(pimpl->ice_a.begin(), n);
     }
 
     // computes n-th moment of the ice polar radius spectrum for the selected particles
     template <typename real_t, backend_t device>
     void particles_t<real_t, device>::diag_ice_c_mom(const int &n)
     {
+      if(pimpl->opts_init.ice_switch == false)
+        throw std::runtime_error("libcloudph++: ice is switched off in opts_init, but diag_ice was called");
       pimpl->moms_calc(pimpl->ice_c.begin(), n);
     }
 
      // computes ice mixing ratio
     template <typename real_t, backend_t device>
       void particles_t<real_t, device>::diag_ice_mix_ratio()
      {
+      if(pimpl->opts_init.ice_switch == false)
+        throw std::runtime_error("libcloudph++: ice is switched off in opts_init, but diag_ice was called");
       pimpl->moms_calc(thrust::make_transform_iterator(
         thrust::make_zip_iterator(thrust::make_tuple(pimpl->ice_a.begin(), pimpl->ice_c.begin(), pimpl->ice_rho.begin())),
         detail::ice_mass<real_t>()
@@ -571,6 +589,9 @@ namespace libcloudphxx
     template <typename real_t, backend_t device>
     void particles_t<real_t, device>::diag_precip_rate_ice_mass()
     {
+      if(pimpl->opts_init.ice_switch == false)
+        throw std::runtime_error("libcloudph++: ice is switched off in opts_init, but diag_ice was called");
+
       // updating terminal velocities
       pimpl->hskpng_vterm_all();
 
@@ -623,7 +644,8 @@ namespace libcloudphxx
     template <typename real_t, backend_t device>
     void particles_t<real_t, device>::diag_chem(const enum chem_species_t &c)
     {
-      if(pimpl->opts_init.chem_switch == false) throw std::runtime_error("libcloudph++: all chemistry was switched off in opts_init");
+      if(pimpl->opts_init.chem_switch == false)
+        throw std::runtime_error("libcloudph++: chemistry is switched off in opts_init, but diag_chem was called");
       pimpl->moms_calc(pimpl->chem_bgn[c], 1.);
     }