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

bx	Box specifying the indices we are initializing
state	Array4 specifying the state data we are to initialize
geomdata	GeometryData object specifying the domain geometry
inputSoundingData	InputSoundingData object we are to initialize from

{
    const Real* z_inp_sound     = inputSoundingData.z_inp_sound_d[0].dataPtr();
    const Real* theta_inp_sound = inputSoundingData.theta_inp_sound_d[0].dataPtr();
    const Real* qv_inp_sound    = inputSoundingData.qv_inp_sound_d[0].dataPtr();
    const int   inp_sound_size  = inputSoundingData.size(0);
 
    // Geometry
    const Real* prob_lo = geomdata.ProbLo();
    const Real* dx = geomdata.CellSize();
    const Real  z_lo = prob_lo[2];
    const Real  dz   = dx[2];
 
    // We want to set the lateral BC values, too
    Box gbx = bx; // Copy constructor
    gbx.grow(0,1); gbx.grow(1,1); // Grow by one in the lateral directions
 
    ParallelFor(gbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
        const Real z = (z_cc_arr) ? z_cc_arr(i,j,k) : z_lo + (k + myhalf) * dz;
 
        Real rho_0 = one;
 
        // Set the density
        state(i, j, k, Rho_comp) = rho_0;
 
        // Initial Rho0*Theta0
        state(i, j, k, RhoTheta_comp) = rho_0 * interpolate_1d(z_inp_sound, theta_inp_sound, z, inp_sound_size);
 
        // total nonprecipitating water (Q1) == water vapor (Qv), i.e., there is no cloud water or cloud ice
        if (l_moist) {
            state(i, j, k, RhoQ1_comp) = rho_0 * interpolate_1d(z_inp_sound, qv_inp_sound, z, inp_sound_size);
        }
    });
}

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 &	,
		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

bx	Box specifying the indices we are initializing
state	Array4 specifying the state data we are to initialize
r_hse_arr	Array4 specifying the density HSE base state data we are to initialize
p_hse_arr	Array4 specifying the pressure HSE base state data we are to initialize
pi_hse_arr	Array4 specifying the Exner pressure HSE base state data we are to initialize
th_hse_arr	Array4 specifying the base state potential temperature we are to initialize
geomdata	GeometryData object specifying the domain geometry
l_gravity	Real number specifying the gravitational acceleration constant
l_rdOcp	Real number specifying the Rhydberg constant ($R_d$) divided by specific heat at constant pressure ($c_p$)
inputSoundingData	InputSoundingData object we are to initialize from

{
    const Real* z_inp_sound     = inputSoundingData.z_inp_sound_d[0].dataPtr();
    const Real* rho_inp_sound   = inputSoundingData.rho_inp_sound_d.dataPtr();
    const Real* theta_inp_sound = inputSoundingData.theta_inp_sound_d[0].dataPtr();
    const Real* qv_inp_sound    = inputSoundingData.qv_inp_sound_d[0].dataPtr();
    const int   inp_sound_size  = inputSoundingData.size(0);
    const bool anel_assume_dry  = inputSoundingData.assume_dry;
 
    // Geometry
    const Real* prob_lo = geomdata.ProbLo();
    const Real* dx = geomdata.CellSize();
    const Real  z_lo = prob_lo[2];
    const Real  dz   = dx[2];
 
    int kbot = geomdata.Domain().smallEnd(2);
    int ktop = geomdata.Domain().bigEnd(2);
 
    // We want to set the lateral BC values, too
    Box gbx = bx; // Copy constructor
    gbx.grow(0,1); gbx.grow(1,1); // Grow by one in the lateral directions
 
    ParallelFor(gbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
        const Real z = (z_cc_arr) ? z_cc_arr(i,j,k)
                                  : z_lo + (k + myhalf) * dz;
 
        Real rho_k   = interpolate_1d(z_inp_sound, rho_inp_sound, z, inp_sound_size);
        Real rhoTh_k = rho_k * interpolate_1d(z_inp_sound, theta_inp_sound, z, inp_sound_size);
 
        Real rho_k_base = rho_k;
        if (l_isentropic) {
            // `rho_inp_sound` previously calculated in calc_rho_p_isentropic()
            // is in HSE and, when multiplied by the specified
            // `theta_input_sound`, give p, T, and theta that are consistent
            // with each other.
            //
            // Here, we do not require thermodynamic consistency between the
            // initial base state rho and the prognostic variables. Instead,
            // we calculate a `rho_hse` that is consistent with the isentropic
            // (constant theta) assumption.
            rho_k_base = rhoTh_k / theta_inp_sound[0];
        }
 
        // Set the density
        state(i, j, k, Rho_comp) = rho_k;
 
        // Initial Rho0*Theta0
        state(i, j, k, RhoTheta_comp) = rhoTh_k;
 
        // Total nonprecipitating water (Q1) == water vapor (Qv), i.e., there
        // is no cloud water or cloud ice
        Real qv_k = zero;
        if (l_moist) {
            qv_k = interpolate_1d(z_inp_sound, qv_inp_sound, z, inp_sound_size);
            state(i, j, k, RhoQ1_comp) = rho_k * qv_k;
        }
 
        // Update hse quantities with values calculated from InputSoundingData.calc_rho_p()
        if (anel_assume_dry) qv_k = 0;
        r_hse_arr (i,j,k) = rho_k_base;
        p_hse_arr (i,j,k) = getPgivenRTh(rhoTh_k, qv_k);
        pi_hse_arr(i,j,k) = getExnergivenRTh(rhoTh_k, l_rdOcp, qv_k);
        th_hse_arr(i,j,k) = getRhoThetagivenP(p_hse_arr(i,j,k), qv_k) / rho_k_base;
        qv_hse_arr(i,j,k) = qv_k;
 
        if (l_isentropic) {
#if 0
            if (i==0 && j==0) {
                Print() << "HSE rho,p,T=pi*th,th at " << IntVect(i,j,k) << " : "
                    << r_hse_arr(i,j,k) << " "
                    << p_hse_arr(i,j,k) << " "
                    << pi_hse_arr(i,j,k)*th_hse_arr(i,j,k) << " "
                    << th_hse_arr(i,j,k)
                    << " with rho,rhotheta=" << rho_k << " " << rhoTh_k
                    << std::endl;
            }
#endif
            // If everything above is thermodynamically consistent, this should be constant
            AMREX_ALWAYS_ASSERT(std::abs(th_hse_arr(i,j,k) - theta_inp_sound[0]) < 1e-12);
        }
 
        // FOEXTRAP hse arrays
        if (k==kbot)
        {
            for (int kk = 1; kk <= ngz; kk++) {
                 r_hse_arr(i, j, k-kk) =  r_hse_arr(i,j,k);
                 p_hse_arr(i, j, k-kk) =  p_hse_arr(i,j,k);
                pi_hse_arr(i, j, k-kk) = pi_hse_arr(i,j,k);
                th_hse_arr(i, j, k-kk) = th_hse_arr(i,j,k);
                qv_hse_arr(i, j, k-kk) = qv_hse_arr(i,j,k);
            }
        }
        else if (k==ktop)
        {
            for (int kk = 1; kk <= ngz; kk++) {
                 r_hse_arr(i, j, k+kk) =  r_hse_arr(i,j,k);
                 p_hse_arr(i, j, k+kk) =  p_hse_arr(i,j,k);
                pi_hse_arr(i, j, k+kk) = pi_hse_arr(i,j,k);
                th_hse_arr(i, j, k+kk) = th_hse_arr(i,j,k);
                qv_hse_arr(i, j, k+kk) = qv_hse_arr(i,j,k);
            }
        }
    });
}

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

bx	Box specifying the indices we are initializing
x_vel	Array4 specifying the x-velocity data we are to initialize
y_vel	Array4 specifying the y-velocity data we are to initialize
z_vel	Array4 specifying the z-velocity data we are to initialize
geomdata	GeometryData object specifying the domain geometry
inputSoundingData	InputSoundingData object we are to initialize from

{
    const Real* z_inp_sound     = inputSoundingData.z_inp_sound_d[0].dataPtr();
    const Real* U_inp_sound     = inputSoundingData.U_inp_sound_d[0].dataPtr();
    const Real* V_inp_sound     = inputSoundingData.V_inp_sound_d[0].dataPtr();
    const int   inp_sound_size  = inputSoundingData.size(0);
 
    // Geometry
    const Real* prob_lo = geomdata.ProbLo();
    const Real* dx = geomdata.CellSize();
    const Real  z_lo = prob_lo[2];
    const Real  dz   = dx[2];
 
    // We want to set the lateral BC values, too
    Box gbx = bx; // Copy constructor
    gbx.grow(0,1); gbx.grow(1,1); // Grow by one in the lateral directions
 
    // Construct a box that is on x-faces
    const Box& xbx = surroundingNodes(gbx,0);
    // Construct a box that is on y-faces
    const Box& ybx = surroundingNodes(gbx,1);
    // Construct a box that is on z-faces
    const Box& zbx = surroundingNodes(gbx,2);
 
    // Set the x,y,z-velocities
    ParallelFor(xbx, ybx, zbx,
    [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
        // Note that this is called on a box of x-faces
        const Real z = (z_nd_arr) ? fourth*( z_nd_arr(i,j  ,k  )
                                         + z_nd_arr(i,j+1,k  )
                                         + z_nd_arr(i,j  ,k+1)
                                         + z_nd_arr(i,j+1,k+1))
                                  : z_lo + (k + myhalf) * dz;
 
        // Set the x-velocity
        x_vel(i, j, k) = interpolate_1d(z_inp_sound, U_inp_sound, z, inp_sound_size);
    },
    [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
        // Note that this is called on a box of y-faces
        const Real z = (z_nd_arr) ? fourth*( z_nd_arr(i  ,j,k  )
                                         + z_nd_arr(i+1,j,k  )
                                         + z_nd_arr(i  ,j,k+1)
                                         + z_nd_arr(i+1,j,k+1))
                                  : z_lo + (k + myhalf) * dz;
 
        // Set the y-velocity
        y_vel(i, j, k) = interpolate_1d(z_inp_sound, V_inp_sound, z, inp_sound_size);
    },
    [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
        // Note that this is called on a box of z-faces
        // Set the z-velocity
        z_vel(i, j, k) = zero;
    });
}

Referenced by ERF::init_from_input_sounding().

Here is the call graph for this function:

Here is the caller graph for this function: