ProblemBase Class Reference

#include <ERF_ProbCommon.H>

Collaboration diagram for ProblemBase:

Public Member Functions
virtual	~ProblemBase ()=default
virtual void	erf_init_dens_hse (amrex::MultiFab &, std::unique_ptr< amrex::MultiFab > &, std::unique_ptr< amrex::MultiFab > &, amrex::Geometry const &)
virtual void	erf_init_dens_hse_moist (amrex::MultiFab &, std::unique_ptr< amrex::MultiFab > &, amrex::Geometry const &)
virtual void	init_custom_pert (const amrex::Box &, const amrex::Box &, const amrex::Box &, const amrex::Box &, amrex::Array4< amrex::Real const > const &, amrex::Array4< amrex::Real > const &, amrex::Array4< amrex::Real > const &, amrex::Array4< amrex::Real > const &, amrex::Array4< amrex::Real > const &, amrex::Array4< amrex::Real > const &, amrex::Array4< amrex::Real > const &, amrex::Array4< amrex::Real const > const &, amrex::Array4< amrex::Real const > const &, amrex::GeometryData const &, amrex::Array4< amrex::Real const > const &, amrex::Array4< amrex::Real const > const &, amrex::Array4< amrex::Real const > const &, const SolverChoice &)
virtual void	update_rhotheta_sources (const amrex::Real &, amrex::Vector< amrex::Real > &src, amrex::Gpu::DeviceVector< amrex::Real > &d_src, const amrex::Geometry &geom, std::unique_ptr< amrex::MultiFab > &z_phys_cc)
virtual void	update_rhoqt_sources (const amrex::Real &, amrex::Vector< amrex::Real > &qsrc, amrex::Gpu::DeviceVector< amrex::Real > &d_qsrc, const amrex::Geometry &geom, std::unique_ptr< amrex::MultiFab > &z_phys_cc)
virtual void	update_w_subsidence (const amrex::Real &, amrex::Vector< amrex::Real > &wbar, amrex::Gpu::DeviceVector< amrex::Real > &d_wbar, const amrex::Geometry &geom, std::unique_ptr< amrex::MultiFab > &z_phys_cc)
virtual void	update_geostrophic_profile (const amrex::Real &, amrex::Vector< amrex::Real > &u_geos, amrex::Gpu::DeviceVector< amrex::Real > &d_u_geos, amrex::Vector< amrex::Real > &v_geos, amrex::Gpu::DeviceVector< amrex::Real > &d_v_geos, const amrex::Geometry &geom, std::unique_ptr< amrex::MultiFab > &z_phys_cc)
void	init_terrain_surface (const amrex::Geometry &geom, amrex::FArrayBox &terrain_fab, const amrex::Real &time)
void	read_custom_terrain (const std::string &fname, const bool is_usgs, const amrex::Geometry &geom, amrex::FArrayBox &terrain_fab, const amrex::Real &)
virtual void	init_custom_terrain (const amrex::Geometry &, amrex::FArrayBox &terrain_fab, const amrex::Real &)
virtual void	erf_init_rayleigh (amrex::Vector< amrex::Vector< amrex::Real > > &, amrex::Geometry const &, std::unique_ptr< amrex::MultiFab > &, amrex::Real)
void	init_uniform (const amrex::Box &bx, amrex::Array4< amrex::Real > const &state)

Protected Member Functions
void	init_base_parms (amrex::Real rho_0, amrex::Real T_0)
virtual std::string	name ()=0

Protected Attributes
ProbParmDefaults	base_parms

Detailed Description

Class to hold problem-specific routines

Constructor & Destructor Documentation

~ProblemBase()

virtual ProblemBase::~ProblemBase ( )

virtualdefault

Virtual destructor to avoid data leakage with derived class

Member Function Documentation

erf_init_dens_hse()

virtual void ProblemBase::erf_init_dens_hse ( amrex::MultiFab &  , std::unique_ptr< amrex::MultiFab > &  , std::unique_ptr< amrex::MultiFab > &  , amrex::Geometry const &    )

inlinevirtual

Function to initialize the hydrostatic reference density

Parameters

[out]	rho_hse	hydrostatic reference density
[in]	z_phys_nd	height coordinate at nodes
[in]	z_phys_cc	height coordinate at cell centers
[in]	geom	container for geometric information

    {
        amrex::Print() << "Hydrostatically balanced density was NOT set"
            << " -- an appropriate init_type should probably have been specified"
            << " (e.g., input_sounding, ideal, real, or metgrid)"
            << std::endl;
        amrex::Error("Should never call erf_init_dens_hse for "+name()+" problem");
    }

Here is the call graph for this function:

erf_init_dens_hse_moist()

virtual void ProblemBase::erf_init_dens_hse_moist ( amrex::MultiFab &  , std::unique_ptr< amrex::MultiFab > &  , amrex::Geometry const &    )

inlinevirtual

    {
 
    }

erf_init_rayleigh()

virtual void ProblemBase::erf_init_rayleigh ( amrex::Vector< amrex::Vector< amrex::Real > > &  , amrex::Geometry const &  , std::unique_ptr< amrex::MultiFab > &  , amrex::Real    )

inlinevirtual

Function to define the quantities needed to impose Rayleigh damping

Parameters

[out]	rayleigh_ptrs	= {strength of Rayleigh damping, reference values for xvel/yvel/zvel/theta used to define Rayleigh damping}
[in]	geom	container for geometric information
[in]	z_phys_cc	height coordinate at cell centers

    {
        amrex::Error("Should never call erf_init_rayleigh for "+name()+" problem");
    }

Here is the call graph for this function:

init_base_parms()

void ProblemBase::init_base_parms ( amrex::Real  rho_0, amrex::Real  T_0  )

inlineprotected

Function to update default base parameters, currently only used for init_type == InitType::Uniform

                                                          {
        base_parms.rho_0 = rho_0;
        base_parms.T_0 = T_0;
    }

init_custom_pert()

virtual void ProblemBase::init_custom_pert ( const amrex::Box &  , const amrex::Box &  , const amrex::Box &  , const amrex::Box &  , amrex::Array4< amrex::Real const > const &  , amrex::Array4< amrex::Real > const &  , amrex::Array4< amrex::Real > const &  , amrex::Array4< amrex::Real > const &  , amrex::Array4< amrex::Real > const &  , amrex::Array4< amrex::Real > const &  , amrex::Array4< amrex::Real > const &  , amrex::Array4< amrex::Real const > const &  , amrex::Array4< amrex::Real const > const &  , amrex::GeometryData const &  , amrex::Array4< amrex::Real const > const &  , amrex::Array4< amrex::Real const > const &  , amrex::Array4< amrex::Real const > const &  , const SolverChoice &    )

inlinevirtual

Function to perform custom initialization of a test problem

Parameters

[in]	bx	cell-centered box on which to initialize scalars
[in]	xbx	box on which to initialize x-component of velocity
[in]	ybx	box on which to initialize y-component of velocity
[in]	zbx	box on which to initialize z-component of velocity
[out]	state	cell-centered variables to be filled in this routine
[out]	x_velocity	x-component of velocity to be filled in this routine
[out]	y_velocity	y-component of velocity to be filled in this routine
[out]	z_velocity	z-component of velocity to be filled in this routine
[out]	r_hse	hydrostatic reference density
[out]	p_hse	hydrostatic reference pressure
[in]	z_nd	height coordinate at nodes
[in]	z_cc	height coordinate at cell centers
[in]	qv	water vapor
[in]	qc	cloud water
[in]	qi	cloud ice
[in]	mf_m	map factor on cell centers
[in]	mf_u	map factor on x-faces
[in]	mf_v	map factor on y-faces
[in]	sc	SolverChoice structure that carries parameters

    {
        amrex::Print() << "No perturbation to background fields supplied for "
            << name() << " problem" << std::endl;
    }

Here is the call graph for this function:

init_custom_terrain()

virtual void ProblemBase::init_custom_terrain ( const amrex::Geometry &  , amrex::FArrayBox &  terrain_fab, const amrex::Real &    )

inlinevirtual

Function to perform custom initialization of terrain

Note: Terrain functionality can also be used to provide grid stretching.

Parameters

[in]	geom	container for geometric information
[out]	terrain_mf	height coordinate at nodes
[in]	time	current time

    {
        amrex::Print() << "Initializing flat terrain" << std::endl;
 
        terrain_fab.template setVal<amrex::RunOn::Device>(0.0);
 
        SolverChoice::set_flat_terrain_flag();
    }

Referenced by init_terrain_surface().

Here is the call graph for this function:

Here is the caller graph for this function:

init_terrain_surface()

void ProblemBase::init_terrain_surface ( const amrex::Geometry &  geom, amrex::FArrayBox &  terrain_fab, const amrex::Real &  time  )

inline

Function to perform custom initialization of terrain

Parameters

[in]	geom	container for geometric information
[out]	z_phys_nd	height coordinate at nodes
[in]	time	current time

    {
        // Check if a valid text file exists for the terrain
        std::string fname, fname_usgs;
        amrex::ParmParse pp("erf");
        auto valid_fname      = pp.query("terrain_file_name",fname);
        auto valid_fname_USGS = pp.query("terrain_file_name_USGS",fname_usgs);
 
        bool is_usgs;
 
        if (valid_fname) {
            is_usgs = false;
            read_custom_terrain(fname,is_usgs,geom,terrain_fab,time);
 
        } else if (valid_fname_USGS) {
            is_usgs = true;
            read_custom_terrain(fname_usgs,is_usgs,geom,terrain_fab,time);
 
        } else {
            init_custom_terrain (geom, terrain_fab, time);
        }
    }

Here is the call graph for this function:

init_uniform()

void ProblemBase::init_uniform ( const amrex::Box &  bx, amrex::Array4< amrex::Real > const &  state  )

inline

Function to set uniform background density and pressure fields

    {
        amrex::Real rho_0 = base_parms.rho_0;
        amrex::Real T_0 = base_parms.T_0;
        amrex::Print() << "Initializing uniform fields"
            << " rho=" << rho_0 << " theta=" << T_0
            << std::endl;
 
        ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
        {
            state(i, j, k, Rho_comp) = rho_0;
            state(i, j, k, RhoTheta_comp) = rho_0 * T_0;
        });
    }

name()

virtual std::string ProblemBase::name ( )

protectedpure virtual

Referenced by erf_init_dens_hse(), erf_init_rayleigh(), init_custom_pert(), update_geostrophic_profile(), update_rhoqt_sources(), update_rhotheta_sources(), and update_w_subsidence().

Here is the caller graph for this function:

read_custom_terrain()

void ProblemBase::read_custom_terrain ( const std::string &  fname, const bool  is_usgs, const amrex::Geometry &  geom, amrex::FArrayBox &  terrain_fab, const amrex::Real &    )

inline

    {
        amrex::Vector<amrex::Real> m_xterrain,m_yterrain,m_zterrain;
 
        int nx = 0; int ny = 0;
 
        if (amrex::ParallelDescriptor::IOProcessor()) {
            // Read terrain file
            amrex::Print()<<"Reading terrain file: "<< fname<< std::endl;
            std::ifstream file(fname);
 
            if (!file.is_open()) {
                amrex::Abort("Error: Could not open the file " + fname+ "\n");
            }
 
            // Check if file is empty
            if (file.peek() == std::ifstream::traits_type::eof()) {
                amrex::Abort("Error: The file  " + fname+ " is empty.\n");
            }
 
            amrex::Real value1,value2,value3;
 
            if (is_usgs) {
                amrex::Real lat_min, lon_min;
 
                file >> lon_min >> lat_min;
                 if(std::fabs(lon_min) > 180.0) {
                    amrex::Error("The value of longitude for entry in the first line in " + fname
                                 + " should not exceed 180. It is " + std::to_string(lon_min));
                }
                 if(std::fabs(lat_min) > 90.0) {
                    amrex::Error("The value of latitude for entry in the first line in " + fname
                                 + " should not exceed 90. It is " + std::to_string(lat_min));
                }
 
                file >> nx >> ny;
 
                int counter = 0;
                while (file >> value1 >> value2 >> value3) {
                    m_xterrain.push_back(value1);
                    if(counter%nx==0) {
                        m_yterrain.push_back(value2);
                    }
                    m_zterrain.push_back(value3);
                    counter += 1;
                }
                AMREX_ASSERT(m_xterrain.size() == static_cast<long int>(nx*ny));
                AMREX_ASSERT(m_yterrain.size() == static_cast<long int>(ny));
                AMREX_ASSERT(m_zterrain.size() == static_cast<long int>(nx*ny));
 
        } else {
            file >> nx; file >> ny;
            AMREX_ALWAYS_ASSERT(nx > 0);
            AMREX_ALWAYS_ASSERT(ny > 0);
            m_xterrain.resize(nx);
            m_yterrain.resize(ny);
            m_zterrain.resize(nx * ny);
            for (int n = 0; n < nx; n++) {
                file >> m_xterrain[n];
            }
            for (int n = 0; n < ny; n++) {
                file >> m_yterrain[n];
            }
            for (int n = 0; n < nx * ny; n++) {
                file >> m_zterrain[n];
            }
        }
 
        amrex::Print() << "NX " << nx << " " << ny << std::endl;
 
            // Close the file after reading
            file.close();
        }
 
        amrex::ParallelDescriptor::Bcast(&nx,1,amrex::ParallelDescriptor::IOProcessorNumber());
        amrex::ParallelDescriptor::Bcast(&ny,1,amrex::ParallelDescriptor::IOProcessorNumber());
 
        int nz = nx * ny;
 
        m_xterrain.resize(nx);
        m_yterrain.resize(ny);
        m_zterrain.resize(nz);
 
        amrex::ParallelDescriptor::Bcast(m_xterrain.data(),nx,amrex::ParallelDescriptor::IOProcessorNumber());
        amrex::ParallelDescriptor::Bcast(m_yterrain.data(),ny,amrex::ParallelDescriptor::IOProcessorNumber());
        amrex::ParallelDescriptor::Bcast(m_zterrain.data(),nz,amrex::ParallelDescriptor::IOProcessorNumber());
 
        // Copy data to the GPU
        amrex::Gpu::DeviceVector<amrex::Real> d_xterrain(nx),d_yterrain(ny),d_zterrain(nz);
        amrex::Gpu::copy(amrex::Gpu::hostToDevice, m_xterrain.begin(), m_xterrain.end(), d_xterrain.begin());
        amrex::Gpu::copy(amrex::Gpu::hostToDevice, m_yterrain.begin(), m_yterrain.end(), d_yterrain.begin());
        amrex::Gpu::copy(amrex::Gpu::hostToDevice, m_zterrain.begin(), m_zterrain.end(), d_zterrain.begin());
 
        amrex::Real* d_xt = d_xterrain.data();
        amrex::Real* d_yt = d_yterrain.data();
        amrex::Real* d_zt = d_zterrain.data();
 
        auto dx = geom.CellSizeArray();
        auto ProbLoArr = geom.ProbLoArray();
 
        int ilo = geom.Domain().smallEnd(0);
        int jlo = geom.Domain().smallEnd(1);
        int klo = geom.Domain().smallEnd(2);
        int ihi = geom.Domain().bigEnd(0) + 1;
        int jhi = geom.Domain().bigEnd(1) + 1;
 
        amrex::Box zbx = terrain_fab.box();
        amrex::Array4<amrex::Real> const& z_arr = terrain_fab.array();
 
        amrex::ParallelFor(zbx, [=] AMREX_GPU_DEVICE (int i, int j, int /*k*/)
        {
            // Clip indices for ghost cells
            int ii = amrex::min(amrex::max(i,ilo),ihi);
            int jj = amrex::min(amrex::max(j,jlo),jhi);
 
            // Location of nodes
            amrex::Real x = ProbLoArr[0] + ii * dx[0] + 1e-3;
            amrex::Real y = ProbLoArr[1] + jj * dx[1] + 1e-3;
 
            int jindex_terrain = -1;
            for (int it = 0; it < ny && jindex_terrain == -1; it++) {
                if (d_yt[it] > y) {
                    jindex_terrain = it-1;
                }
            }
            if (jindex_terrain == -1) {
                jindex_terrain = ny-1;
            }
 
            int ind11, ind12, ind21, ind22;
            amrex::Real x1, x2, y1, y2;
 
            int iindex_terrain=-1;
            if (is_usgs) {
                int gstart = (jindex_terrain  )*nx;
                int gend   = (jindex_terrain+1)*nx-1;
                for (int it = gstart; it <= gend && iindex_terrain == -1; it++) {
                    if (d_xt[it] > x) {
                        iindex_terrain = it-gstart;
                    }
                }
 
                // Define the four values to interpolate between
                ind11 = jindex_terrain*nx + iindex_terrain;
                ind12 = ind11+nx;
                ind21 = ind11+1;
                ind22 = ind12+1;
                x1 = d_xt[ind11];
                x2 = d_xt[ind21];
                y1 = d_yt[jindex_terrain];
                y2 = d_yt[jindex_terrain+1];
 
            } else {
                for (int it = 0; it < nx && iindex_terrain == -1; it++) {
                    if (d_xt[it] > x) {
                        iindex_terrain = it-1;
                    }
                }
                if (iindex_terrain == -1) {
                    iindex_terrain = nx-1;
                }
 
                // Define the four values to interpolate between
                x1 = d_xt[iindex_terrain];
                x2 = d_xt[iindex_terrain+1];
                y1 = d_yt[jindex_terrain];
                y2 = d_yt[jindex_terrain+1];
 
                ind11 =           jindex_terrain         * nx + iindex_terrain;
                ind21 =  std::min(jindex_terrain+1,ny-1) * nx + iindex_terrain;
 
                ind12 =           jindex_terrain         * nx + std::min(iindex_terrain+1,nx-1);
                ind22  = std::min(jindex_terrain+1,ny-1) * nx + std::min(iindex_terrain+1,nx-1);
            }
 
            // Do bilinear interpolation
            amrex::Real denom = (x2-x1)*(y2-y1);
            amrex::Real  w_11 = (x2-x)*(y2-y)/denom;
            amrex::Real  w_12 = (x2-x)*(y-y1)/denom;
            amrex::Real  w_21 = (x-x1)*(y2-y)/denom;
            amrex::Real  w_22 = (x-x1)*(y-y1)/denom;
 
            z_arr(i,j,klo) = w_11*d_zt[ind11] + w_12*d_zt[ind12] + w_21*d_zt[ind21] + w_22*d_zt[ind22];
        });
    }

Referenced by init_terrain_surface().

Here is the caller graph for this function:

update_geostrophic_profile()

virtual void ProblemBase::update_geostrophic_profile ( const amrex::Real &  , amrex::Vector< amrex::Real > &  u_geos, amrex::Gpu::DeviceVector< amrex::Real > &  d_u_geos, amrex::Vector< amrex::Real > &  v_geos, amrex::Gpu::DeviceVector< amrex::Real > &  d_v_geos, const amrex::Geometry &  geom, std::unique_ptr< amrex::MultiFab > &  z_phys_cc  )

inlinevirtual

Function to update user-specified geostrophic wind profile.

Parameters

[in]	time	current time
[out]	u_geos	geostrophic wind profile
[out]	v_geos	geostrophic wind profile
[in]	geom	container for geometric information
[in]	z_phys_cc	height coordinate at cell centers

    {
        if (u_geos.empty()) return;
 
        if (z_phys_cc) {
            // use_terrain=1
            amrex::Error("Geostrophic wind profile not defined for "+name()+" problem with terrain");
        } else {
            const int khi              = geom.Domain().bigEnd()[2];
            // const amrex::Real* prob_lo = geom.ProbLo();
            // const auto dx              = geom.CellSize();
            for (int k = 0; k <= khi; k++)
            {
                // const amrex::Real z_cc = prob_lo[2] + (k+0.5)* dx[2];
                // set RHS term of RhoTheta equation based on time, z_cc here
                u_geos[k] = 0.0;
                v_geos[k] = 0.0;
            }
        }
 
        // Copy from host version to device version
        amrex::Gpu::copy(amrex::Gpu::hostToDevice, u_geos.begin(), u_geos.end(), d_u_geos.begin());
        amrex::Gpu::copy(amrex::Gpu::hostToDevice, v_geos.begin(), v_geos.end(), d_v_geos.begin());
    }

Here is the call graph for this function:

update_rhoqt_sources()

virtual void ProblemBase::update_rhoqt_sources ( const amrex::Real &  , amrex::Vector< amrex::Real > &  qsrc, amrex::Gpu::DeviceVector< amrex::Real > &  d_qsrc, const amrex::Geometry &  geom, std::unique_ptr< amrex::MultiFab > &  z_phys_cc  )

inlinevirtual

Function to update user-specified moisture source terms that can vary with time and height.

Parameters

[in]	time	current time
[out]	rhoqt_source	moisture forcing profile
[in]	geom	container for geometric information
[in]	z_phys_cc	height coordinate at cell centers

    {
        if (qsrc.empty()) return;
 
        if (z_phys_cc) {
            // use_terrain=1
            amrex::Error("Moisture forcing not defined for "+name()+" problem with terrain");
        } else {
            const int khi              = geom.Domain().bigEnd()[2];
            // const amrex::Real* prob_lo = geom.ProbLo();
            // const auto dx              = geom.CellSize();
            for (int k = 0; k <= khi; k++)
            {
                // const amrex::Real z_cc = prob_lo[2] + (k+0.5)* dx[2];
                // set RHS term of RhoQ1 equation based on time, z_cc here
                qsrc[k] = 0.0;
            }
        }
 
        // Copy from host version to device version
        amrex::Gpu::copy(amrex::Gpu::hostToDevice, qsrc.begin(), qsrc.end(), d_qsrc.begin());
    }

Here is the call graph for this function:

update_rhotheta_sources()

virtual void ProblemBase::update_rhotheta_sources ( const amrex::Real &  , amrex::Vector< amrex::Real > &  src, amrex::Gpu::DeviceVector< amrex::Real > &  d_src, const amrex::Geometry &  geom, std::unique_ptr< amrex::MultiFab > &  z_phys_cc  )

inlinevirtual

Function to update user-specified temperature source terms that can vary with time and height.

Parameters

[in]	time	current time
[out]	rhotheta_source	forcing profile
[in]	geom	container for geometric information
[in]	z_phys_cc	height coordinate at cell centers

    {
        if (src.empty()) return;
 
        if (z_phys_cc) {
            // use_terrain=1
            amrex::Error("Temperature forcing not defined for "+name()+" problem with terrain");
        } else {
            const int khi              = geom.Domain().bigEnd()[2];
            // const amrex::Real* prob_lo = geom.ProbLo();
            // const auto dx              = geom.CellSize();
            for (int k = 0; k <= khi; k++)
            {
                // const amrex::Real z_cc = prob_lo[2] + (k+0.5)* dx[2];
                // set RHS term of RhoTheta equation based on time, z_cc here
                src[k] = 0.0;
            }
        }
 
        // Copy from host version to device version
        amrex::Gpu::copy(amrex::Gpu::hostToDevice, src.begin(), src.end(), d_src.begin());
    }

Here is the call graph for this function:

update_w_subsidence()

virtual void ProblemBase::update_w_subsidence ( const amrex::Real &  , amrex::Vector< amrex::Real > &  wbar, amrex::Gpu::DeviceVector< amrex::Real > &  d_wbar, const amrex::Geometry &  geom, std::unique_ptr< amrex::MultiFab > &  z_phys_cc  )

inlinevirtual

Function to update the vertical velocity profile, used to add subsidence source terms for x-mom, y-mom, rho*theta, rho*Q1, and rho*Q2.

TODO: Currently, this is only called by InitData, so there is no time dependence.

Parameters

[in]	time	current time
[out]	wbar	w vel forcing profile
[in]	geom	container for geometric information
[in]	z_phys_cc	height coordinate at cell centers

    {
        if (wbar.empty()) return;
 
        if (z_phys_cc) {
            // use_terrain=1
            amrex::Error("Moisture forcing not defined for "+name()+" problem with terrain");
        } else {
            const int khi              = geom.Domain().bigEnd()[2];
            // const amrex::Real* prob_lo = geom.ProbLo();
            // const auto dx              = geom.CellSize();
            for (int k = 0; k <= khi; k++)
            {
                // const amrex::Real z_cc = prob_lo[2] + (k+0.5)* dx[2];
                // set vertical velocity profile based on time, z_cc here
                wbar[k] = 0.0;
            }
        }
 
        // Copy from host version to device version
        amrex::Gpu::copy(amrex::Gpu::hostToDevice, wbar.begin(), wbar.end(), d_wbar.begin());
    }

Here is the call graph for this function:

Member Data Documentation

base_parms

ProbParmDefaults ProblemBase::base_parms

protected

Referenced by init_base_parms(), and init_uniform().

---

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

• Source/ERF_ProbCommon.H