Port shoc_assumed_pdf functions individually (like in shoc.F90)

components/eamxx/src/physics/shoc/impl/shoc_assumed_pdf_compute_buoyancy_flux_impl.hpp

@@ -0,0 +1,41 @@
+#ifndef SHOC_SHOC_ASSUMED_PDF_COMPUTE_BUOYANCY_FLUX_IMPL_HPP
+#define SHOC_SHOC_ASSUMED_PDF_COMPUTE_BUOYANCY_FLUX_IMPL_HPP
+
+#include "shoc_functions.hpp" // for ETI only but harmless for GPU
+
+#include <iomanip>
+
+namespace scream {
+namespace shoc {
+
+/*
+ * Implementation of shoc_assumed_pdf_compute_buoyancy_flux. Clients should NOT
+ * #include this file, but include shoc_functions.hpp instead.
+ */
+
+template<typename S, typename D>
+KOKKOS_INLINE_FUNCTION
+void Functions<S,D>::shoc_assumed_pdf_compute_buoyancy_flux(
+  const Spack& wthlsec,
+  const Spack& wqwsec,
+  const Spack& pval,
+  const Spack& wqls,
+  Spack&       wthv_sec)
+{
+  const Scalar basepres = C::P0;
+  const Scalar rair = C::Rair;
+  const Scalar rv = C::RV;
+  const Scalar cp = C::CP;
+  const Scalar lcond = C::LatVap;
+  const Scalar basetemp = C::basetemp;
+  const Scalar epsterm = rair/rv;
+
+  wthv_sec = wthlsec + ((1 - epsterm)/epsterm)*basetemp*wqwsec
+             + ((lcond/cp)*ekat::pow(basepres/pval, (rair/cp))
+             - (1/epsterm)*basetemp)*wqls;
+}
+
+} // namespace shoc
+} // namespace scream
+
+#endif

components/eamxx/src/physics/shoc/impl/shoc_assumed_pdf_compute_cloud_liquid_variance_impl.hpp

@@ -0,0 +1,39 @@
+#ifndef SHOC_SHOC_ASSUMED_PDF_COMPUTE_CLOUD_LIQUID_VARIANCE_IMPL_HPP
+#define SHOC_SHOC_ASSUMED_PDF_COMPUTE_CLOUD_LIQUID_VARIANCE_IMPL_HPP
+
+#include "shoc_functions.hpp" // for ETI only but harmless for GPU
+
+#include <iomanip>
+
+namespace scream {
+namespace shoc {
+
+/*
+ * Implementation of shoc_assumed_pdf_compute_cloud_liquid_variance. Clients should NOT
+ * #include this file, but include shoc_functions.hpp instead.
+ */
+
+template<typename S, typename D>
+KOKKOS_INLINE_FUNCTION
+void Functions<S,D>::shoc_assumed_pdf_compute_cloud_liquid_variance(
+  const Spack& a,
+  const Spack& s1,
+  const Spack& ql1,
+  const Spack& C1,
+  const Spack& std_s1,
+  const Spack& s2,
+  const Spack& ql2,
+  const Spack& C2,
+  const Spack& std_s2,
+  const Spack& shoc_ql,
+  Spack&       shoc_ql2)
+{
+  shoc_ql2 = ekat::max(0, a*(s1*ql1 + C1*ekat::square(std_s1))
+                          + (1 - a)*(s2*ql2 + C2*ekat::square(std_s2))
+                          - ekat::square(shoc_ql));
+}
+
+} // namespace shoc
+} // namespace scream
+
+#endif

components/eamxx/src/physics/shoc/impl/shoc_assumed_pdf_compute_liquid_water_flux_impl.hpp

@@ -0,0 +1,33 @@
+#ifndef SHOC_SHOC_ASSUMED_PDF_COMPUTE_LIQUID_WATER_FLUX_IMPL_HPP
+#define SHOC_SHOC_ASSUMED_PDF_COMPUTE_LIQUID_WATER_FLUX_IMPL_HPP
+
+#include "shoc_functions.hpp" // for ETI only but harmless for GPU
+
+#include <iomanip>
+
+namespace scream {
+namespace shoc {
+
+/*
+ * Implementation of shoc_assumed_pdf_compute_liquid_water_flux. Clients should NOT
+ * #include this file, but include shoc_functions.hpp instead.
+ */
+
+template<typename S, typename D>
+KOKKOS_INLINE_FUNCTION
+void Functions<S,D>::shoc_assumed_pdf_compute_liquid_water_flux(
+  const Spack& a,
+  const Spack& w1_1,
+  const Spack& w_first,
+  const Spack& ql1,
+  const Spack& w1_2,
+  const Spack& ql2,
+  Spack&       wqls)
+{
+  wqls = a*((w1_1 - w_first)*ql1) + (1 - a)*((w1_2 - w_first)*ql2);
+}
+
+} // namespace shoc
+} // namespace scream
+
+#endif

components/eamxx/src/physics/shoc/impl/shoc_assumed_pdf_compute_qs_impl.hpp

@@ -0,0 +1,55 @@
+#ifndef SHOC_SHOC_ASSUMED_PDF_COMPUTE_QS_IMPL_HPP
+#define SHOC_SHOC_ASSUMED_PDF_COMPUTE_QS_IMPL_HPP
+
+#include "shoc_functions.hpp" // for ETI only but harmless for GPU
+#include "physics_functions.hpp"
+
+#include <iomanip>
+
+namespace scream {
+namespace shoc {
+
+/*
+ * Implementation of shoc_assumed_pdf_compute_qs. Clients should NOT
+ * #include this file, but include shoc_functions.hpp instead.
+ */
+
+template<typename S, typename D>
+KOKKOS_INLINE_FUNCTION
+void Functions<S,D>::shoc_assumed_pdf_compute_qs(
+  const Spack& Tl1_1,
+  const Spack& Tl1_2,
+  const Spack& pval,
+  const Smask& active_entries,
+  Spack&       qs1,
+  Spack&       beta1,
+  Spack&       qs2,
+  Spack&       beta2)
+{
+  const Scalar rair = C::Rair;
+  const Scalar rv = C::RV;
+  const Scalar cp = C::CP;
+  const Scalar lcond = C::LatVap;
+
+  // Compute MurphyKoop_svp
+  const int liquid = 0;
+  const Spack esval1_1 = scream::physics::Functions<S,D>::MurphyKoop_svp(Tl1_1,liquid,active_entries,"shoc::shoc_assumed_pdf (Tl1_1)");
+  const Spack esval1_2 = scream::physics::Functions<S,D>::MurphyKoop_svp(Tl1_2,liquid,active_entries,"shoc::shoc_assumed_pdf (Tl1_2)");
+  const Spack lstarn(lcond);
+
+  qs1 = sp(0.622)*esval1_1/ekat::max(esval1_1, pval - esval1_1);
+  beta1 = (rair/rv)*(lstarn/(rair*Tl1_1))*(lstarn/(cp*Tl1_1));
+
+  // Only compute qs2 and beta2 if the two plumes are not equal
+  const Smask condition = (Tl1_1 != Tl1_2);
+  qs2 = qs1;
+  beta2 = beta1;
+
+  qs2.set(condition, sp(0.622)*esval1_2/ekat::max(esval1_2, pval - esval1_2));
+  beta2.set(condition, (rair/rv)*(lstarn/(rair*Tl1_2))*(lstarn/(cp*Tl1_2)));
+}
+
+} // namespace shoc
+} // namespace scream
+
+#endif

components/eamxx/src/physics/shoc/impl/shoc_assumed_pdf_compute_s_impl.hpp

@@ -0,0 +1,69 @@
+#ifndef SHOC_SHOC_ASSUMED_PDF_COMPUTE_S_IMPL_HPP
+#define SHOC_SHOC_ASSUMED_PDF_COMPUTE_S_IMPL_HPP
+
+#include "shoc_functions.hpp" // for ETI only but harmless for GPU
+
+#include <iomanip>
+
+namespace scream {
+namespace shoc {
+
+/*
+ * Implementation of shoc_assumed_pdf_compute_s. Clients should NOT
+ * #include this file, but include shoc_functions.hpp instead.
+ */
+
+template<typename S, typename D>
+KOKKOS_INLINE_FUNCTION
+void Functions<S,D>::shoc_assumed_pdf_compute_s(
+  const Spack& qw1,
+  const Spack& qs,
+  const Spack& beta,
+  const Spack& pval,
+  const Spack& thl2,
+  const Spack& qw2,
+  const Spack& sqrtthl2,
+  const Spack& sqrtqw2,
+  const Spack& r_qwthl,
+  Spack&       s,
+  Spack&       std_s,
+  Spack&       qn,
+  Spack&       C)
+{
+  const Scalar rair = C::Rair;
+  const Scalar basepres = C::P0;
+  const Scalar cp = C::CP;
+  const Scalar lcond = C::LatVap;
+  const Scalar pi = C::Pi;
+
+  const Scalar sqrt2(std::sqrt(Scalar(2.0))), sqrt2pi(std::sqrt(2*pi));
+
+  const Spack cthl=((1 + beta*qw1)/ekat::square(1 + beta*qs))*(cp/lcond)*
+                    beta*qs*ekat::pow(pval/basepres, (rair/cp));
+  const Spack cqt = 1/(1 + beta*qs);
+
+  std_s = ekat::sqrt(ekat::max(0,
+                               ekat::square(cthl)*thl2
+                               + ekat::square(cqt)*qw2 - 2*cthl*sqrtthl2*cqt*sqrtqw2*r_qwthl));
+  const auto std_s_not_small = std_s > std::sqrt(std::numeric_limits<Scalar>::min()) * 100;
+  s = qw1-qs*((1 + beta*qw1)/(1 + beta*qs));
+  if (std_s_not_small.any()) {
+    C.set(std_s_not_small, sp(0.5)*(1 + ekat::erf(s/(sqrt2*std_s))));
+  }
+  C.set(!std_s_not_small && s > 0, 1);
+  const auto std_s_C_not_small = std_s_not_small && C != 0;
+  if (std_s_C_not_small.any()) {
+    qn.set(std_s_C_not_small, s*C+(std_s/sqrt2pi)*ekat::exp(-sp(0.5)*ekat::square(s/std_s)));
+  }
+  qn.set(!std_s_not_small && s > 0, s);
+
+  // Checking to prevent empty clouds
+  const auto qn_le_zero = qn <= 0;
+  C.set(qn_le_zero,0);
+  qn.set(qn_le_zero,0);
+}
+
+} // namespace shoc
+} // namespace scream
+
+#endif

components/eamxx/src/physics/shoc/impl/shoc_assumed_pdf_compute_sgs_liquid_impl.hpp

@@ -0,0 +1,30 @@
+#ifndef SHOC_SHOC_ASSUMED_PDF_COMPUTE_SGS_LIQUID_IMPL_HPP
+#define SHOC_SHOC_ASSUMED_PDF_COMPUTE_SGS_LIQUID_IMPL_HPP
+
+#include "shoc_functions.hpp" // for ETI only but harmless for GPU
+
+#include <iomanip>
+
+namespace scream {
+namespace shoc {
+
+/*
+ * Implementation of shoc_assumed_pdf_compute_sgs_liquid. Clients should NOT
+ * #include this file, but include shoc_functions.hpp instead.
+ */
+
+template<typename S, typename D>
+KOKKOS_INLINE_FUNCTION
+void Functions<S,D>::shoc_assumed_pdf_compute_sgs_liquid(
+  const Spack& a,
+  const Spack& ql1,
+  const Spack& ql2,
+  Spack&       shoc_ql)
+{
+  shoc_ql = ekat::max(0, a*ql1 + (1 - a)*ql2);
+}
+
+} // namespace shoc
+} // namespace scream
+
+#endif

components/eamxx/src/physics/shoc/impl/shoc_assumed_pdf_compute_temperature_impl.hpp

@@ -0,0 +1,32 @@
+#ifndef SHOC_SHOC_ASSUMED_PDF_COMPUTE_TEMPERATURE_IMPL_HPP
+#define SHOC_SHOC_ASSUMED_PDF_COMPUTE_TEMPERATURE_IMPL_HPP
+
+#include "shoc_functions.hpp" // for ETI only but harmless for GPU
+
+#include <iomanip>
+
+namespace scream {
+namespace shoc {
+
+/*
+ * Implementation of shoc_assumed_pdf_compute_temperature. Clients should NOT
+ * #include this file, but include shoc_functions.hpp instead.
+ */
+
+template<typename S, typename D>
+KOKKOS_INLINE_FUNCTION
+void Functions<S,D>::shoc_assumed_pdf_compute_temperature(
+  const Spack& thl1,
+  const Spack& pval,
+  Spack&       Tl1)
+{
+  constexpr Scalar basepres = C::P0;
+  constexpr Scalar rair = C::Rair;
+  constexpr Scalar cp = C::CP;
+  Tl1 = thl1/(ekat::pow(basepres/pval,(rair/cp)));
+}
+
+} // namespace shoc
+} // namespace scream
+
+#endif