ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
ERF_InitFromInputSounding.cpp File Reference
#include <ERF.H>
#include <ERF_EOS.H>
#include <ERF_Constants.H>
#include <ERF_Utils.H>
#include <ERF_ProbCommon.H>
Include dependency graph for ERF_InitFromInputSounding.cpp:

Functions

void init_state_from_input_sounding (const Box &bx, Array4< Real > const &state, GeometryData const &geomdata, Array4< const Real > const &z_cc_arr, const bool &l_moist, InputSoundingData const &inputSoundingData)
 
void init_state_from_input_sounding_hse (const Box &bx, Array4< Real > const &state, Array4< Real > const &r_hse_arr, Array4< Real > const &p_hse_arr, Array4< Real > const &pi_hse_arr, Array4< Real > const &th_hse_arr, Array4< Real > const &qv_hse_arr, GeometryData const &geomdata, Array4< const Real > const &z_cc_arr, const Real &l_gravity, const Real &l_rdOcp, const bool &l_moist, InputSoundingData const &inputSoundingData, const bool &l_isentropic, const int &ngz)
 
void init_velocities_from_input_sounding (const Box &bx, Array4< Real > const &x_vel, Array4< Real > const &y_vel, Array4< Real > const &z_vel, GeometryData const &geomdata, Array4< const Real > const &z_nd_arr, InputSoundingData const &inputSoundingData)
 

Function Documentation

◆ init_state_from_input_sounding()

void init_state_from_input_sounding ( const Box &  bx,
Array4< Real > const &  state,
GeometryData const &  geomdata,
Array4< const Real > const &  z_cc_arr,
const bool &  l_moist,
InputSoundingData const &  inputSoundingData 
)

Box level wrapper for initializing scalar data from input sounding data.

Parameters
[in]bxBox specifying the indices we are initializing
[out]stateArray4 specifying the state data we are to initialize
[in]geomdataGeometryData object specifying the domain geometry
[in]z_cc_arrArray4 specifying cell-centered terrain heights, if present
[in]l_moistWhether moisture variables should be initialized
[in]inputSoundingDataInputSoundingData object we are to initialize from

Box level wrapper for initializing scalar data from input sounding data.

Parameters
bxBox specifying the indices we are initializing
stateArray4 specifying the state data we are to initialize
geomdataGeometryData object specifying the domain geometry
z_cc_arrArray4 specifying cell-centered terrain heights, if present
l_moistWhether moisture variables should be initialized
inputSoundingDataInputSoundingData object we are to initialize from
327 {
328  const Real* z_inp_sound = inputSoundingData.z_inp_sound_d[0].dataPtr();
329  const Real* theta_inp_sound = inputSoundingData.theta_inp_sound_d[0].dataPtr();
330  const Real* qv_inp_sound = inputSoundingData.qv_inp_sound_d[0].dataPtr();
331  const int inp_sound_size = inputSoundingData.size(0);
332 
333  // Geometry
334  const Real* prob_lo = geomdata.ProbLo();
335  const Real* dx = geomdata.CellSize();
336  const Real z_lo = prob_lo[2];
337  const Real dz = dx[2];
338 
339  // We want to set the lateral BC values, too
340  Box gbx = bx; // Copy constructor
341  gbx.grow(0,1); gbx.grow(1,1); // Grow by one in the lateral directions
342 
343  ParallelFor(gbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
344  const Real z = (z_cc_arr) ? z_cc_arr(i,j,k) : z_lo + (k + myhalf) * dz;
345 
346  Real rho_0 = one;
347 
348  // Set the density
349  state(i, j, k, Rho_comp) = rho_0;
350 
351  // Initial Rho0*Theta0
352  state(i, j, k, RhoTheta_comp) = rho_0 * interpolate_1d(z_inp_sound, theta_inp_sound, z, inp_sound_size);
353 
354  // total nonprecipitating water (Q1) == water vapor (Qv), i.e., there is no cloud water or cloud ice
355  if (l_moist) {
356  state(i, j, k, RhoQ1_comp) = rho_0 * interpolate_1d(z_inp_sound, qv_inp_sound, z, inp_sound_size);
357  }
358  });
359 }
#define Rho_comp
Definition: ERF_IndexDefines.H:42
#define RhoTheta_comp
Definition: ERF_IndexDefines.H:43
#define RhoQ1_comp
Definition: ERF_IndexDefines.H:48
Real rho_0
Definition: ERF_InitCustomPert_ABL.H:18
const Real dx
Definition: ERF_InitCustomPert_ABL.H:44
const GpuArray< Real, AMREX_SPACEDIM > prob_lo
Definition: ERF_InitCustomPert_CloudChamber.H:33
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real interpolate_1d(const amrex::Real *alpha, const amrex::Real *beta, const amrex::Real alpha_interp, const int alpha_size)
Definition: ERF_Interpolation_1D.H:14
ParallelFor(fab_box, [=] AMREX_GPU_DEVICE(int i, int j, int k) { qrcuten_arr(i, j, k)=Real(0);qscuten_arr(i, j, k)=Real(0);qicuten_arr(i, j, k)=Real(0);})
constexpr amrex::Real one
Definition: ERF_NumericalConstants.H:37
constexpr amrex::Real myhalf
Definition: ERF_NumericalConstants.H:41
amrex::Real Real
Definition: ERF_ShocInterface.H:19
@ dz
Definition: ERF_AdvanceWDM6.cpp:272

Referenced by ERF::init_from_input_sounding().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ init_state_from_input_sounding_hse()

void init_state_from_input_sounding_hse ( const Box &  bx,
Array4< Real > const &  state,
Array4< Real > const &  r_hse_arr,
Array4< Real > const &  p_hse_arr,
Array4< Real > const &  pi_hse_arr,
Array4< Real > const &  th_hse_arr,
Array4< Real > const &  qv_hse_arr,
GeometryData const &  geomdata,
Array4< const Real > const &  z_cc_arr,
const Real &  l_gravity,
const Real &  l_rdOcp,
const bool &  l_moist,
InputSoundingData const &  inputSoundingData,
const bool &  l_isentropic,
const int &  ngz 
)

Box level wrapper for initializing scalar and hydrostatic base state data from input sounding data.

Parameters
[in]bxBox specifying the indices we are initializing
[out]stateArray4 specifying the state data we are to initialize
[out]r_hse_arrArray4 specifying the density HSE base state data we are to initialize
[out]p_hse_arrArray4 specifying the pressure HSE base state data we are to initialize
[out]pi_hse_arrArray4 specifying the Exner pressure HSE base state data we are to initialize
[out]th_hse_arrArray4 specifying the base state potential temperature we are to initialize
[out]qv_hse_arrArray4 specifying the base state water vapor mixing ratio we are to initialize
[in]geomdataGeometryData object specifying the domain geometry
[in]z_cc_arrArray4 specifying cell-centered terrain heights, if present
[in]l_gravityGravity constant (retained for interface compatibility).
[in]l_rdOcpReal number specifying the Rhydberg constant ($R_d$) divided by specific heat at constant pressure ($c_p$)
[in]l_moistWhether moisture variables should be initialized
[in]inputSoundingDataInputSoundingData object we are to initialize from
[in]l_isentropicWhether to construct an isentropic HSE base state
[in]ngzNumber of vertical ghost cells to fill in HSE arrays

Box level wrapper for initializing scalar and hydrostatic base state data from input sounding data.

Parameters
bxBox specifying the indices we are initializing
stateArray4 specifying the state data we are to initialize
r_hse_arrArray4 specifying the density HSE base state data we are to initialize
p_hse_arrArray4 specifying the pressure HSE base state data we are to initialize
pi_hse_arrArray4 specifying the Exner pressure HSE base state data we are to initialize
th_hse_arrArray4 specifying the base state potential temperature we are to initialize
qv_hse_arrArray4 specifying the base state water vapor mixing ratio we are to initialize
geomdataGeometryData object specifying the domain geometry The unused gravity argument is retained for interface compatibility.
z_cc_arrArray4 specifying cell-centered terrain heights, if present
l_rdOcpReal number specifying the Rhydberg constant ($R_d$) divided by specific heat at constant pressure ($c_p$)
l_moistWhether moisture variables should be initialized
inputSoundingDataInputSoundingData object we are to initialize from
l_isentropicWhether to construct an isentropic HSE base state
ngzNumber of vertical ghost cells to fill in HSE arrays
397 {
398  const Real* z_inp_sound = inputSoundingData.z_inp_sound_d[0].dataPtr();
399  const Real* rho_inp_sound = inputSoundingData.rho_inp_sound_d.dataPtr();
400  const Real* theta_inp_sound = inputSoundingData.theta_inp_sound_d[0].dataPtr();
401  const Real* qv_inp_sound = inputSoundingData.qv_inp_sound_d[0].dataPtr();
402  const int inp_sound_size = inputSoundingData.size(0);
403  const bool anel_assume_dry = inputSoundingData.assume_dry;
404 
405  // Geometry
406  const Real* prob_lo = geomdata.ProbLo();
407  const Real* dx = geomdata.CellSize();
408  const Real z_lo = prob_lo[2];
409  const Real dz = dx[2];
410 
411  int kbot = geomdata.Domain().smallEnd(2);
412  int ktop = geomdata.Domain().bigEnd(2);
413 
414  // We want to set the lateral BC values, too
415  Box gbx = bx; // Copy constructor
416  gbx.grow(0,1); gbx.grow(1,1); // Grow by one in the lateral directions
417 
418  ParallelFor(gbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
419  const Real z = (z_cc_arr) ? z_cc_arr(i,j,k)
420  : z_lo + (k + myhalf) * dz;
421 
422  Real rho_k = interpolate_1d(z_inp_sound, rho_inp_sound, z, inp_sound_size);
423  Real rhoTh_k = rho_k * interpolate_1d(z_inp_sound, theta_inp_sound, z, inp_sound_size);
424 
425  Real rho_k_base = rho_k;
426  if (l_isentropic) {
427  // `rho_inp_sound` previously calculated in calc_rho_p_isentropic()
428  // is in HSE and, when multiplied by the specified
429  // `theta_input_sound`, give p, T, and theta that are consistent
430  // with each other.
431  //
432  // Here, we do not require thermodynamic consistency between the
433  // initial base state rho and the prognostic variables. Instead,
434  // we calculate a `rho_hse` that is consistent with the isentropic
435  // (constant theta) assumption.
436  rho_k_base = rhoTh_k / theta_inp_sound[0];
437  }
438 
439  // Set the density
440  state(i, j, k, Rho_comp) = rho_k;
441 
442  // Initial Rho0*Theta0
443  state(i, j, k, RhoTheta_comp) = rhoTh_k;
444 
445  // Total nonprecipitating water (Q1) == water vapor (Qv), i.e., there
446  // is no cloud water or cloud ice
447  Real qv_k = zero;
448  if (l_moist) {
449  qv_k = interpolate_1d(z_inp_sound, qv_inp_sound, z, inp_sound_size);
450  state(i, j, k, RhoQ1_comp) = rho_k * qv_k;
451  }
452 
453  // Update hse quantities with values calculated from InputSoundingData.calc_rho_p()
454  if (anel_assume_dry) qv_k = 0;
455  r_hse_arr (i,j,k) = rho_k_base;
456  p_hse_arr (i,j,k) = getPgivenRTh(rhoTh_k, qv_k);
457  pi_hse_arr(i,j,k) = getExnergivenRTh(rhoTh_k, l_rdOcp, qv_k);
458  th_hse_arr(i,j,k) = getRhoThetagivenP(p_hse_arr(i,j,k), qv_k) / rho_k_base;
459  qv_hse_arr(i,j,k) = qv_k;
460 
461  if (l_isentropic) {
462 #if 0
463  if (i==0 && j==0) {
464  Print() << "HSE rho,p,T=pi*th,th at " << IntVect(i,j,k) << " : "
465  << r_hse_arr(i,j,k) << " "
466  << p_hse_arr(i,j,k) << " "
467  << pi_hse_arr(i,j,k)*th_hse_arr(i,j,k) << " "
468  << th_hse_arr(i,j,k)
469  << " with rho,rhotheta=" << rho_k << " " << rhoTh_k
470  << std::endl;
471  }
472 #endif
473  // If everything above is thermodynamically consistent, this should be constant
474  AMREX_ALWAYS_ASSERT(std::abs(th_hse_arr(i,j,k) - theta_inp_sound[0]) < 1e-12);
475  }
476 
477  // FOEXTRAP hse arrays
478  if (k==kbot)
479  {
480  for (int kk = 1; kk <= ngz; kk++) {
481  r_hse_arr(i, j, k-kk) = r_hse_arr(i,j,k);
482  p_hse_arr(i, j, k-kk) = p_hse_arr(i,j,k);
483  pi_hse_arr(i, j, k-kk) = pi_hse_arr(i,j,k);
484  th_hse_arr(i, j, k-kk) = th_hse_arr(i,j,k);
485  qv_hse_arr(i, j, k-kk) = qv_hse_arr(i,j,k);
486  }
487  }
488  else if (k==ktop)
489  {
490  for (int kk = 1; kk <= ngz; kk++) {
491  r_hse_arr(i, j, k+kk) = r_hse_arr(i,j,k);
492  p_hse_arr(i, j, k+kk) = p_hse_arr(i,j,k);
493  pi_hse_arr(i, j, k+kk) = pi_hse_arr(i,j,k);
494  th_hse_arr(i, j, k+kk) = th_hse_arr(i,j,k);
495  qv_hse_arr(i, j, k+kk) = qv_hse_arr(i,j,k);
496  }
497  }
498  });
499 }
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getRhoThetagivenP(const amrex::Real p, const amrex::Real qv=amrex::Real(0))
Definition: ERF_EOS.H:172
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getExnergivenRTh(const amrex::Real rhotheta, const amrex::Real rdOcp, const amrex::Real qv=amrex::Real(0))
Definition: ERF_EOS.H:156
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getPgivenRTh(const amrex::Real rhotheta, const amrex::Real qv=amrex::Real(0))
Definition: ERF_EOS.H:81
AMREX_ALWAYS_ASSERT(bx.length()[2]==khi+1)
constexpr amrex::Real zero
Definition: ERF_NumericalConstants.H:36

Referenced by ERF::init_from_input_sounding().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ init_velocities_from_input_sounding()

void init_velocities_from_input_sounding ( const Box &  bx,
Array4< Real > const &  x_vel,
Array4< Real > const &  y_vel,
Array4< Real > const &  z_vel,
GeometryData const &  geomdata,
Array4< const Real > const &  z_nd_arr,
InputSoundingData const &  inputSoundingData 
)

Box level wrapper for initializing velocities from input sounding data.

Parameters
[in]bxBox specifying the indices we are initializing
[out]x_velArray4 specifying the x-velocity data we are to initialize
[out]y_velArray4 specifying the y-velocity data we are to initialize
[out]z_velArray4 specifying the z-velocity data we are to initialize
[in]geomdataGeometryData object specifying the domain geometry
[in]z_nd_arrArray4 specifying node-centered terrain heights, if present
[in]inputSoundingDataInputSoundingData object we are to initialize from

Box level wrapper for initializing velocities from input sounding data.

Parameters
bxBox specifying the indices we are initializing
x_velArray4 specifying the x-velocity data we are to initialize
y_velArray4 specifying the y-velocity data we are to initialize
z_velArray4 specifying the z-velocity data we are to initialize
geomdataGeometryData object specifying the domain geometry
z_nd_arrArray4 specifying node-centered terrain heights, if present
inputSoundingDataInputSoundingData object we are to initialize from
520 {
521  const Real* z_inp_sound = inputSoundingData.z_inp_sound_d[0].dataPtr();
522  const Real* U_inp_sound = inputSoundingData.U_inp_sound_d[0].dataPtr();
523  const Real* V_inp_sound = inputSoundingData.V_inp_sound_d[0].dataPtr();
524  const int inp_sound_size = inputSoundingData.size(0);
525 
526  // Geometry
527  const Real* prob_lo = geomdata.ProbLo();
528  const Real* dx = geomdata.CellSize();
529  const Real z_lo = prob_lo[2];
530  const Real dz = dx[2];
531 
532  // We want to set the lateral BC values, too
533  Box gbx = bx; // Copy constructor
534  gbx.grow(0,1); gbx.grow(1,1); // Grow by one in the lateral directions
535 
536  // Construct a box that is on x-faces
537  const Box& xbx = surroundingNodes(gbx,0);
538  // Construct a box that is on y-faces
539  const Box& ybx = surroundingNodes(gbx,1);
540  // Construct a box that is on z-faces
541  const Box& zbx = surroundingNodes(gbx,2);
542 
543  // Set the x,y,z-velocities
544  ParallelFor(xbx, ybx, zbx,
545  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
546  // Note that this is called on a box of x-faces
547  const Real z = (z_nd_arr) ? fourth*( z_nd_arr(i,j ,k )
548  + z_nd_arr(i,j+1,k )
549  + z_nd_arr(i,j ,k+1)
550  + z_nd_arr(i,j+1,k+1))
551  : z_lo + (k + myhalf) * dz;
552 
553  // Set the x-velocity
554  x_vel(i, j, k) = interpolate_1d(z_inp_sound, U_inp_sound, z, inp_sound_size);
555  },
556  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
557  // Note that this is called on a box of y-faces
558  const Real z = (z_nd_arr) ? fourth*( z_nd_arr(i ,j,k )
559  + z_nd_arr(i+1,j,k )
560  + z_nd_arr(i ,j,k+1)
561  + z_nd_arr(i+1,j,k+1))
562  : z_lo + (k + myhalf) * dz;
563 
564  // Set the y-velocity
565  y_vel(i, j, k) = interpolate_1d(z_inp_sound, V_inp_sound, z, inp_sound_size);
566  },
567  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
568  // Note that this is called on a box of z-faces
569  // Set the z-velocity
570  z_vel(i, j, k) = zero;
571  });
572 }
constexpr amrex::Real fourth
Definition: ERF_NumericalConstants.H:42

Referenced by ERF::init_from_input_sounding().

Here is the call graph for this function:
Here is the caller graph for this function: