ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
ERF_ShocInterface.H
Go to the documentation of this file.
1 #ifndef ERF_SHOCINTERFACE_H
2 #define ERF_SHOCINTERFACE_H
3 
4 /*
5  * SHOC model interface to ERF -- we use this interface to link to the SHOC code
6  * as cloned from the E3SM git repository: https://github.com/E3SM-Project/E3SM/components/eamxx/src/physics/shoc
7  */
8 
9 #include "physics/shoc/shoc_constants.hpp"
10 #include "physics/shoc/shoc_functions.hpp"
11 #include "share/util/eamxx_common_physics_functions.hpp"
12 
13 #include "ekat_subview_utils.hpp"
14 #include "ekat_parameter_list.hpp"
15 
17 using Int=int;
18 
19 #include <string>
20 
21 #include <Kokkos_Core.hpp>
22 #include <AMReX_ParmParse.H>
23 #include <AMReX_MultiFabUtil.H>
24 #include <ERF_DataStruct.H>
25 #include <ERF_Kokkos.H>
26 #include <ERF_IndexDefines.H>
27 #include <ERF_EOS.H>
28 
29 //#include "share/atm_process/atmosphere_process.hpp"
30 //#include "share/atm_process/ATMBufferManager.hpp"
31 
33 {
34  using SHF = scream::shoc::Functions<Real, KokkosDefaultDevice>;
35  using PF = scream::PhysicsFunctions<KokkosDefaultDevice>;
36  using C = scream::physics::Constants<Real>;
37  using KT = ekat::KokkosTypes<KokkosDefaultDevice>;
38  using SC = scream::shoc::Constants<Real>;
39 
40  using Spack = typename SHF::Spack;
41  using IntSmallPack = typename SHF::IntSmallPack;
42  using Smask = typename SHF::Smask;
43  using view_1d_int = typename KT::template view_1d<Int>;
44  using view_1d = typename SHF::view_1d<Real>;
45  using view_1d_const = typename SHF::view_1d<const Real>;
46  using view_2d = typename SHF::view_2d<SHF::Spack>;
47  using view_2d_const = typename SHF::view_2d<const Spack>;
48  using sview_2d = typename ekat::KokkosTypes<KokkosDefaultDevice>::template view_2d<Real>;
49  using sview_2d_const = typename ekat::KokkosTypes<KokkosDefaultDevice>::template view_2d<const Real>;
50  using view_3d = typename SHF::view_3d<Spack>;
51  using view_3d_const = typename SHF::view_3d<const Spack>;
52  using view_3d_strided = typename SHF::view_3d_strided<Spack>;
53 
54  using WSM = ekat::WorkspaceManager<Spack, KT::Device>;
55 
56  template<typename ScalarT>
57  using uview_1d = ekat::Unmanaged<typename KT::template view_1d<ScalarT>>;
58  template<typename ScalarT>
59  using uview_2d = ekat::Unmanaged<typename KT::template view_2d<ScalarT>>;
60 
61 public:
62 
63  // Constructor
64  SHOCInterface (const int& lev,
65  SolverChoice& sc);
66 
67  // Set the grid
68  void
69  set_grids (int& level,
70  const amrex::BoxArray& ba,
71  amrex::Geometry& geom,
72  amrex::MultiFab* cons_in,
73  amrex::MultiFab* z_phys,
74  amrex::MultiFab* lat,
75  amrex::MultiFab* lon);
76 
77  // Initialize the implementation
78  void
79  initialize_impl ();
80 
81  // Run the implementation
82  void
83  run_impl (const Real dt);
84 
85  // Finalize the implementation
86  void
87  finalize_impl ();
88 
89  // Allocate buffer space
90  void
91  alloc_buffers ();
92 
93  // De-allocate buffer space
94  void
95  dealloc_buffers ();
96 
97  // Fill KOKKOS Views from AMReX MultiFabs
98  void
100 
101  // Fill AMReX MultiFabs from KOKKOS Views
102  void
104 
105  // The name of the subcomponent
106  std::string name () const { return "shoc"; }
107 
108 #ifndef KOKKOS_ENABLE_CUDA
109  // Cuda requires methods enclosing __device__ lambda's to be public
110 protected:
111 #endif
112 
113  /*--------------------------------------------------------------------------------------------*/
114  // Most individual processes have a pre-processing step that constructs needed variables from
115  // the set of fields stored in the field manager. A structure like this defines those operations,
116  // which can then be called during run_impl in the main .cpp code.
117  // Structure to handle the local generation of data needed by shoc_main in run_impl
118  struct SHOCPreprocess {
119  SHOCPreprocess () = default;
120 
121  KOKKOS_INLINE_FUNCTION
122  void operator()(const Kokkos::TeamPolicy<KT::ExeSpace>::member_type& team) const
123  {
124  const int i = team.league_rank();
125 
126  const Real zvir = C::ZVIR;
127  const Real cpair = C::Cpair;
128  const Real ggr = C::gravit;
129  const Real inv_ggr = 1/ggr;
130  const Real mintke = SC::mintke;
131 
132  const int nlev_packs = ekat::npack<Spack>(nlev);
133 
134  Kokkos::parallel_for(Kokkos::TeamVectorRange(team, nlev_packs), [&] (const Int& k)
135  {
136  cldfrac_liq_prev(i,k)=cldfrac_liq(i,k);
137 
138  // Inverse of Exner. In non-rel builds, assert that exner != 0 when in range before computing.
139  const Spack exner = PF::exner_function(p_mid(i,k));
140  const Smask nonzero = (exner != 0);
141  EKAT_KERNEL_ASSERT((nonzero || !(ekat::range<IntSmallPack>(k*Spack::n) < nlev)).all());
142  inv_exner(i,k).set(nonzero, 1/exner);
143 
144  tke(i,k) = ekat::max(mintke, tke(i,k));
145 
146  // Tracers are updated as a group. The tracers tke and qc act as separate inputs to shoc_main()
147  // and are therefore updated differently to the tracers group's monolithic field. Here, we make
148  // a copy if each of these tracers and pass to shoc_main() so that changes to the tracer group
149  // does not alter tke or qc values. Then during post processing, we copy back correct values of
150  // tke and qc to tracer group in postprocessing.
151  // TODO: remove *_copy views once SHOC can request a subset of tracers.
152  tke_copy(i,k) = tke(i,k);
153  qc_copy(i,k) = qc(i,k);
154 
155  qw(i,k) = qv(i,k) + qc(i,k);
156 
157  // Temperature
158  // NOTE: theta_v (thv) is intentionally different from one in HOMME
159  const auto theta_zt = PF::calculate_theta_from_T(T_mid(i,k),p_mid(i,k));
160  thlm(i,k) = PF::calculate_thetal_from_theta(theta_zt,T_mid(i,k),qc(i,k));
161  thv(i,k) = theta_zt*(1 + zvir*qv(i,k) - qc(i,k));
162 
163  // Vertical layer thickness
164  dz(i,k) = PF::calculate_dz(pseudo_density(i,k), p_mid(i,k), T_mid(i,k), qv(i,k));
165 
166  rrho(i,k) = inv_ggr*(pseudo_density(i,k)/dz(i,k));
167  wm_zt(i,k) = -1*omega(i,k)/(rrho(i,k)*ggr);
168  });
169  team.team_barrier();
170 
171  // Compute vertical layer heights
172  const auto dz_s = ekat::subview(dz, i);
173  const auto z_int_s = ekat::subview(z_int, i);
174  const auto z_mid_s = ekat::subview(z_mid, i);
175  PF::calculate_z_int(team,nlev,dz_s,z_surf,z_int_s);
176  team.team_barrier();
177  PF::calculate_z_mid(team,nlev,z_int_s,z_mid_s);
178  team.team_barrier();
179 
180  const int nlevi_v = nlev/Spack::n;
181  const int nlevi_p = nlev%Spack::n;
182  Kokkos::parallel_for(Kokkos::TeamVectorRange(team, nlev_packs), [&] (const Int& k)
183  {
184  zt_grid(i,k) = z_mid(i,k) - z_int(i, nlevi_v)[nlevi_p];
185  zi_grid(i,k) = z_int(i,k) - z_int(i, nlevi_v)[nlevi_p];
186 
187  // Dry static energy
188  shoc_s(i,k) = PF::calculate_dse(T_mid(i,k),z_mid(i,k),phis(i));
189 
190  if (k+1 == nlev_packs) zi_grid(i,nlevi_v)[nlevi_p] = 0;
191  });
192  team.team_barrier();
193 
194  const auto zt_grid_s = ekat::subview(zt_grid, i);
195  const auto zi_grid_s = ekat::subview(zi_grid, i);
196  const auto rrho_s = ekat::subview(rrho, i);
197  const auto rrho_i_s = ekat::subview(rrho_i, i);
198  SHF::linear_interp(team,zt_grid_s,zi_grid_s,rrho_s,rrho_i_s,nlev,nlev+1,0);
199  team.team_barrier();
200 
201  const auto exner_int = PF::exner_function(p_int(i,nlevi_v)[nlevi_p]);
202  const auto inv_exner_int_surf = 1/exner_int;
203 
204  wpthlp_sfc(i) = (surf_sens_flux(i)/(cpair*rrho_i(i,nlevi_v)[nlevi_p]))*inv_exner_int_surf;
205  wprtp_sfc(i) = surf_evap(i)/rrho_i(i,nlevi_v)[nlevi_p];
206  upwp_sfc(i) = surf_mom_flux(i,0)/rrho_i(i,nlevi_v)[nlevi_p];
207  vpwp_sfc(i) = surf_mom_flux(i,1)/rrho_i(i,nlevi_v)[nlevi_p];
208  } // operator
209 
210  // Local variables
211  int ncol, nlev;
249 
250  // Assigning local variables
251  void set_variables(const int ncol_, const int nlev_,
252  const Real z_surf_,
253  const view_2d_const& T_mid_, const view_2d_const& p_mid_, const view_2d_const& p_int_, const view_2d_const& pseudo_density_,
254  const view_2d_const& omega_,
255  const view_1d_const& phis_, const view_1d_const& surf_sens_flux_, const view_1d_const& surf_evap_,
256  const sview_2d_const& surf_mom_flux_,
257  const view_3d_strided& qtracers_,
258  const view_2d& qv_, const view_2d_const& qc_, const view_2d& qc_copy_,
259  const view_2d& tke_, const view_2d& tke_copy_,
260  const view_2d& z_mid_, const view_2d& z_int_,
261  const view_2d& dse_, const view_2d& rrho_, const view_2d& rrho_i_,
262  const view_2d& thv_, const view_2d& dz_,const view_2d& zt_grid_,const view_2d& zi_grid_, const view_1d& wpthlp_sfc_,
263  const view_1d& wprtp_sfc_,const view_1d& upwp_sfc_,const view_1d& vpwp_sfc_, const view_2d& wtracer_sfc_,
264  const view_2d& wm_zt_,const view_2d& inv_exner_,const view_2d& thlm_,const view_2d& qw_,
265  const view_2d& cldfrac_liq_, const view_2d& cldfrac_liq_prev_)
266  {
267  ncol = ncol_;
268  nlev = nlev_;
269  z_surf = z_surf_;
270  // IN
271  T_mid = T_mid_;
272  p_mid = p_mid_;
273  p_int = p_int_;
274  pseudo_density = pseudo_density_;
275  omega = omega_;
276  phis = phis_;
277  surf_sens_flux = surf_sens_flux_;
278  surf_evap = surf_evap_;
279  surf_mom_flux = surf_mom_flux_;
280  qv = qv_;
281  // OUT
282  qtracers = qtracers_;
283  qc = qc_;
284  qc_copy = qc_copy_;
285  shoc_s = dse_;
286  tke = tke_;
287  tke_copy = tke_copy_;
288  z_mid = z_mid_;
289  z_int = z_int_;
290  rrho = rrho_;
291  rrho_i = rrho_i_;
292  thv = thv_;
293  dz = dz_;
294  zt_grid = zt_grid_;
295  zi_grid = zi_grid_;
296  wpthlp_sfc = wpthlp_sfc_;
297  wprtp_sfc = wprtp_sfc_;
298  upwp_sfc = upwp_sfc_;
299  vpwp_sfc = vpwp_sfc_;
300  wtracer_sfc = wtracer_sfc_;
301  wm_zt = wm_zt_;
302  inv_exner = inv_exner_;
303  thlm = thlm_;
304  qw = qw_;
305  cldfrac_liq=cldfrac_liq_;
306  cldfrac_liq_prev=cldfrac_liq_prev_;
307  } // set_variables
308  }; // SHOCPreprocess
309  /* --------------------------------------------------------------------------------------------*/
310 
311  /*--------------------------------------------------------------------------------------------*/
312  // Structure to handle the generation of data needed by the rest of the model based on output from
313  // shoc_main.
315  SHOCPostprocess () = default;
316 
317  KOKKOS_INLINE_FUNCTION
318  void operator()(const Kokkos::TeamPolicy<KT::ExeSpace>::member_type& team) const
319  {
320  const int i = team.league_rank();
321 
322  const Real inv_qc_relvar_max = 10;
323  const Real inv_qc_relvar_min = 0.001;
324 
325  const int nlev_packs = ekat::npack<Spack>(nlev);
326  Kokkos::parallel_for(Kokkos::TeamVectorRange(team, nlev_packs), [&] (const Int& k)
327  {
328  // See comment in SHOCPreprocess::operator() about the necessity of *_copy views
329  tke(i,k) = tke_copy(i,k);
330  qc(i,k) = qc_copy(i,k);
331 
332  qv(i,k) = qw(i,k) - qc(i,k);
333 
334  cldfrac_liq(i,k) = ekat::min(cldfrac_liq(i,k), 1);
335 
336  //P3 uses inv_qc_relvar, P3 is using dry mmrs, but
337  //wet<->dry conversion is a constant factor that cancels out in mean(qc)^2/mean(qc'*qc').
338  inv_qc_relvar(i,k) = 1;
339  const auto condition = (qc(i,k) != 0 && qc2(i,k) != 0);
340  if (condition.any()) {
341  inv_qc_relvar(i,k).set(condition,
342  ekat::min(inv_qc_relvar_max,
343  ekat::max(inv_qc_relvar_min,
344  ekat::square(qc(i,k))/qc2(i,k))));
345  }
346 
347  // Temperature
348  const Spack dse_ik(dse(i,k));
349  const Spack z_mid_ik(z_mid(i,k));
350  const Real phis_i(phis(i));
351  T_mid(i,k) = PF::calculate_temperature_from_dse(dse_ik,z_mid_ik,phis_i);
352 
353  });
354 
355  // If necessary, set appropriate boundary fluxes for energy and mass conservation checks.
356  // Any boundary fluxes not included in SHOC interface are set to 0.
358  vapor_flux(i) = surf_evap(i);
359  water_flux(i) = 0.0;
360  ice_flux(i) = 0.0;
361  heat_flux(i) = surf_sens_flux(i);
362  }
363  } // operator
364 
365  // Local variables
366  int ncol, nlev;
384 
385  // Assigning local variables
386  void set_variables (const int ncol_, const int nlev_,
387  const view_2d_const& rrho_,
388  const view_2d& qv_, const view_2d_const& qw_, const view_2d& qc_, const view_2d_const& qc_copy_,
389  const view_2d& tke_, const view_2d_const& tke_copy_, const view_3d_strided& qtracers_, const view_2d_const& qc2_,
390  const view_2d& cldfrac_liq_, const view_2d& inv_qc_relvar_,
391  const view_2d& T_mid_, const view_2d_const& dse_, const view_2d_const& z_mid_, const view_1d_const phis_)
392  {
393  ncol = ncol_;
394  nlev = nlev_;
395  rrho = rrho_;
396  qv = qv_;
397  qw = qw_;
398  qc = qc_;
399  qc_copy = qc_copy_;
400  tke = tke_;
401  tke_copy = tke_copy_;
402  qtracers = qtracers_;
403  qc2 = qc2_;
404  cldfrac_liq = cldfrac_liq_;
405  inv_qc_relvar = inv_qc_relvar_;
406  T_mid = T_mid_;
407  dse = dse_;
408  z_mid = z_mid_;
409  phis = phis_;
410  } // set_variables
411 
412  void set_mass_and_energy_fluxes (const view_1d_const& surf_evap_, const view_1d_const& surf_sens_flux_,
413  const view_1d& vapor_flux_, const view_1d& water_flux_,
414  const view_1d& ice_flux_, const view_1d& heat_flux_)
415  {
417  surf_evap = surf_evap_;
418  surf_sens_flux = surf_sens_flux_;
419  vapor_flux = vapor_flux_;
420  water_flux = water_flux_;
421  ice_flux = ice_flux_;
422  heat_flux = heat_flux_;
423  }
424  }; // SHOCPostprocess
425  /* --------------------------------------------------------------------------------------------*/
426 
427  // Structure for storing local variables initialized using the ATMBufferManager
428  struct Buffer {
429 #ifndef SCREAM_SHOC_SMALL_KERNELS
430  static constexpr int num_1d_scalar_ncol = 4;
431 #else
432  static constexpr int num_1d_scalar_ncol = 15;
433 #endif
434  static constexpr int num_1d_scalar_nlev = 1;
435 #ifndef SCREAM_SHOC_SMALL_KERNELS
436  static constexpr int num_2d_vector_mid = 19;
437  static constexpr int num_2d_vector_int = 12;
438 #else
439  static constexpr int num_2d_vector_mid = 23;
440  static constexpr int num_2d_vector_int = 13;
441 #endif
442  static constexpr int num_2d_vector_tr = 1;
443 
448 #ifdef SCREAM_SHOC_SMALL_KERNELS
449  uview_1d<Real> se_b;
450  uview_1d<Real> ke_b;
451  uview_1d<Real> wv_b;
452  uview_1d<Real> wl_b;
453  uview_1d<Real> se_a;
454  uview_1d<Real> ke_a;
455  uview_1d<Real> wv_a;
456  uview_1d<Real> wl_a;
457  uview_1d<Real> kbfs;
458  uview_1d<Real> ustar2;
459  uview_1d<Real> wstar;
460 #endif
462 
463  uview_2d<Spack> unused; // Placeholder for unused views
495 #ifdef SCREAM_SHOC_SMALL_KERNELS
496  uview_2d<Spack> rho_zt;
497  uview_2d<Spack> shoc_qv;
499  uview_2d<Spack> dz_zt;
500  uview_2d<Spack> dz_zi;
502 #endif
503 
505  }; // BUFFER
506  /* --------------------------------------------------------------------------------------------*/
507 
508 #ifndef KOKKOS_ENABLE_CUDA
509  // Cuda requires methods enclosing __device__ lambda's to be public
510 protected:
511 #endif
512 
513  // Update flux (if necessary)
514  void check_flux_state_consistency (const double dt);
515 
516  // Apply TMS drag coeff to shoc_main inputs (if necessary)
518 
519 protected:
520 
521  // SHOC updates the 'tracers' group.
523 
524  // Computes total number of bytes needed for local variables
526 
527  // Set local variables
528  void init_buffers ();
529 
530  // Offsets for MultiFab <-> KOKKOS transfer
531  amrex::Vector<int> m_col_offsets;
532 
533  // Keep track of field dimensions and other scalar values
534  // needed in shoc_main
542 
543  // Grid level
544  int m_lev;
545 
546  // Step number
547  int m_step;
548 
549  // Geometry at given level
550  amrex::Geometry m_geom;
551 
552  // Boxarray at given level
553  amrex::BoxArray m_ba;
554 
555  // Pointer to the CC conserved vars
556  amrex::MultiFab* m_cons_in = nullptr;
557 
558  // Pointer to the terrain heights
559  amrex::MultiFab* m_z_phys = nullptr;
560 
561  // Flag for first step
562  bool m_first_step = true;
563 
564  // Struct which contains local variables
566 
567  // Store the structures for each argument to shoc_main;
568  SHF::SHOCInput input;
569  SHF::SHOCInputOutput input_output;
570  SHF::SHOCOutput output;
571  SHF::SHOCHistoryOutput history_output;
572  SHF::SHOCRuntime runtime_options;
573 #ifdef SCREAM_SHOC_SMALL_KERNELS
574  SHF::SHOCTemporaries temporaries;
575 #endif
576 
577  // Structures which compute pre/post process
580 
581  // WSM for internal local variables
582  ekat::WorkspaceManager<Spack, KT::Device> workspace_mgr;
583 
584  // SHOCPreprocess data structures
585  //=======================================================
594 
623 
625 
626  // SHOCPostprocess data structures
627  //=======================================================
630 
631  // SHOC InputOutput data structures
632  //=======================================================
636 
637  // SHOC Output data structures
638  //=======================================================
643 
644  // SHOC HistoryOutput data structures
645  //=======================================================
662 
663  // SHOC Miscellaneous data structures
664  //=======================================================
665 
666 }; // SHOCInterface
667 #endif
int Int
Definition: ERF_ShocInterface.H:17
amrex::Real Real
Definition: ERF_ShocInterface.H:16
Definition: ERF_ShocInterface.H:33
view_2d p_int
Definition: ERF_ShocInterface.H:597
typename SHF::view_1d< const Real > view_1d_const
Definition: ERF_ShocInterface.H:45
SHOCPreprocess shoc_preprocess
Definition: ERF_ShocInterface.H:578
view_2d wqw_sec
Definition: ERF_ShocInterface.H:652
scream::shoc::Constants< Real > SC
Definition: ERF_ShocInterface.H:38
view_2d inv_exner
Definition: ERF_ShocInterface.H:617
amrex::Geometry m_geom
Definition: ERF_ShocInterface.H:550
void set_computed_group_impl()
void alloc_buffers()
Definition: ERF_ShocInterface.cpp:106
view_2d wtke_sec
Definition: ERF_ShocInterface.H:653
view_2d thl_sec
Definition: ERF_ShocInterface.H:648
view_2d tkh
Definition: ERF_ShocInterface.H:642
view_2d cldfrac_liq
Definition: ERF_ShocInterface.H:621
amrex::MultiFab * m_z_phys
Definition: ERF_ShocInterface.H:559
void init_buffers()
Int m_num_layers
Definition: ERF_ShocInterface.H:536
view_2d vw_sec
Definition: ERF_ShocInterface.H:655
view_2d thv
Definition: ERF_ShocInterface.H:610
size_t requested_buffer_size_in_bytes() const
view_2d shoc_cond
Definition: ERF_ShocInterface.H:660
typename ekat::KokkosTypes< KokkosDefaultDevice >::template view_2d< Real > sview_2d
Definition: ERF_ShocInterface.H:48
view_1d vpwp_sfc
Definition: ERF_ShocInterface.H:590
view_2d brunt
Definition: ERF_ShocInterface.H:658
Int m_num_vel_comp
Definition: ERF_ShocInterface.H:540
view_2d dens
Definition: ERF_ShocInterface.H:598
view_1d surf_sens_flux
Definition: ERF_ShocInterface.H:591
view_1d pblh
Definition: ERF_ShocInterface.H:639
Int m_num_cols
Definition: ERF_ShocInterface.H:535
view_2d sgs_buoy_flux
Definition: ERF_ShocInterface.H:633
SHF::SHOCOutput output
Definition: ERF_ShocInterface.H:570
typename KT::template view_1d< Int > view_1d_int
Definition: ERF_ShocInterface.H:43
ekat::WorkspaceManager< Spack, KT::Device > workspace_mgr
Definition: ERF_ShocInterface.H:582
view_2d p_mid
Definition: ERF_ShocInterface.H:596
view_2d uw_sec
Definition: ERF_ShocInterface.H:654
view_2d dse
Definition: ERF_ShocInterface.H:612
Int m_nadv
Definition: ERF_ShocInterface.H:538
view_2d qc
Definition: ERF_ShocInterface.H:601
view_2d dz
Definition: ERF_ShocInterface.H:611
view_2d qwthl_sec
Definition: ERF_ShocInterface.H:650
void run_impl(const Real dt)
Definition: ERF_ShocInterface.cpp:772
Int m_npbl
Definition: ERF_ShocInterface.H:537
view_2d shoc_mix
Definition: ERF_ShocInterface.H:646
view_2d qw_sec
Definition: ERF_ShocInterface.H:649
view_2d cldfrac_liq_prev
Definition: ERF_ShocInterface.H:622
scream::physics::Constants< Real > C
Definition: ERF_ShocInterface.H:36
view_2d qw
Definition: ERF_ShocInterface.H:619
SHF::SHOCInput input
Definition: ERF_ShocInterface.H:568
SHOCInterface(const int &lev, SolverChoice &sc)
Definition: ERF_ShocInterface.cpp:6
view_2d cloud_frac
Definition: ERF_ShocInterface.H:620
view_2d omega
Definition: ERF_ShocInterface.H:599
typename SHF::view_2d< const Spack > view_2d_const
Definition: ERF_ShocInterface.H:47
ekat::Unmanaged< typename KT::template view_1d< ScalarT > > uview_1d
Definition: ERF_ShocInterface.H:57
typename SHF::Smask Smask
Definition: ERF_ShocInterface.H:42
typename SHF::view_3d_strided< Spack > view_3d_strided
Definition: ERF_ShocInterface.H:52
view_2d rrho_i
Definition: ERF_ShocInterface.H:609
view_1d obklen
Definition: ERF_ShocInterface.H:641
void kokkos_buffers_to_mf()
Definition: ERF_ShocInterface.cpp:447
view_2d z_mid
Definition: ERF_ShocInterface.H:603
view_2d qc_copy
Definition: ERF_ShocInterface.H:602
scream::shoc::Functions< Real, KokkosDefaultDevice > SHF
Definition: ERF_ShocInterface.H:34
SHF::SHOCRuntime runtime_options
Definition: ERF_ShocInterface.H:572
view_2d rrho
Definition: ERF_ShocInterface.H:608
view_2d inv_qc_relvar
Definition: ERF_ShocInterface.H:629
view_1d surf_evap
Definition: ERF_ShocInterface.H:592
view_2d w_sec
Definition: ERF_ShocInterface.H:647
view_1d phis
Definition: ERF_ShocInterface.H:586
view_2d shoc_evap
Definition: ERF_ShocInterface.H:661
typename ekat::KokkosTypes< KokkosDefaultDevice >::template view_2d< const Real > sview_2d_const
Definition: ERF_ShocInterface.H:49
view_3d horiz_wind
Definition: ERF_ShocInterface.H:635
void check_flux_state_consistency(const double dt)
typename SHF::view_1d< Real > view_1d
Definition: ERF_ShocInterface.H:44
view_2d wqls_sec
Definition: ERF_ShocInterface.H:657
view_3d_strided qtracers
Definition: ERF_ShocInterface.H:624
ekat::KokkosTypes< KokkosDefaultDevice > KT
Definition: ERF_ShocInterface.H:37
view_1d upwp_sfc
Definition: ERF_ShocInterface.H:589
view_1d ustar
Definition: ERF_ShocInterface.H:640
view_2d tke_copy
Definition: ERF_ShocInterface.H:607
view_2d wm_zt
Definition: ERF_ShocInterface.H:616
scream::PhysicsFunctions< KokkosDefaultDevice > PF
Definition: ERF_ShocInterface.H:35
void finalize_impl()
Definition: ERF_ShocInterface.cpp:863
view_2d w3
Definition: ERF_ShocInterface.H:656
void initialize_impl()
Definition: ERF_ShocInterface.cpp:452
ekat::Unmanaged< typename KT::template view_2d< ScalarT > > uview_2d
Definition: ERF_ShocInterface.H:59
ekat::WorkspaceManager< Spack, KT::Device > WSM
Definition: ERF_ShocInterface.H:54
typename SHF::view_3d< const Spack > view_3d_const
Definition: ERF_ShocInterface.H:51
void set_grids(int &level, const amrex::BoxArray &ba, amrex::Geometry &geom, amrex::MultiFab *cons_in, amrex::MultiFab *z_phys, amrex::MultiFab *lat, amrex::MultiFab *lon)
Definition: ERF_ShocInterface.cpp:62
bool m_first_step
Definition: ERF_ShocInterface.H:562
Int m_num_tracers
Definition: ERF_ShocInterface.H:539
std::string name() const
Definition: ERF_ShocInterface.H:106
int m_step
Definition: ERF_ShocInterface.H:547
sview_2d surf_mom_flux
Definition: ERF_ShocInterface.H:593
typename SHF::view_2d< SHF::Spack > view_2d
Definition: ERF_ShocInterface.H:46
view_2d zi_grid
Definition: ERF_ShocInterface.H:614
view_2d thlm
Definition: ERF_ShocInterface.H:618
int m_lev
Definition: ERF_ShocInterface.H:544
view_2d wtracer_sfc
Definition: ERF_ShocInterface.H:615
SHOCPostprocess shoc_postprocess
Definition: ERF_ShocInterface.H:579
typename SHF::Spack Spack
Definition: ERF_ShocInterface.H:40
SHF::SHOCInputOutput input_output
Definition: ERF_ShocInterface.H:569
Buffer m_buffer
Definition: ERF_ShocInterface.H:565
void mf_to_kokkos_buffers()
Definition: ERF_ShocInterface.cpp:290
view_2d shoc_ql2
Definition: ERF_ShocInterface.H:628
view_2d tk
Definition: ERF_ShocInterface.H:634
view_2d qv
Definition: ERF_ShocInterface.H:600
view_2d wthl_sec
Definition: ERF_ShocInterface.H:651
amrex::BoxArray m_ba
Definition: ERF_ShocInterface.H:553
view_2d tke
Definition: ERF_ShocInterface.H:606
view_2d zt_grid
Definition: ERF_ShocInterface.H:613
view_2d T_mid
Definition: ERF_ShocInterface.H:595
view_2d shoc_s
Definition: ERF_ShocInterface.H:605
void dealloc_buffers()
Definition: ERF_ShocInterface.cpp:201
Int hdtime
Definition: ERF_ShocInterface.H:541
view_2d z_int
Definition: ERF_ShocInterface.H:604
void apply_turbulent_mountain_stress()
typename SHF::view_3d< Spack > view_3d
Definition: ERF_ShocInterface.H:50
view_2d isotropy
Definition: ERF_ShocInterface.H:659
typename SHF::IntSmallPack IntSmallPack
Definition: ERF_ShocInterface.H:41
amrex::Vector< int > m_col_offsets
Definition: ERF_ShocInterface.H:531
SHF::SHOCHistoryOutput history_output
Definition: ERF_ShocInterface.H:571
view_1d wpthlp_sfc
Definition: ERF_ShocInterface.H:587
view_1d wprtp_sfc
Definition: ERF_ShocInterface.H:588
amrex::MultiFab * m_cons_in
Definition: ERF_ShocInterface.H:556
@ tabs
Definition: ERF_Kessler.H:24
Definition: ERF_ShocInterface.H:428
uview_2d< Spack > w3
Definition: ERF_ShocInterface.H:492
uview_2d< Spack > qc_copy
Definition: ERF_ShocInterface.H:479
uview_2d< Spack > dse
Definition: ERF_ShocInterface.H:477
uview_1d< Real > upwp_sfc
Definition: ERF_ShocInterface.H:446
uview_2d< Spack > vw_sec
Definition: ERF_ShocInterface.H:491
uview_1d< Real > wpthlp_sfc
Definition: ERF_ShocInterface.H:444
uview_2d< Spack > thv
Definition: ERF_ShocInterface.H:468
uview_2d< Spack > wqls_sec
Definition: ERF_ShocInterface.H:493
uview_2d< Spack > wthl_sec
Definition: ERF_ShocInterface.H:487
uview_2d< Spack > thlm
Definition: ERF_ShocInterface.H:475
uview_2d< Spack > uw_sec
Definition: ERF_ShocInterface.H:490
uview_2d< Spack > z_mid
Definition: ERF_ShocInterface.H:464
uview_1d< Real > vpwp_sfc
Definition: ERF_ShocInterface.H:447
uview_2d< Spack > rrho
Definition: ERF_ShocInterface.H:466
uview_2d< Spack > wtke_sec
Definition: ERF_ShocInterface.H:489
uview_2d< Spack > rrho_i
Definition: ERF_ShocInterface.H:467
uview_2d< Spack > unused
Definition: ERF_ShocInterface.H:463
Spack * wsm_data
Definition: ERF_ShocInterface.H:504
static constexpr int num_2d_vector_int
Definition: ERF_ShocInterface.H:437
uview_2d< Spack > wtracer_sfc
Definition: ERF_ShocInterface.H:472
static constexpr int num_2d_vector_tr
Definition: ERF_ShocInterface.H:442
static constexpr int num_1d_scalar_nlev
Definition: ERF_ShocInterface.H:434
uview_2d< Spack > inv_exner
Definition: ERF_ShocInterface.H:474
uview_2d< Spack > qwthl_sec
Definition: ERF_ShocInterface.H:486
uview_2d< Spack > qw
Definition: ERF_ShocInterface.H:476
uview_2d< Spack > brunt
Definition: ERF_ShocInterface.H:494
uview_2d< Spack > shoc_mix
Definition: ERF_ShocInterface.H:481
uview_2d< Spack > tke_copy
Definition: ERF_ShocInterface.H:478
static constexpr int num_1d_scalar_ncol
Definition: ERF_ShocInterface.H:430
uview_2d< Spack > isotropy
Definition: ERF_ShocInterface.H:482
uview_2d< Spack > zt_grid
Definition: ERF_ShocInterface.H:470
uview_2d< Spack > w_sec
Definition: ERF_ShocInterface.H:483
uview_2d< Spack > thl_sec
Definition: ERF_ShocInterface.H:484
uview_2d< Spack > wm_zt
Definition: ERF_ShocInterface.H:473
uview_2d< Spack > qw_sec
Definition: ERF_ShocInterface.H:485
uview_2d< Spack > wqw_sec
Definition: ERF_ShocInterface.H:488
uview_2d< Spack > z_int
Definition: ERF_ShocInterface.H:465
uview_2d< Spack > shoc_ql2
Definition: ERF_ShocInterface.H:480
uview_1d< Real > wprtp_sfc
Definition: ERF_ShocInterface.H:445
uview_2d< Spack > zi_grid
Definition: ERF_ShocInterface.H:471
static constexpr int num_2d_vector_mid
Definition: ERF_ShocInterface.H:436
uview_1d< Spack > pref_mid
Definition: ERF_ShocInterface.H:461
uview_2d< Spack > dz
Definition: ERF_ShocInterface.H:469
Definition: ERF_ShocInterface.H:314
view_1d_const surf_sens_flux
Definition: ERF_ShocInterface.H:379
view_1d ice_flux
Definition: ERF_ShocInterface.H:382
view_1d vapor_flux
Definition: ERF_ShocInterface.H:380
view_1d_const phis
Definition: ERF_ShocInterface.H:376
int nlev
Definition: ERF_ShocInterface.H:366
view_1d heat_flux
Definition: ERF_ShocInterface.H:383
view_2d cldfrac_liq
Definition: ERF_ShocInterface.H:372
view_2d_const rrho
Definition: ERF_ShocInterface.H:367
KOKKOS_INLINE_FUNCTION void operator()(const Kokkos::TeamPolicy< KT::ExeSpace >::member_type &team) const
Definition: ERF_ShocInterface.H:318
view_3d_strided qtracers
Definition: ERF_ShocInterface.H:370
void set_variables(const int ncol_, const int nlev_, const view_2d_const &rrho_, const view_2d &qv_, const view_2d_const &qw_, const view_2d &qc_, const view_2d_const &qc_copy_, const view_2d &tke_, const view_2d_const &tke_copy_, const view_3d_strided &qtracers_, const view_2d_const &qc2_, const view_2d &cldfrac_liq_, const view_2d &inv_qc_relvar_, const view_2d &T_mid_, const view_2d_const &dse_, const view_2d_const &z_mid_, const view_1d_const phis_)
Definition: ERF_ShocInterface.H:386
view_2d_const qw
Definition: ERF_ShocInterface.H:369
bool compute_mass_and_energy_fluxes
Definition: ERF_ShocInterface.H:377
view_1d water_flux
Definition: ERF_ShocInterface.H:381
view_1d_const surf_evap
Definition: ERF_ShocInterface.H:378
view_2d_const qc_copy
Definition: ERF_ShocInterface.H:369
view_2d_const z_mid
Definition: ERF_ShocInterface.H:375
view_2d_const dse
Definition: ERF_ShocInterface.H:375
view_2d tke
Definition: ERF_ShocInterface.H:368
view_2d T_mid
Definition: ERF_ShocInterface.H:374
view_2d inv_qc_relvar
Definition: ERF_ShocInterface.H:373
view_2d qc
Definition: ERF_ShocInterface.H:368
void set_mass_and_energy_fluxes(const view_1d_const &surf_evap_, const view_1d_const &surf_sens_flux_, const view_1d &vapor_flux_, const view_1d &water_flux_, const view_1d &ice_flux_, const view_1d &heat_flux_)
Definition: ERF_ShocInterface.H:412
view_2d qv
Definition: ERF_ShocInterface.H:368
int ncol
Definition: ERF_ShocInterface.H:366
view_2d_const qc2
Definition: ERF_ShocInterface.H:371
view_2d_const tke_copy
Definition: ERF_ShocInterface.H:369
Definition: ERF_ShocInterface.H:118
int ncol
Definition: ERF_ShocInterface.H:211
view_3d_strided qtracers
Definition: ERF_ShocInterface.H:222
int nlev
Definition: ERF_ShocInterface.H:211
view_2d cldfrac_liq_prev
Definition: ERF_ShocInterface.H:248
view_2d rrho
Definition: ERF_ShocInterface.H:231
view_2d_const pseudo_density
Definition: ERF_ShocInterface.H:216
sview_2d_const surf_mom_flux
Definition: ERF_ShocInterface.H:221
view_2d_const qc
Definition: ERF_ShocInterface.H:224
view_2d thlm
Definition: ERF_ShocInterface.H:244
view_1d_const phis
Definition: ERF_ShocInterface.H:218
view_1d_const surf_evap
Definition: ERF_ShocInterface.H:220
view_2d zt_grid
Definition: ERF_ShocInterface.H:235
view_2d cldfrac_liq
Definition: ERF_ShocInterface.H:247
view_2d inv_exner
Definition: ERF_ShocInterface.H:243
view_1d_const surf_sens_flux
Definition: ERF_ShocInterface.H:219
view_2d wm_zt
Definition: ERF_ShocInterface.H:242
KOKKOS_INLINE_FUNCTION void operator()(const Kokkos::TeamPolicy< KT::ExeSpace >::member_type &team) const
Definition: ERF_ShocInterface.H:122
view_2d qc_copy
Definition: ERF_ShocInterface.H:225
view_2d wtracer_sfc
Definition: ERF_ShocInterface.H:241
view_2d z_mid
Definition: ERF_ShocInterface.H:226
view_2d rrho_i
Definition: ERF_ShocInterface.H:232
view_2d_const omega
Definition: ERF_ShocInterface.H:217
view_2d_const T_mid
Definition: ERF_ShocInterface.H:213
view_2d cloud_frac
Definition: ERF_ShocInterface.H:246
Real z_surf
Definition: ERF_ShocInterface.H:212
view_2d_const p_int
Definition: ERF_ShocInterface.H:215
view_2d tke
Definition: ERF_ShocInterface.H:229
view_1d wprtp_sfc
Definition: ERF_ShocInterface.H:238
view_2d qw
Definition: ERF_ShocInterface.H:245
view_1d upwp_sfc
Definition: ERF_ShocInterface.H:239
view_1d wpthlp_sfc
Definition: ERF_ShocInterface.H:237
view_2d z_int
Definition: ERF_ShocInterface.H:227
view_2d_const p_mid
Definition: ERF_ShocInterface.H:214
view_1d vpwp_sfc
Definition: ERF_ShocInterface.H:240
view_2d dz
Definition: ERF_ShocInterface.H:234
view_2d shoc_s
Definition: ERF_ShocInterface.H:228
view_2d tke_copy
Definition: ERF_ShocInterface.H:230
view_2d zi_grid
Definition: ERF_ShocInterface.H:236
view_2d qv
Definition: ERF_ShocInterface.H:223
view_2d thv
Definition: ERF_ShocInterface.H:233
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_)
Definition: ERF_ShocInterface.H:251
Definition: ERF_DataStruct.H:123