ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
SHOCInterface::SHOCPreprocess Struct Reference

#include <ERF_ShocInterface.H>

Collaboration diagram for SHOCInterface::SHOCPreprocess:

Public Member Functions

 SHOCPreprocess ()=default
 
KOKKOS_INLINE_FUNCTION void operator() (const Kokkos::TeamPolicy< KT::ExeSpace >::member_type &team) const
 
void set_variables (const int ncol_, const int nlev_, const Real z_surf_, const view_2d_const &T_mid_, const view_2d_const &p_mid_, const view_2d_const &p_int_, const view_2d_const &pseudo_density_, const view_2d_const &omega_, const view_1d_const &phis_, const view_1d_const &surf_sens_flux_, const view_1d_const &surf_evap_, const sview_2d_const &surf_mom_flux_, const view_3d_strided &qtracers_, const view_2d &qv_, const view_2d_const &qc_, const view_2d &qc_copy_, const view_2d &tke_, const view_2d &tke_copy_, const view_2d &z_mid_, const view_2d &z_int_, const view_2d &dse_, const view_2d &rrho_, const view_2d &rrho_i_, const view_2d &thv_, const view_2d &dz_, const view_2d &zt_grid_, const view_2d &zi_grid_, const view_1d &wpthlp_sfc_, const view_1d &wprtp_sfc_, const view_1d &upwp_sfc_, const view_1d &vpwp_sfc_, const view_2d &wtracer_sfc_, const view_2d &wm_zt_, const view_2d &inv_exner_, const view_2d &thlm_, const view_2d &qw_, const view_2d &cldfrac_liq_, const view_2d &cldfrac_liq_prev_)
 

Public Attributes

int ncol
 
int nlev
 
Real z_surf
 
view_2d_const T_mid
 
view_2d_const p_mid
 
view_2d_const p_int
 
view_2d_const pseudo_density
 
view_2d_const omega
 
view_1d_const phis
 
view_1d_const surf_sens_flux
 
view_1d_const surf_evap
 
sview_2d_const surf_mom_flux
 
view_3d_strided qtracers
 
view_2d qv
 
view_2d_const qc
 
view_2d qc_copy
 
view_2d z_mid
 
view_2d z_int
 
view_2d shoc_s
 
view_2d tke
 
view_2d tke_copy
 
view_2d rrho
 
view_2d rrho_i
 
view_2d thv
 
view_2d dz
 
view_2d zt_grid
 
view_2d zi_grid
 
view_1d wpthlp_sfc
 
view_1d wprtp_sfc
 
view_1d upwp_sfc
 
view_1d vpwp_sfc
 
view_2d wtracer_sfc
 
view_2d wm_zt
 
view_2d inv_exner
 
view_2d thlm
 
view_2d qw
 
view_2d cloud_frac
 
view_2d cldfrac_liq
 
view_2d cldfrac_liq_prev
 

Constructor & Destructor Documentation

◆ SHOCPreprocess()

SHOCInterface::SHOCPreprocess::SHOCPreprocess ( )
default

Member Function Documentation

◆ operator()()

KOKKOS_INLINE_FUNCTION void SHOCInterface::SHOCPreprocess::operator() ( const Kokkos::TeamPolicy< KT::ExeSpace >::member_type &  team) const
inline
124  {
125  const int i = team.league_rank();
126 
127  const Real zvir = C::ZVIR;
128  const Real cpair = C::Cpair;
129  const Real ggr = C::gravit;
130  const Real inv_ggr = 1/ggr;
131  const Real mintke = SC::mintke;
132 
133  const int nlev_packs = ekat::npack<Spack>(nlev);
134 
135  Kokkos::parallel_for(Kokkos::TeamVectorRange(team, nlev_packs), [&] (const Int& k)
136  {
137  cldfrac_liq_prev(i,k)=cldfrac_liq(i,k);
138 
139  // Inverse of Exner. In non-rel builds, assert that exner != 0 when in range before computing.
140  const Spack exner = PF::exner_function(p_mid(i,k));
141  const Smask nonzero = (exner != 0);
142  EKAT_KERNEL_ASSERT((nonzero || !(ekat::range<IntSmallPack>(k*Spack::n) < nlev)).all());
143  inv_exner(i,k).set(nonzero, 1/exner);
144 
145  tke(i,k) = ekat::max(mintke, tke(i,k));
146 
147  // Tracers are updated as a group. The tracers tke and qc act as separate inputs to shoc_main()
148  // and are therefore updated differently to the tracers group's monolithic field. Here, we make
149  // a copy if each of these tracers and pass to shoc_main() so that changes to the tracer group
150  // does not alter tke or qc values. Then during post processing, we copy back correct values of
151  // tke and qc to tracer group in postprocessing.
152  // TODO: remove *_copy views once SHOC can request a subset of tracers.
153  tke_copy(i,k) = tke(i,k);
154  qc_copy(i,k) = qc(i,k);
155 
156  qw(i,k) = qv(i,k) + qc(i,k);
157 
158  // Temperature
159  // NOTE: theta_v (thv) is intentionally different from one in HOMME
160  const auto theta_zt = PF::calculate_theta_from_T(T_mid(i,k),p_mid(i,k));
161  thlm(i,k) = PF::calculate_thetal_from_theta(theta_zt,T_mid(i,k),qc(i,k));
162  thv(i,k) = theta_zt*(1 + zvir*qv(i,k) - qc(i,k));
163 
164  // Vertical layer thickness
165  dz(i,k) = PF::calculate_dz(pseudo_density(i,k), p_mid(i,k), T_mid(i,k), qv(i,k));
166 
167  rrho(i,k) = inv_ggr*(pseudo_density(i,k)/dz(i,k));
168  wm_zt(i,k) = -1*omega(i,k)/(rrho(i,k)*ggr);
169  });
170  team.team_barrier();
171 
172  // Compute vertical layer heights
173  const auto dz_s = ekat::subview(dz, i);
174  const auto z_int_s = ekat::subview(z_int, i);
175  const auto z_mid_s = ekat::subview(z_mid, i);
176  PF::calculate_z_int(team,nlev,dz_s,z_surf,z_int_s);
177  team.team_barrier();
178  PF::calculate_z_mid(team,nlev,z_int_s,z_mid_s);
179  team.team_barrier();
180 
181  const int nlevi_v = nlev/Spack::n;
182  const int nlevi_p = nlev%Spack::n;
183  Kokkos::parallel_for(Kokkos::TeamVectorRange(team, nlev_packs), [&] (const Int& k)
184  {
185  zt_grid(i,k) = z_mid(i,k) - z_int(i, nlevi_v)[nlevi_p];
186  zi_grid(i,k) = z_int(i,k) - z_int(i, nlevi_v)[nlevi_p];
187 
188  // Dry static energy
189  shoc_s(i,k) = PF::calculate_dse(T_mid(i,k),z_mid(i,k),phis(i));
190 
191  if (k+1 == nlev_packs) zi_grid(i,nlevi_v)[nlevi_p] = 0;
192  });
193  team.team_barrier();
194 
195  const auto zt_grid_s = ekat::subview(zt_grid, i);
196  const auto zi_grid_s = ekat::subview(zi_grid, i);
197  const auto rrho_s = ekat::subview(rrho, i);
198  const auto rrho_i_s = ekat::subview(rrho_i, i);
199  SHF::linear_interp(team,zt_grid_s,zi_grid_s,rrho_s,rrho_i_s,nlev,nlev+1,0);
200  team.team_barrier();
201 
202  const auto exner_int = PF::exner_function(p_int(i,nlevi_v)[nlevi_p]);
203  const auto inv_exner_int_surf = 1/exner_int;
204 
205  wpthlp_sfc(i) = (surf_sens_flux(i)/(cpair*rrho_i(i,nlevi_v)[nlevi_p]))*inv_exner_int_surf;
206  wprtp_sfc(i) = surf_evap(i)/rrho_i(i,nlevi_v)[nlevi_p];
207  upwp_sfc(i) = surf_mom_flux(i,0)/rrho_i(i,nlevi_v)[nlevi_p];
208  vpwp_sfc(i) = surf_mom_flux(i,1)/rrho_i(i,nlevi_v)[nlevi_p];
209  } // operator
int Int
Definition: ERF_ShocInterface.H:20
amrex::Real Real
Definition: ERF_ShocInterface.H:19
typename SHF::Smask Smask
Definition: ERF_ShocInterface.H:45
typename SHF::Spack Spack
Definition: ERF_ShocInterface.H:43
int nlev
Definition: ERF_ShocInterface.H:212
view_2d cldfrac_liq_prev
Definition: ERF_ShocInterface.H:249
view_2d rrho
Definition: ERF_ShocInterface.H:232
view_2d_const pseudo_density
Definition: ERF_ShocInterface.H:217
sview_2d_const surf_mom_flux
Definition: ERF_ShocInterface.H:222
view_2d_const qc
Definition: ERF_ShocInterface.H:225
view_2d thlm
Definition: ERF_ShocInterface.H:245
view_1d_const phis
Definition: ERF_ShocInterface.H:219
view_1d_const surf_evap
Definition: ERF_ShocInterface.H:221
view_2d zt_grid
Definition: ERF_ShocInterface.H:236
view_2d cldfrac_liq
Definition: ERF_ShocInterface.H:248
view_2d inv_exner
Definition: ERF_ShocInterface.H:244
view_1d_const surf_sens_flux
Definition: ERF_ShocInterface.H:220
view_2d wm_zt
Definition: ERF_ShocInterface.H:243
view_2d qc_copy
Definition: ERF_ShocInterface.H:226
view_2d z_mid
Definition: ERF_ShocInterface.H:227
view_2d rrho_i
Definition: ERF_ShocInterface.H:233
view_2d_const omega
Definition: ERF_ShocInterface.H:218
view_2d_const T_mid
Definition: ERF_ShocInterface.H:214
Real z_surf
Definition: ERF_ShocInterface.H:213
view_2d_const p_int
Definition: ERF_ShocInterface.H:216
view_2d tke
Definition: ERF_ShocInterface.H:230
view_1d wprtp_sfc
Definition: ERF_ShocInterface.H:239
view_2d qw
Definition: ERF_ShocInterface.H:246
view_1d upwp_sfc
Definition: ERF_ShocInterface.H:240
view_1d wpthlp_sfc
Definition: ERF_ShocInterface.H:238
view_2d z_int
Definition: ERF_ShocInterface.H:228
view_2d_const p_mid
Definition: ERF_ShocInterface.H:215
view_1d vpwp_sfc
Definition: ERF_ShocInterface.H:241
view_2d dz
Definition: ERF_ShocInterface.H:235
view_2d shoc_s
Definition: ERF_ShocInterface.H:229
view_2d tke_copy
Definition: ERF_ShocInterface.H:231
view_2d zi_grid
Definition: ERF_ShocInterface.H:237
view_2d qv
Definition: ERF_ShocInterface.H:224
view_2d thv
Definition: ERF_ShocInterface.H:234

◆ set_variables()

void SHOCInterface::SHOCPreprocess::set_variables ( const int  ncol_,
const int  nlev_,
const Real  z_surf_,
const view_2d_const T_mid_,
const view_2d_const p_mid_,
const view_2d_const p_int_,
const view_2d_const pseudo_density_,
const view_2d_const omega_,
const view_1d_const phis_,
const view_1d_const surf_sens_flux_,
const view_1d_const surf_evap_,
const sview_2d_const surf_mom_flux_,
const view_3d_strided qtracers_,
const view_2d qv_,
const view_2d_const qc_,
const view_2d qc_copy_,
const view_2d tke_,
const view_2d tke_copy_,
const view_2d z_mid_,
const view_2d z_int_,
const view_2d dse_,
const view_2d rrho_,
const view_2d rrho_i_,
const view_2d thv_,
const view_2d dz_,
const view_2d zt_grid_,
const view_2d zi_grid_,
const view_1d wpthlp_sfc_,
const view_1d wprtp_sfc_,
const view_1d upwp_sfc_,
const view_1d vpwp_sfc_,
const view_2d wtracer_sfc_,
const view_2d wm_zt_,
const view_2d inv_exner_,
const view_2d thlm_,
const view_2d qw_,
const view_2d cldfrac_liq_,
const view_2d cldfrac_liq_prev_ 
)
inline
267  {
268  ncol = ncol_;
269  nlev = nlev_;
270  z_surf = z_surf_;
271  // IN
272  T_mid = T_mid_;
273  p_mid = p_mid_;
274  p_int = p_int_;
275  pseudo_density = pseudo_density_;
276  omega = omega_;
277  phis = phis_;
278  surf_sens_flux = surf_sens_flux_;
279  surf_evap = surf_evap_;
280  surf_mom_flux = surf_mom_flux_;
281  qv = qv_;
282  // OUT
283  qtracers = qtracers_;
284  qc = qc_;
285  qc_copy = qc_copy_;
286  shoc_s = dse_;
287  tke = tke_;
288  tke_copy = tke_copy_;
289  z_mid = z_mid_;
290  z_int = z_int_;
291  rrho = rrho_;
292  rrho_i = rrho_i_;
293  thv = thv_;
294  dz = dz_;
295  zt_grid = zt_grid_;
296  zi_grid = zi_grid_;
297  wpthlp_sfc = wpthlp_sfc_;
298  wprtp_sfc = wprtp_sfc_;
299  upwp_sfc = upwp_sfc_;
300  vpwp_sfc = vpwp_sfc_;
301  wtracer_sfc = wtracer_sfc_;
302  wm_zt = wm_zt_;
303  inv_exner = inv_exner_;
304  thlm = thlm_;
305  qw = qw_;
306  cldfrac_liq=cldfrac_liq_;
307  cldfrac_liq_prev=cldfrac_liq_prev_;
308  } // set_variables
int ncol
Definition: ERF_ShocInterface.H:212
view_3d_strided qtracers
Definition: ERF_ShocInterface.H:223
view_2d wtracer_sfc
Definition: ERF_ShocInterface.H:242

Member Data Documentation

◆ cldfrac_liq

view_2d SHOCInterface::SHOCPreprocess::cldfrac_liq

Referenced by operator()(), and set_variables().

◆ cldfrac_liq_prev

view_2d SHOCInterface::SHOCPreprocess::cldfrac_liq_prev

Referenced by operator()(), and set_variables().

◆ cloud_frac

view_2d SHOCInterface::SHOCPreprocess::cloud_frac

◆ dz

view_2d SHOCInterface::SHOCPreprocess::dz

Referenced by operator()(), and set_variables().

◆ inv_exner

view_2d SHOCInterface::SHOCPreprocess::inv_exner

Referenced by operator()(), and set_variables().

◆ ncol

int SHOCInterface::SHOCPreprocess::ncol

Referenced by set_variables().

◆ nlev

int SHOCInterface::SHOCPreprocess::nlev

Referenced by operator()(), and set_variables().

◆ omega

view_2d_const SHOCInterface::SHOCPreprocess::omega

Referenced by operator()(), and set_variables().

◆ p_int

view_2d_const SHOCInterface::SHOCPreprocess::p_int

Referenced by operator()(), and set_variables().

◆ p_mid

view_2d_const SHOCInterface::SHOCPreprocess::p_mid

Referenced by operator()(), and set_variables().

◆ phis

view_1d_const SHOCInterface::SHOCPreprocess::phis

Referenced by operator()(), and set_variables().

◆ pseudo_density

view_2d_const SHOCInterface::SHOCPreprocess::pseudo_density

Referenced by operator()(), and set_variables().

◆ qc

view_2d_const SHOCInterface::SHOCPreprocess::qc

Referenced by operator()(), and set_variables().

◆ qc_copy

view_2d SHOCInterface::SHOCPreprocess::qc_copy

Referenced by operator()(), and set_variables().

◆ qtracers

view_3d_strided SHOCInterface::SHOCPreprocess::qtracers

Referenced by set_variables().

◆ qv

view_2d SHOCInterface::SHOCPreprocess::qv

Referenced by operator()(), and set_variables().

◆ qw

view_2d SHOCInterface::SHOCPreprocess::qw

Referenced by operator()(), and set_variables().

◆ rrho

view_2d SHOCInterface::SHOCPreprocess::rrho

Referenced by operator()(), and set_variables().

◆ rrho_i

view_2d SHOCInterface::SHOCPreprocess::rrho_i

Referenced by operator()(), and set_variables().

◆ shoc_s

view_2d SHOCInterface::SHOCPreprocess::shoc_s

Referenced by operator()(), and set_variables().

◆ surf_evap

view_1d_const SHOCInterface::SHOCPreprocess::surf_evap

Referenced by operator()(), and set_variables().

◆ surf_mom_flux

sview_2d_const SHOCInterface::SHOCPreprocess::surf_mom_flux

Referenced by operator()(), and set_variables().

◆ surf_sens_flux

view_1d_const SHOCInterface::SHOCPreprocess::surf_sens_flux

Referenced by operator()(), and set_variables().

◆ T_mid

view_2d_const SHOCInterface::SHOCPreprocess::T_mid

Referenced by operator()(), and set_variables().

◆ thlm

view_2d SHOCInterface::SHOCPreprocess::thlm

Referenced by operator()(), and set_variables().

◆ thv

view_2d SHOCInterface::SHOCPreprocess::thv

Referenced by operator()(), and set_variables().

◆ tke

view_2d SHOCInterface::SHOCPreprocess::tke

Referenced by operator()(), and set_variables().

◆ tke_copy

view_2d SHOCInterface::SHOCPreprocess::tke_copy

Referenced by operator()(), and set_variables().

◆ upwp_sfc

view_1d SHOCInterface::SHOCPreprocess::upwp_sfc

Referenced by operator()(), and set_variables().

◆ vpwp_sfc

view_1d SHOCInterface::SHOCPreprocess::vpwp_sfc

Referenced by operator()(), and set_variables().

◆ wm_zt

view_2d SHOCInterface::SHOCPreprocess::wm_zt

Referenced by operator()(), and set_variables().

◆ wprtp_sfc

view_1d SHOCInterface::SHOCPreprocess::wprtp_sfc

Referenced by operator()(), and set_variables().

◆ wpthlp_sfc

view_1d SHOCInterface::SHOCPreprocess::wpthlp_sfc

Referenced by operator()(), and set_variables().

◆ wtracer_sfc

view_2d SHOCInterface::SHOCPreprocess::wtracer_sfc

Referenced by set_variables().

◆ z_int

view_2d SHOCInterface::SHOCPreprocess::z_int

Referenced by operator()(), and set_variables().

◆ z_mid

view_2d SHOCInterface::SHOCPreprocess::z_mid

Referenced by operator()(), and set_variables().

◆ z_surf

Real SHOCInterface::SHOCPreprocess::z_surf

Referenced by operator()(), and set_variables().

◆ zi_grid

view_2d SHOCInterface::SHOCPreprocess::zi_grid

Referenced by operator()(), and set_variables().

◆ zt_grid

view_2d SHOCInterface::SHOCPreprocess::zt_grid

Referenced by operator()(), and set_variables().


The documentation for this struct was generated from the following file: