Merge branch 'cuda_math_opt' into cache_mfp

Author: pdziekan

include/libcloudph++/blk_2m/activation_formulae.hpp

@@ -36,8 +36,8 @@ namespace libcloudphxx
         const quantity<si::dimensionless, real_t> &chem_b,
         const quantity<si::dimensionless, real_t> beta = beta_default<real_t>()
       ) {
-        return pow(mean_rd / si::metres, -(1+beta)) * sqrt(4 * pow(common::kelvin::A<real_t>(T) / si::metres, 3) / 27 / chem_b);
-      }           // can't have pow<3/2>(si:: )
+        return std::pow(mean_rd / si::metres, -(1+beta)) * std::sqrt(4 * std::pow(common::kelvin::A<real_t>(T) / si::metres, 3) / 27 / chem_b);
+      }           // can't have std::pow<3/2>(si:: )
 
       // helper for activation formulae (see eq. 10 in Morrison and Grabowski 2007)
       template <typename real_t>
@@ -55,7 +55,7 @@ namespace libcloudphxx
         const quantity<si::dimensionless, real_t> &sdev_rd,
         const quantity<si::dimensionless, real_t> beta = beta_default<real_t>()
       ) {
-        return pow(sdev_rd, 1+beta);
+        return std::pow(sdev_rd, 1+beta);
       }
 
       //helper for activation formulae (see eq. 10 in Morrison and Grabowski 2007)
@@ -73,7 +73,7 @@ namespace libcloudphxx
         return log(
           s_0(T, mean_rd, chem_b) / 
           std::min(real_t(s(p, T, rv)), real_t(RH_max - 1))
-        ) / sqrt(2) / log(sdev_rd_s(sdev_rd));
+        ) / std::sqrt(2) / log(sdev_rd_s(sdev_rd));
       }
 
       // helper for activation formulae (see eq. 10 in Morrison and Grabowski 2007)

include/libcloudph++/blk_2m/common_formulae.hpp

@@ -61,7 +61,7 @@ namespace libcloudphxx
          const quantity<si::dimensionless, real_t> &rc,
         const quantity<si::mass_density, real_t> &rhod
       ) {
-        auto tmp = pow(
+        auto tmp = std::pow(
 	  c_md<real_t>() * n * std::tgamma(miu_c(n*rhod) + d_md<real_t>() + real_t(1))
           /
 	  (rc * std::tgamma(miu_c(n*rhod) + real_t(1))) * si::cubic_metres
@@ -83,7 +83,7 @@ namespace libcloudphxx
          const quantity<si::mass_density, real_t> &rhod
       ) {
         auto tmp = nc
-          * pow(lambda_c(nc, rc, rhod) * si::metres, miu_c(nc * rhod) + real_t(1))
+          * std::pow(lambda_c(nc, rc, rhod) * si::metres, miu_c(nc * rhod) + real_t(1))
           / std::tgamma(miu_c(nc * rhod) + real_t(1))
           / si::metres;
         assert(finite(tmp * si::metres * si::cubic_metres) && "intercept N0_c is finite failed");
@@ -98,7 +98,7 @@ namespace libcloudphxx
          const quantity<divide_typeof_helper<si::dimensionless, si::mass>::type, real_t> &nr,
          const quantity<si::dimensionless, real_t> &rr
       ) {
-        auto tmp = pow(
+        auto tmp = std::pow(
 	  c_md<real_t>() * nr * std::tgamma(d_md<real_t>() + real_t(1)) / rr * si::cubic_metres
           ,
 	  real_t(1) / d_md<real_t>()

include/libcloudph++/blk_2m/terminal_vel_formulae.hpp

@@ -60,8 +60,8 @@ namespace libcloudphxx
           const quantity<divide_typeof_helper<si::dimensionless, si::length>::type, real_t> &lbd, //slope of assumed exponential size distribution
           const quantity<si::length, real_t> &D // diameter
         ) {
-          auto tmp =  - pow(lbd * si::metres, -6) * exp(-lbd * D) *
-            (pow(real_t(lbd * D), 5) + 5 * pow(real_t(lbd * D), 4) * 20 * pow(real_t(lbd * D), 3) + 60 * pow(real_t(lbd * D), 2) + 120 * real_t(lbd * D) + 120);
+          auto tmp =  - std::pow(lbd * si::metres, -6) * exp(-lbd * D) *
+            (pow(real_t(lbd * D), 5) + 5 * std::pow(real_t(lbd * D), 4) * 20 * std::pow(real_t(lbd * D), 3) + 60 * std::pow(real_t(lbd * D), 2) + 120 * real_t(lbd * D) + 120);
 
           assert(finite(tmp) && "mint_1 is finite failed");
           return tmp;
@@ -72,8 +72,8 @@ namespace libcloudphxx
           const quantity<divide_typeof_helper<si::dimensionless, si::length>::type, real_t> &lbd, //slope of assumed exponential size distribution
           const quantity<si::length, real_t> &D // diameter
         ) {
-          auto tmp = - pow(real_t(lbd * si::metres), -5) * exp(-lbd * D) *
-            (pow(real_t(lbd * D), 4) + 4 * pow(real_t(lbd * D), 3) * 12 * pow(real_t(lbd * D), 2) + 24 * real_t(lbd * D) + 24);
+          auto tmp = - std::pow(real_t(lbd * si::metres), -5) * exp(-lbd * D) *
+            (pow(real_t(lbd * D), 4) + 4 * std::pow(real_t(lbd * D), 3) * 12 * std::pow(real_t(lbd * D), 2) + 24 * real_t(lbd * D) + 24);
 
           assert(finite(tmp) && "mint_2 is finite failed");
           return tmp;
@@ -84,9 +84,9 @@ namespace libcloudphxx
           const quantity<divide_typeof_helper<si::dimensionless, si::length>::type, real_t> &lbd, // slope of assumed exponential size distribution
           const quantity<si::length, real_t> &D // diameter
         ) {
-          auto tmp = real_t(1./16) / pow(real_t(lbd * si::metres), real_t(9./2))
-                   * (105 * sqrt(pi<real_t>()) * std::erf(sqrt(lbd * D))
-                      - 2 * sqrt(lbd * D) * exp(-lbd * D) * (8 * pow(real_t(lbd * D), 3) + 28 * pow(real_t(lbd * D), 2) + 70 * real_t(lbd * D) + 105)
+          auto tmp = real_t(1./16) / std::pow(real_t(lbd * si::metres), real_t(9./2))
+                   * (105 * std::sqrt(pi<real_t>()) * std::erf(sqrt(lbd * D))
+                      - 2 * std::sqrt(lbd * D) * exp(-lbd * D) * (8 * std::pow(real_t(lbd * D), 3) + 28 * std::pow(real_t(lbd * D), 2) + 70 * real_t(lbd * D) + 105)
                    );
 
           assert(finite(tmp) && "mint_3 is finite failed");
@@ -99,7 +99,7 @@ namespace libcloudphxx
           const quantity<divide_typeof_helper<si::dimensionless, si::length>::type, real_t> &lbd, //slope of assumed exponential size distribution
           const quantity<si::length, real_t> &D // diameter
         ) {
-          auto tmp = - pow(real_t(lbd * si::metres), -1) * exp(-lbd * D);
+          auto tmp = - std::pow(real_t(lbd * si::metres), -1) * exp(-lbd * D);
 
           assert(finite(tmp) && "int_4 is finite failed");
           return tmp;
@@ -112,7 +112,7 @@ namespace libcloudphxx
          const quantity<divide_typeof_helper<si::dimensionless, si::length>::type, real_t> &lbd, // slope of assumed exponential size distribution
          const quantity<si::length, real_t> &D // diameter
        ) {
-         auto tmp = pow(real_t(lbd * si::metres), -3) * exp(-lbd * D) * (- real_t(lbd * D) * (real_t(lbd * D) + 2) - 2);
+         auto tmp = std::pow(real_t(lbd * si::metres), -3) * exp(-lbd * D) * (- real_t(lbd * D) * (real_t(lbd * D) + 2) - 2);
 
          assert(finite(tmp) && "nint_1 is finite failed");
          return tmp;
@@ -124,7 +124,7 @@ namespace libcloudphxx
          const quantity<divide_typeof_helper<si::dimensionless, si::length>::type, real_t> &lbd, //slope of assumed exponential size distribution
          const quantity<si::length, real_t> &D // diameter
        ) {
-         auto tmp =  - pow(real_t(lbd * si::metres), -2) * exp(-lbd * D) * (lbd * D + 1);
+         auto tmp =  - std::pow(real_t(lbd * si::metres), -2) * exp(-lbd * D) * (lbd * D + 1);
 
          assert(finite(tmp) && "nint_2 is finite failed");
          return tmp;
@@ -136,8 +136,8 @@ namespace libcloudphxx
          const quantity<divide_typeof_helper<si::dimensionless, si::length>::type, real_t> &lbd, //slope of assumed exponential size distribution
          const quantity<si::length, real_t> &D //diameter
        ) {
-         auto tmp = sqrt(pi<real_t>()) * std::erf(sqrt(lbd * D)) / 2 / pow(real_t(lbd * si::metres), 3./2)
-           - sqrt(D / si::metres) * exp(-lbd*D) / (lbd * si::metres);
+         auto tmp = std::sqrt(pi<real_t>()) * std::erf(sqrt(lbd * D)) / 2 / std::pow(real_t(lbd * si::metres), 3./2)
+           - std::sqrt(D / si::metres) * exp(-lbd*D) / (lbd * si::metres);
 
          assert(finite(tmp) && "nint_3 is finite failed");
          return tmp;
@@ -161,15 +161,15 @@ namespace libcloudphxx
           * real_t(1000) // mass of the drop in grams
           * (
 	    alpha_fall(d1<real_t>() / real_t(2))
-            * pow(c_md<real_t>() * si::cubic_metres / si::kilograms * real_t(1000), beta_fall(d1<real_t>() / real_t(2)))
+            * std::pow(c_md<real_t>() * si::cubic_metres / si::kilograms * real_t(1000), beta_fall(d1<real_t>() / real_t(2)))
 	    * (mint_1(lbd, d1<real_t>()) - mint_1(lbd, real_t(0) * si::metres))
 	    +
 	    alpha_fall(d1<real_t>() + d2<real_t>() / real_t(2))
-            * pow(c_md<real_t>() * si::cubic_metres / si::kilograms * real_t(1000), beta_fall(d1<real_t>() + d2<real_t>() / real_t(2)))
+            * std::pow(c_md<real_t>() * si::cubic_metres / si::kilograms * real_t(1000), beta_fall(d1<real_t>() + d2<real_t>() / real_t(2)))
 	    * (mint_2(lbd, d2<real_t>()) - mint_2(lbd, d1<real_t>()))
 	    +
 	    alpha_fall(d2<real_t>() + d3<real_t>() / real_t(2))
-            * pow(c_md<real_t>() * si::cubic_metres / si::kilograms * real_t(1000), beta_fall(d2<real_t>() + d3<real_t>() / real_t(2)))
+            * std::pow(c_md<real_t>() * si::cubic_metres / si::kilograms * real_t(1000), beta_fall(d2<real_t>() + d3<real_t>() / real_t(2)))
 	    * (mint_3(lbd, d3<real_t>()) - mint_3(lbd, d2<real_t>()))
 	    +
 	    alpha_fall(real_t(2) * d3<real_t>())
@@ -197,15 +197,15 @@ namespace libcloudphxx
         auto tmp = rho_stp<real_t>() / rhod
           * (
 	   alpha_fall(d1<real_t>() / real_t(2))
-           * pow(c_md<real_t>() * si::cubic_metres / si::kilograms * real_t(1000), beta_fall(d1<real_t>() / real_t(2)))
+           * std::pow(c_md<real_t>() * si::cubic_metres / si::kilograms * real_t(1000), beta_fall(d1<real_t>() / real_t(2)))
            * (nint_1(lbd, d1<real_t>()) - nint_1(lbd, real_t(0) * si::metres))
 	   +
            alpha_fall(d1<real_t>() + d2<real_t>() / real_t(2))
-           * pow(c_md<real_t>() * si::cubic_metres / si::kilograms * real_t(1000), beta_fall(d1<real_t>() + d2<real_t>() / real_t(2)))
+           * std::pow(c_md<real_t>() * si::cubic_metres / si::kilograms * real_t(1000), beta_fall(d1<real_t>() + d2<real_t>() / real_t(2)))
            * (nint_2(lbd, d2<real_t>()) - nint_2(lbd, d1<real_t>()))
            +
            alpha_fall(d2<real_t>() + d3<real_t>() / real_t(2))
-           * pow(c_md<real_t>() * si::cubic_metres / si::kilograms * real_t(1000), beta_fall(d2<real_t>() + d3<real_t>() / real_t(2)))
+           * std::pow(c_md<real_t>() * si::cubic_metres / si::kilograms * real_t(1000), beta_fall(d2<real_t>() + d3<real_t>() / real_t(2)))
            * (nint_3(lbd, d3<real_t>()) - nint_3(lbd, d2<real_t>()))
            +
            alpha_fall(real_t(2) * d3<real_t>()) * (real_t(0) - int_4(lbd, d3<real_t>()))

include/libcloudph++/common/GA17_turbulence.hpp

@@ -7,6 +7,13 @@ namespace libcloudphxx
     // Grabowski and Abade 2017
     namespace GA17_turbulence
     {
+
+#if !defined(__NVCC__)
+      using std::pow;
+      using std::cbrt;
+      using std::sqrt;
+#endif
+
       libcloudphxx_const_derived(si::dimensionless, C_E, real_t(0.845));
       libcloudphxx_const_derived(si::dimensionless, C_tau, real_t(1.5));
       libcloudphxx_const_derived(si::dimensionless, cube_root_of_two_pi, pow(real_t(2) *
@@ -18,13 +25,6 @@ namespace libcloudphxx
         , real_t(1./3.))
       );
 
-/*
-#if !defined(__NVCC__)
-      //using std::pow;
-      using std::cbrt;
-      //using std::sqrt;
-#endif
-*/
       typedef divide_typeof_helper<
         si::dimensionless,  
         si::length

include/libcloudph++/common/detail/fake_units.hpp

@@ -48,7 +48,11 @@ namespace libcloudphxx
           const quantity<qntt_t, real_t> &b
         )
         {
+#if !defined(__NVCC__)
+          return quantity<qntt_t, real_t>(std::pow(a.value, b.value)); // pow from the parent namespace, hopefuly this will be built-in CUDA pow
+#else
           return quantity<qntt_t, real_t>(::pow(a.value, b.value)); // pow from the parent namespace, hopefuly this will be built-in CUDA pow
+#endif
         }
 
         // ...

include/libcloudph++/common/henry.hpp

@@ -73,6 +73,9 @@ namespace libcloudphxx
         const quantity<si::temperature, real_t> &T,
         const quantity<mass_over_amount, real_t> &M
       ) {
+#if !defined(__NVCC__)
+        using std::sqrt;
+#endif
         return (
           real_t(  //bug in boost #6957
             sqrt(
@@ -99,6 +102,9 @@ namespace libcloudphxx
         const quantity<si::temperature, real_t> &T,
         const quantity<mass_over_amount, real_t> &M
       ) {
+#if !defined(__NVCC__)
+        using std::sqrt;
+#endif
         return(
           real_t(1.) / (
             rw2 / real_t(3.) / D 

include/libcloudph++/common/hydrostatic.hpp

@@ -11,6 +11,9 @@ namespace libcloudphxx
   {
     namespace hydrostatic
     {
+#if !defined(__NVCC__)
+      using std::pow;
+#endif
       using theta_std::p_1000;
       using moist_air::R; 
       using moist_air::c_pd;