ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
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_bx_scalars_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_bx_scalars_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)
 
void init_bx_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_bx_scalars_from_input_sounding()

void init_bx_scalars_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
bxBox specifying the indices we are initializing
stateArray4 specifying the state data we are to initialize
geomdataGeometryData object specifying the domain geometry
inputSoundingDataInputSoundingData object we are to initialize from
182 {
183  const Real* z_inp_sound = inputSoundingData.z_inp_sound_d[0].dataPtr();
184  const Real* theta_inp_sound = inputSoundingData.theta_inp_sound_d[0].dataPtr();
185  const Real* qv_inp_sound = inputSoundingData.qv_inp_sound_d[0].dataPtr();
186  const int inp_sound_size = inputSoundingData.size(0);
187 
188  // Geometry
189  const Real* prob_lo = geomdata.ProbLo();
190  const Real* dx = geomdata.CellSize();
191  const Real z_lo = prob_lo[2];
192  const Real dz = dx[2];
193 
194  // We want to set the lateral BC values, too
195  Box gbx = bx; // Copy constructor
196  gbx.grow(0,1); gbx.grow(1,1); // Grow by one in the lateral directions
197 
198  ParallelFor(gbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
199  const Real z = (z_cc_arr) ? z_cc_arr(i,j,k)
200  : z_lo + (k + 0.5) * dz;
201 
202  Real rho_0 = 1.0;
203 
204  // Set the density
205  state(i, j, k, Rho_comp) = rho_0;
206 
207  // Initial Rho0*Theta0
208  state(i, j, k, RhoTheta_comp) = rho_0 * interpolate_1d(z_inp_sound, theta_inp_sound, z, inp_sound_size);
209 
210  // Initialize all scalars to 0.
211  for (int n = 0; n < NSCALARS; n++) {
212  state(i, j, k, RhoScalar_comp+n) = 0;
213  }
214 
215  // total nonprecipitating water (Q1) == water vapor (Qv), i.e., there is no cloud water or cloud ice
216  if (l_moist) {
217  state(i, j, k, RhoQ1_comp) = rho_0 * interpolate_1d(z_inp_sound, qv_inp_sound, z, inp_sound_size);
218  }
219  });
220 }
#define RhoScalar_comp
Definition: ERF_IndexDefines.H:40
#define Rho_comp
Definition: ERF_IndexDefines.H:36
#define RhoTheta_comp
Definition: ERF_IndexDefines.H:37
#define NSCALARS
Definition: ERF_IndexDefines.H:16
#define RhoQ1_comp
Definition: ERF_IndexDefines.H:42
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:12

Referenced by ERF::init_from_input_sounding().

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

◆ init_bx_scalars_from_input_sounding_hse()

void init_bx_scalars_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 &  ,
const Real &  l_rdOcp,
const bool &  l_moist,
InputSoundingData const &  inputSoundingData,
const bool &  l_isentropic 
)

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
geomdataGeometryData object specifying the domain geometry
l_gravityReal number specifying the gravitational acceleration constant
l_rdOcpReal number specifying the Rhydberg constant ($R_d$) divided by specific heat at constant pressure ($c_p$)
inputSoundingDataInputSoundingData object we are to initialize from
252 {
253  const Real* z_inp_sound = inputSoundingData.z_inp_sound_d[0].dataPtr();
254  const Real* rho_inp_sound = inputSoundingData.rho_inp_sound_d.dataPtr();
255  const Real* theta_inp_sound = inputSoundingData.theta_inp_sound_d[0].dataPtr();
256  const Real* qv_inp_sound = inputSoundingData.qv_inp_sound_d[0].dataPtr();
257  const int inp_sound_size = inputSoundingData.size(0);
258  const bool anel_assume_dry = inputSoundingData.assume_dry;
259 
260  // Geometry
261  const Real* prob_lo = geomdata.ProbLo();
262  const Real* dx = geomdata.CellSize();
263  const Real z_lo = prob_lo[2];
264  const Real dz = dx[2];
265 
266  int kbot = geomdata.Domain().smallEnd(2);
267  int ktop = geomdata.Domain().bigEnd(2);
268 
269  // We want to set the lateral BC values, too
270  Box gbx = bx; // Copy constructor
271  gbx.grow(0,1); gbx.grow(1,1); // Grow by one in the lateral directions
272 
273  ParallelFor(gbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
274  const Real z = (z_cc_arr) ? z_cc_arr(i,j,k)
275  : z_lo + (k + 0.5) * dz;
276 
277  Real rho_k = interpolate_1d(z_inp_sound, rho_inp_sound, z, inp_sound_size);
278  Real rhoTh_k = rho_k * interpolate_1d(z_inp_sound, theta_inp_sound, z, inp_sound_size);
279 
280  Real rho_k_base = rho_k;
281  if (l_isentropic) {
282  // `rho_inp_sound` previously calculated in calc_rho_p_isentropic()
283  // is in HSE and, when multiplied by the specified
284  // `theta_input_sound`, give p, T, and theta that are consistent
285  // with each other.
286  //
287  // Here, we do not require thermodynamic consistency between the
288  // initial base state rho and the prognostic variables. Instead,
289  // we calculate a `rho_hse` that is consistent with the isentropic
290  // (constant theta) assumption.
291  rho_k_base = rhoTh_k / theta_inp_sound[0];
292  }
293 
294  // Set the density
295  state(i, j, k, Rho_comp) = rho_k;
296 
297  // Initial Rho0*Theta0
298  state(i, j, k, RhoTheta_comp) = rhoTh_k;
299 
300  // Initialize all scalars to 0.
301  for (int n = 0; n < NSCALARS; n++) {
302  state(i, j, k, RhoScalar_comp+n) = 0;
303  }
304 
305  // total nonprecipitating water (Q1) == water vapor (Qv), i.e., there
306  // is no cloud water or cloud ice
307  Real qv_k = 0.0;
308  if (l_moist) {
309  qv_k = interpolate_1d(z_inp_sound, qv_inp_sound, z, inp_sound_size);
310  state(i, j, k, RhoQ1_comp) = rho_k * qv_k;
311  }
312 
313  // Update hse quantities with values calculated from InputSoundingData.calc_rho_p()
314  if (anel_assume_dry) qv_k = 0;
315  r_hse_arr (i,j,k) = rho_k_base * (1.0 + qv_k);
316  p_hse_arr (i,j,k) = getPgivenRTh(rhoTh_k, qv_k);
317  pi_hse_arr(i,j,k) = getExnergivenRTh(rhoTh_k, l_rdOcp);
318  th_hse_arr(i,j,k) = getRhoThetagivenP(p_hse_arr(i,j,k), qv_k) / rho_k_base;
319  qv_hse_arr(i,j,k) = qv_k;
320  if (l_isentropic) {
321 #if 0
322  if (i==0 && j==0) {
323  Print() << "HSE rho,p,T=pi*th,th at " << IntVect(i,j,k) << " : "
324  << r_hse_arr(i,j,k) << " "
325  << p_hse_arr(i,j,k) << " "
326  << pi_hse_arr(i,j,k)*th_hse_arr(i,j,k) << " "
327  << th_hse_arr(i,j,k)
328  << " with rho,rhotheta=" << rho_k << " " << rhoTh_k
329  << std::endl;
330  }
331 #endif
332  // If everything above is thermodynamically consistent, this should be constant
333  AMREX_ALWAYS_ASSERT(std::abs(th_hse_arr(i,j,k) - theta_inp_sound[0]) < 1e-12);
334  }
335 
336  // TODO: we should be setting this to the number of ghost cells of base_state[lev]
337  // instead of hard-wiring it here!
338  int ng = 3;
339 
340  // FOEXTRAP hse arrays
341  if (k==kbot)
342  {
343  for (int kk = 1; kk <= ng; kk++) {
344  r_hse_arr(i, j, k-kk) = r_hse_arr(i,j,k);
345  p_hse_arr(i, j, k-kk) = p_hse_arr(i,j,k);
346  pi_hse_arr(i, j, k-kk) = pi_hse_arr(i,j,k);
347  th_hse_arr(i, j, k-kk) = th_hse_arr(i,j,k);
348  qv_hse_arr(i, j, k-kk) = qv_hse_arr(i,j,k);
349  }
350  }
351  else if (k==ktop)
352  {
353  for (int kk = 1; kk <= ng; kk++) {
354  r_hse_arr(i, j, k+kk) = r_hse_arr(i,j,k);
355  p_hse_arr(i, j, k+kk) = p_hse_arr(i,j,k);
356  pi_hse_arr(i, j, k+kk) = pi_hse_arr(i,j,k);
357  qv_hse_arr(i, j, k+kk) = qv_hse_arr(i,j,k);
358  }
359  }
360  });
361 }
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getPgivenRTh(const amrex::Real rhotheta, const amrex::Real qv=0.)
Definition: ERF_EOS.H:81
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getExnergivenRTh(const amrex::Real rhotheta, const amrex::Real rdOcp, const amrex::Real qv=0.0)
Definition: ERF_EOS.H:156
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getRhoThetagivenP(const amrex::Real p, const amrex::Real qv=0.0)
Definition: ERF_EOS.H:172
@ ng
Definition: ERF_Morrison.H:48

Referenced by ERF::init_from_input_sounding().

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

◆ init_bx_velocities_from_input_sounding()

void init_bx_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
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
inputSoundingDataInputSoundingData object we are to initialize from
381 {
382  const Real* z_inp_sound = inputSoundingData.z_inp_sound_d[0].dataPtr();
383  const Real* U_inp_sound = inputSoundingData.U_inp_sound_d[0].dataPtr();
384  const Real* V_inp_sound = inputSoundingData.V_inp_sound_d[0].dataPtr();
385  const int inp_sound_size = inputSoundingData.size(0);
386 
387  // Geometry
388  const Real* prob_lo = geomdata.ProbLo();
389  const Real* dx = geomdata.CellSize();
390  const Real z_lo = prob_lo[2];
391  const Real dz = dx[2];
392 
393  // We want to set the lateral BC values, too
394  Box gbx = bx; // Copy constructor
395  gbx.grow(0,1); gbx.grow(1,1); // Grow by one in the lateral directions
396 
397  // Construct a box that is on x-faces
398  const Box& xbx = surroundingNodes(gbx,0);
399  // Construct a box that is on y-faces
400  const Box& ybx = surroundingNodes(gbx,1);
401  // Construct a box that is on z-faces
402  const Box& zbx = surroundingNodes(gbx,2);
403 
404  // Set the x,y,z-velocities
405  ParallelFor(xbx, ybx, zbx,
406  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
407  // Note that this is called on a box of x-faces
408  const Real z = (z_nd_arr) ? 0.25*( z_nd_arr(i,j ,k )
409  + z_nd_arr(i,j+1,k )
410  + z_nd_arr(i,j ,k+1)
411  + z_nd_arr(i,j+1,k+1))
412  : z_lo + (k + 0.5) * dz;
413 
414  // Set the x-velocity
415  x_vel(i, j, k) = interpolate_1d(z_inp_sound, U_inp_sound, z, inp_sound_size);
416  },
417  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
418  // Note that this is called on a box of y-faces
419  const Real z = (z_nd_arr) ? 0.25*( z_nd_arr(i ,j,k )
420  + z_nd_arr(i+1,j,k )
421  + z_nd_arr(i ,j,k+1)
422  + z_nd_arr(i+1,j,k+1))
423  : z_lo + (k + 0.5) * dz;
424 
425  // Set the y-velocity
426  y_vel(i, j, k) = interpolate_1d(z_inp_sound, V_inp_sound, z, inp_sound_size);
427  },
428  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
429  // Note that this is called on a box of z-faces
430  // Set the z-velocity
431  z_vel(i, j, k) = 0.0;
432  });
433 }
const Box zbx
Definition: ERF_DiffSetup.H:23
const Box xbx
Definition: ERF_DiffSetup.H:21
const Box ybx
Definition: ERF_DiffSetup.H:22

Referenced by ERF::init_from_input_sounding().

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