ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
ERF_RRTMGP_Utils.H
Go to the documentation of this file.
1 #ifndef ERF_RRTMGP_UTILS_H
2 #define ERF_RRTMGP_UTILS_H
3 
4 #include "rrtmgp_const.h"
5 //#include "rrtmgp_conversion.h"
6 
7 #include "YAKL.h"
8 #include "YAKL_Bounds_fortran.h"
9 
10 #include "ERF_Constants.H"
11 
12 namespace rrtmgp {
13 
14 // Things we need from YAKL
15 using yakl::intrinsics::maxval;
16 using yakl::intrinsics::minval;
17 using yakl::intrinsics::count;
18 using yakl::intrinsics::sum;
19 
20 
21 template<class T, int myMem, int myStyle>
22 void
23 mixing_ratio_to_cloud_mass (yakl::Array<T,2,myMem,myStyle> const &mixing_ratio ,
24  yakl::Array<T,2,myMem,myStyle> const &cloud_fraction,
25  yakl::Array<T,2,myMem,myStyle> const &rho ,
26  yakl::Array<T,2,myMem,myStyle> const &dz ,
27  yakl::Array<T,2,myMem,myStyle> &cloud_mass)
28 {
29  int ncol = mixing_ratio.dimension[0];
30  int nlay = mixing_ratio.dimension[1];
31  yakl::fortran::parallel_for(yakl::fortran::SimpleBounds<2>(nlay, ncol), YAKL_LAMBDA(int ilay, int icol)
32  {
33  // Compute in-cloud mixing ratio (mixing ratio of the cloudy part of the layer)
34  // NOTE: these thresholds (from E3SM) seem arbitrary, but included here for consistency
35  // This limits in-cloud mixing ratio to 0.005 kg/kg. According to note in cloud_diagnostics
36  // in EAM, this is consistent with limits in MG2. Is this true for P3?
37  if (cloud_fraction(icol,ilay) > 0) {
38  // Compute layer-integrated cloud mass (per unit area)
39  auto incloud_mixing_ratio = std::min(mixing_ratio(icol,ilay) / std::max(0.0001, cloud_fraction(icol,ilay)), 0.005);
40  cloud_mass(icol,ilay) = incloud_mixing_ratio * rho(icol,ilay) * dz(icol,ilay);
41  } else {
42  cloud_mass(icol,ilay) = 0;
43  }
44  });
45 }
46 
47 
48 template<class S, class T>
49 void
50 limit_to_bounds (S const &arr_in,
51  T const lower ,
52  T const upper ,
53  S &arr_out)
54 {
55  yakl::c::parallel_for(arr_in.totElems(), YAKL_LAMBDA(int i)
56  {
57  arr_out.data()[i] = std::min(std::max(arr_in.data()[i], lower), upper);
58  });
59 }
60 
61 
62 // Provide a routine to compute heating due to radiative fluxes. This is
63 // computed as net flux into a layer, converted to a heating rate. It is
64 // the responsibility of the user to ensure fields are passed with the
65 // proper units. I.e., pressure at level interfaces should be in Pa,
66 // fluxes in W m-2, Cpair in J kg-1 K-1, gravit in m s-2. This will give
67 // heating in units of K s-1.
68 // TODO: we should probably update this to use the pseudo-density pdel instead
69 // of approximating pdel by differencing the level interface pressures.
70 // We are leaving this for the time being for consistency with SCREAMv0,
71 // from which this code was directly ported.
72 template <class T, int myMem, int myStyle>
73 void
74 compute_heating_rate (yakl::Array<T,2,myMem,myStyle> const &flux_up,
75  yakl::Array<T,2,myMem,myStyle> const &flux_dn,
76  yakl::Array<T,2,myMem,myStyle> const &rho ,
77  yakl::Array<T,2,myMem,myStyle> const &dz ,
78  yakl::Array<T,2,myMem,myStyle> &heating_rate)
79 {
80  auto ncol = flux_up.dimension[0];
81  auto nlay = flux_up.dimension[1]-1;
82  yakl::fortran::parallel_for(yakl::fortran::SimpleBounds<2>(nlay,ncol), YAKL_LAMBDA(int ilay, int icol)
83  {
84  heating_rate(icol,ilay) = ( flux_up(icol,ilay+1) - flux_up(icol,ilay)
85  - flux_dn(icol,ilay+1) + flux_dn(icol,ilay) )
86  / ( -rho(icol,ilay) * dz(icol,ilay) * Cp_d );
87  });
88 }
89 
90 
91 inline
92 bool
93 radiation_do (const int irad,
94  const int nstep)
95 {
96  // If irad == 0, then never do radiation;
97  // Otherwise, we always call radiation at the first step,
98  // and afterwards we do radiation if the timestep is divisible
99  // by irad
100  if (irad == 0) {
101  return false;
102  } else {
103  return ( (nstep == 0) || (nstep % irad == 0) );
104  }
105 }
106 
107 } // namespace rrtmgp
108 #endif
constexpr amrex::Real Cp_d
Definition: ERF_Constants.H:12
@ rho
Definition: ERF_Kessler.H:22
@ T
Definition: ERF_IndexDefines.H:110
Definition: ERF_RRTMGP_Interface.cpp:15
void mixing_ratio_to_cloud_mass(yakl::Array< T, 2, myMem, myStyle > const &mixing_ratio, yakl::Array< T, 2, myMem, myStyle > const &cloud_fraction, yakl::Array< T, 2, myMem, myStyle > const &rho, yakl::Array< T, 2, myMem, myStyle > const &dz, yakl::Array< T, 2, myMem, myStyle > &cloud_mass)
Definition: ERF_RRTMGP_Utils.H:23
bool radiation_do(const int irad, const int nstep)
Definition: ERF_RRTMGP_Utils.H:93
void compute_heating_rate(yakl::Array< T, 2, myMem, myStyle > const &flux_up, yakl::Array< T, 2, myMem, myStyle > const &flux_dn, yakl::Array< T, 2, myMem, myStyle > const &rho, yakl::Array< T, 2, myMem, myStyle > const &dz, yakl::Array< T, 2, myMem, myStyle > &heating_rate)
Definition: ERF_RRTMGP_Utils.H:74
void limit_to_bounds(S const &arr_in, T const lower, T const upper, S &arr_out)
Definition: ERF_RRTMGP_Utils.H:50