Radiation Class Reference

#include <ERF_Radiation.H>

Inheritance diagram for Radiation:

Collaboration diagram for Radiation:

Public Member Functions
	Radiation (const int &lev, SolverChoice &sc)
	~Radiation ()
virtual void	Init (const amrex::Geometry &geom, const amrex::BoxArray &ba, amrex::MultiFab *cons_in) override
virtual void	Run (int &level, int &step, double &time, const double &dt, const amrex::BoxArray &ba, amrex::Geometry &geom, amrex::MultiFab *cons_in, amrex::iMultiFab *lmask, amrex::MultiFab *t_surf, amrex::Vector< amrex::MultiFab * > &lsm_input_ptrs, amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs, amrex::MultiFab *qheating_rates, amrex::MultiFab *rad_fluxes, amrex::MultiFab *z_phys, amrex::MultiFab *lat_ptr, amrex::MultiFab *lon_ptr, const bool updated_lsm) override
void	set_grids (int &level, int &step, double &time, const double &dt, const amrex::BoxArray &ba, amrex::Geometry &geom, amrex::MultiFab *cons_in, amrex::iMultiFab *lmask, amrex::MultiFab *t_surf, amrex::Vector< amrex::MultiFab * > &lsm_input_ptrs, amrex::MultiFab *qheating_rates, amrex::MultiFab *rad_fluxes, amrex::MultiFab *z_phys, amrex::MultiFab *lat, amrex::MultiFab *lon, const bool updated_lsm)
void	alloc_buffers ()
void	dealloc_buffers ()
void	mf_to_kokkos_buffers (amrex::iMultiFab *lmask, amrex::MultiFab *t_surf, amrex::Vector< amrex::MultiFab * > &lsm_input_ptrs)
void	kokkos_buffers_to_mf (amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
void	write_rrtmgp_fluxes ()
void	initialize_impl ()
void	run_impl ()
void	finalize_impl (amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
void	rad_run_impl (amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
virtual amrex::Vector< std::string >	get_lsm_input_varnames () override
virtual amrex::Vector< std::string >	get_lsm_output_varnames () override
void	populateDatalogMF ()
virtual void	WriteDataLog (const double &time) override
Public Member Functions inherited from IRadiation
virtual	~IRadiation ()=default
void	setupDataLog ()
void	setDataLogFrequency (const int nstep)
bool	hasDatalog ()

Private Attributes
int	m_lev
int	m_step
double	m_time
double	m_dt
amrex::Geometry	m_geom
amrex::BoxArray	m_ba
bool	m_update_rad = false
bool	m_rad_write_fluxes = false
bool	m_moist = false
bool	m_ice = false
int	m_qi_comp = -1
bool	m_lsm = false
amrex::Vector< std::string >	m_lsm_input_names
amrex::Vector< std::string >	m_lsm_output_names
amrex::Real	m_rad_t_sfc = -1
amrex::MultiFab *	m_cons_in = nullptr
amrex::MultiFab *	m_qheating_rates = nullptr
amrex::MultiFab *	m_rad_fluxes = nullptr
amrex::MultiFab *	m_z_phys = nullptr
amrex::MultiFab *	m_lat = nullptr
amrex::MultiFab *	m_lon = nullptr
amrex::Real	m_lat_cons = amrex::Real(39.809860)
amrex::Real	m_lon_cons = -amrex::Real(98.555183)
amrex::MultiFab	datalog_mf
std::string	rrtmgp_file_path = "."
std::string	rrtmgp_coeffs_sw = "rrtmgp-data-sw-g224-2018-12-04.nc"
std::string	rrtmgp_coeffs_lw = "rrtmgp-data-lw-g256-2018-12-04.nc"
std::string	rrtmgp_cloud_optics_sw = "rrtmgp-cloud-optics-coeffs-sw.nc"
std::string	rrtmgp_cloud_optics_lw = "rrtmgp-cloud-optics-coeffs-lw.nc"
std::string	rrtmgp_coeffs_file_sw
std::string	rrtmgp_coeffs_file_lw
std::string	rrtmgp_cloud_optics_file_sw
std::string	rrtmgp_cloud_optics_file_lw
int	m_ngas = 8
const std::vector< std::string >	m_gas_names
const std::vector< amrex::Real >	m_mol_weight_gas
amrex::Real	m_co2vmr = amrex::Real(388.717e-6)
amrex::Vector< amrex::Real >	m_o3vmr
amrex::Real	m_n2ovmr = amrex::Real(323.141e-9)
amrex::Real	m_covmr = amrex::Real(1.0e-7)
amrex::Real	m_ch4vmr = amrex::Real(1807.851e-9)
amrex::Real	m_o2vmr = amrex::Real(0.209448)
amrex::Real	m_n2vmr = amrex::Real(0.7906)
int	m_o3_size
real1d_k	m_gas_mol_weights
std::vector< std::string >	gas_names_offset
GasConcsK< amrex::Real, layout_t, KokkosDefaultDevice >	m_gas_concs
int	m_ncol
int	m_nlay
amrex::Vector< int >	m_col_offsets
bool	m_do_aerosol_rad = false
bool	m_extra_clnsky_diag = false
bool	m_extra_clnclrsky_diag = false
int	m_orbital_year = -9999
int	m_orbital_mon = -9999
int	m_orbital_day = -9999
int	m_orbital_sec = -9999
bool	m_fixed_orbital_year = false
amrex::Real	m_orbital_eccen = -amrex::Real(9999.)
amrex::Real	m_orbital_obliq = -amrex::Real(9999.)
amrex::Real	m_orbital_mvelp = -amrex::Real(9999.)
amrex::Real	m_fixed_total_solar_irradiance = -amrex::Real(9999.)
amrex::Real	m_fixed_solar_zenith_angle = -amrex::Real(9999.)
int	m_nswbands
int	m_nlwbands
int	m_nswgpts
int	m_nlwgpts
int	m_rad_freq_in_steps = 1
int	m_ncol_chunk_requested = 1024
int	m_ncol_chunk = 1024
int	m_rad_nvar = 12
bool	m_do_subcol_sampling = true
real1d_k	o3_lay
real1d_k	mu0
real1d_k	sfc_alb_dir_vis
real1d_k	sfc_alb_dir_nir
real1d_k	sfc_alb_dif_vis
real1d_k	sfc_alb_dif_nir
real1d_k	sfc_flux_dir_vis
real1d_k	sfc_flux_dir_nir
real1d_k	sfc_flux_dif_vis
real1d_k	sfc_flux_dif_nir
real1d_k	lat
real1d_k	lon
real1d_k	sfc_emis
real1d_k	t_sfc
real1d_k	lw_src
real2d_k	r_lay
real2d_k	p_lay
real2d_k	t_lay
real2d_k	z_del
real2d_k	qv_lay
real2d_k	qc_lay
real2d_k	qi_lay
real2d_k	cldfrac_tot
real2d_k	eff_radius_qc
real2d_k	eff_radius_qi
real2d_k	lwp
real2d_k	iwp
real2d_k	sw_heating
real2d_k	lw_heating
real2d_k	sw_clrsky_heating
real2d_k	lw_clrsky_heating
real2d_k	d_tint
real2d_k	p_lev
real2d_k	t_lev
real2d_k	sw_flux_up
real2d_k	sw_flux_dn
real2d_k	sw_flux_dn_dir
real2d_k	lw_flux_up
real2d_k	lw_flux_dn
real2d_k	sw_clnclrsky_flux_up
real2d_k	sw_clnclrsky_flux_dn
real2d_k	sw_clnclrsky_flux_dn_dir
real2d_k	sw_clrsky_flux_up
real2d_k	sw_clrsky_flux_dn
real2d_k	sw_clrsky_flux_dn_dir
real2d_k	sw_clnsky_flux_up
real2d_k	sw_clnsky_flux_dn
real2d_k	sw_clnsky_flux_dn_dir
real2d_k	lw_clnclrsky_flux_up
real2d_k	lw_clnclrsky_flux_dn
real2d_k	lw_clrsky_flux_up
real2d_k	lw_clrsky_flux_dn
real2d_k	lw_clnsky_flux_up
real2d_k	lw_clnsky_flux_dn
real3d_k	sw_bnd_flux_up
real3d_k	sw_bnd_flux_dn
real3d_k	sw_bnd_flux_dir
real3d_k	sw_bnd_flux_dif
real3d_k	lw_bnd_flux_up
real3d_k	lw_bnd_flux_dn
real2d_k	sfc_alb_dir
real2d_k	sfc_alb_dif
real3d_k	aero_tau_sw
real3d_k	aero_ssa_sw
real3d_k	aero_g_sw
real3d_k	aero_tau_lw

Additional Inherited Members
Protected Attributes inherited from IRadiation
std::unique_ptr< std::fstream >	datalog = nullptr
std::string	datalogname
int	datalog_int = -1

Constructor & Destructor Documentation

Radiation()

Radiation::Radiation	(	const int &	lev,
		SolverChoice &	sc
	)

{
    // Note that Kokkos is now initialized in main.cpp
 
    // Check if we have a valid moisture model
    if (sc.moisture_type != MoistureType::None) { m_moist = true; }
 
    // Cloud-ice support follows the configured moisture-component mapping.
    m_qi_comp = sc.moisture_indices.qi;
    m_ice = (m_qi_comp >= 0);
 
    // Check if we have a land surface model enabled
    if (sc.lsm_type != LandSurfaceType::None) { m_lsm = true; }
 
    // Construct parser object for following reads
    ParmParse pp("erf");
 
    // Must specify a surface temp (LSM can overwrite)
    pp.get("rad_t_sfc", m_rad_t_sfc);
 
    // Radiation timestep, as a number of atm steps
    pp.queryAdd("rad_freq_in_steps", m_rad_freq_in_steps);
 
    // Get nvar if specified
    pp.queryAdd("rad_nvar", m_rad_nvar);
    AMREX_ALWAYS_ASSERT_WITH_MESSAGE(m_rad_nvar >= 0,
        "erf.rad_nvar must be greater than 0. "
        "It controls the amount of memory allocated for temporaries with RRTMGP; "
        "a value of 0 would allocate no memory.");
 
    // Number of columns per RRTMGP chunk (controls peak GPU memory)
    pp.queryAdd("rad_ncol_chunk", m_ncol_chunk_requested);
    AMREX_ALWAYS_ASSERT_WITH_MESSAGE(m_ncol_chunk_requested > 0,
        "erf.rad_ncol_chunk must be a positive integer (default 5000). "
        "It controls the number of columns processed per RRTMGP kernel launch; "
        "a value of 0 or negative would produce an infinite loop.");
    m_ncol_chunk = m_ncol_chunk_requested;
 
    // Flag to write fluxes to plt file
    pp.queryAdd("rad_write_fluxes", m_rad_write_fluxes);
 
    // Do MCICA subcolumn sampling
    pp.queryAdd("rad_do_subcol_sampling", m_do_subcol_sampling);
 
    // Determine orbital year. If orbital_year is negative, use current year
    // from timestamp for orbital year; if non-negative, use provided orbital year
    // for duration of simulation.  Note that this is keyed off the value itself,
    // not off whether the input was present, so that a negative value defers to
    // the timestamp as documented.
    pp.queryAdd("rad_orbital_year", m_orbital_year);
    m_fixed_orbital_year = (m_orbital_year >= 0);
 
    // Get orbital parameters from inputs file
    pp.queryAdd("rad_orbital_eccentricity", m_orbital_eccen);
    pp.queryAdd("rad_orbital_obliquity"   , m_orbital_obliq);
    pp.queryAdd("rad_orbital_mvelp"       , m_orbital_mvelp);
 
    // Get a constant lat/lon for idealized simulations
    pp.queryAdd("rad_cons_lat", m_lat_cons);
    pp.queryAdd("rad_cons_lon", m_lon_cons);
 
    // Value for prescribing an invariant solar constant (i.e. total solar irradiance at
    // TOA).  Used for idealized experiments such as RCE. Disabled when value is less than zero
    pp.queryAdd("fixed_total_solar_irradiance", m_fixed_total_solar_irradiance);
 
    // Determine whether or not we are using a fixed solar zenith angle (positive value)
    pp.queryAdd("fixed_solar_zenith_angle", m_fixed_solar_zenith_angle);
 
    // Get prescribed surface values of greenhouse gases
    pp.queryAdd("co2vmr", m_co2vmr);
    pp.queryarr("o3vmr" , m_o3vmr );
    pp.queryAdd("n2ovmr", m_n2ovmr);
    pp.queryAdd("covmr" , m_covmr );
    pp.queryAdd("ch4vmr", m_ch4vmr);
    pp.queryAdd("o2vmr" , m_o2vmr );
    pp.queryAdd("n2vmr" , m_n2vmr );
 
    // Aerosol forcing hook (not implemented). The aerosol arrays that used to be
    // passed through rrtmgp_main were never populated with real data, so enabling
    // this flag only ever multiplied radiation by zero aerosol optics. The hook is
    // kept so a future SPA/prescribed-aerosol scheme can wire in without touching
    // the ParmParse surface.
    pp.queryAdd("rad_do_aerosol", m_do_aerosol_rad);
    if (m_do_aerosol_rad) {
        amrex::Abort("erf.rad_do_aerosol = true is not supported: aerosol forcing is "
                     "currently not implemented in the ERF RRTMGP interface. The hook "
                     "is retained for a future aerosol coupling; set rad_do_aerosol = "
                     "false (or remove it) to continue.");
    }
 
    // Whether we do extra clean/clear sky calculations
    pp.queryAdd("rad_extra_clnclrsky_diag", m_extra_clnclrsky_diag);
    pp.queryAdd("rad_extra_clnsky_diag"   , m_extra_clnsky_diag);
 
    // Parse path and file names
    pp.queryAdd("rrtmgp_file_path"      , rrtmgp_file_path);
    pp.queryAdd("rrtmgp_coeffs_sw"      , rrtmgp_coeffs_sw  );
    pp.queryAdd("rrtmgp_coeffs_lw"      , rrtmgp_coeffs_lw  );
    pp.queryAdd("rrtmgp_cloud_optics_sw", rrtmgp_cloud_optics_sw);
    pp.queryAdd("rrtmgp_cloud_optics_lw", rrtmgp_cloud_optics_lw);
 
    // Fail early, and with a message that names the missing files and the inputs that
    // control them, rather than letting netCDF report a bare "No such file or directory"
    // when the coefficients are opened just below (or, for the cloud optics files, much
    // later inside rrtmgp_initialize). Every rank does this so that all of them abort.
    check_rrtmgp_data_files(rrtmgp_file_path,
                            {{"rrtmgp_coeffs_sw"      , rrtmgp_coeffs_sw      },
                             {"rrtmgp_coeffs_lw"      , rrtmgp_coeffs_lw      },
                             {"rrtmgp_cloud_optics_sw", rrtmgp_cloud_optics_sw},
                             {"rrtmgp_cloud_optics_lw", rrtmgp_cloud_optics_lw}});
 
    // Append file names to path
    rrtmgp_coeffs_file_sw       = rrtmgp_file_path + "/" + rrtmgp_coeffs_sw;
    rrtmgp_coeffs_file_lw       = rrtmgp_file_path + "/" + rrtmgp_coeffs_lw;
    rrtmgp_cloud_optics_file_sw = rrtmgp_file_path + "/" + rrtmgp_cloud_optics_sw;
    rrtmgp_cloud_optics_file_lw = rrtmgp_file_path + "/" + rrtmgp_cloud_optics_lw;
 
    // Get dimensions from lookup data
    if (ParallelDescriptor::IOProcessor()) {
        auto ncf_sw = ncutils::NCFile::open(rrtmgp_coeffs_file_sw, NC_CLOBBER | NC_NETCDF4);
        m_nswbands = ncf_sw.dim("bnd").len();
        m_nswgpts  = ncf_sw.dim("gpt").len();
        ncf_sw.close();
 
        auto ncf_lw = ncutils::NCFile::open(rrtmgp_coeffs_file_lw, NC_CLOBBER | NC_NETCDF4);
        m_nlwbands = ncf_lw.dim("bnd").len();
        m_nlwgpts  = ncf_lw.dim("gpt").len();
        ncf_lw.close();
    }
    int ioproc = ParallelDescriptor::IOProcessorNumber();  // I/O rank
    ParallelDescriptor::Bcast(&m_nswbands, 1, ioproc);
    ParallelDescriptor::Bcast(&m_nlwbands, 1, ioproc);
    ParallelDescriptor::Bcast(&m_nswgpts,  1, ioproc);
    ParallelDescriptor::Bcast(&m_nlwgpts,  1, ioproc);
 
    // Output for user
    if (lev == 0) {
        Print() << "Radiation interface constructed:\n";
        Print() << "========================================================\n";
        Print() << "Coeff SW file: " << rrtmgp_coeffs_file_sw << "\n";
        Print() << "Coeff LW file: " << rrtmgp_coeffs_file_lw << "\n";
        Print() << "Cloud SW file: " << rrtmgp_cloud_optics_file_sw << "\n";
        Print() << "Cloud LW file: " << rrtmgp_cloud_optics_file_lw << "\n";
        Print() << "Number of short/longwave bands: "
                << m_nswbands << " " << m_nlwbands << "\n";
        Print() << "Number of short/longwave gauss points: "
                << m_nswgpts  << " " << m_nlwgpts  << "\n";
        Print() << "========================================================\n";
    }
}

Here is the call graph for this function:

~Radiation()

Radiation::~Radiation ( )

inline

    {
        // Release k-distribution data and memory pool
        if (rrtmgp::initialized) {
            rrtmgp::rrtmgp_finalize();
        }
        // Note that Kokkos is now finalized in main.cpp
    }

Here is the call graph for this function:

Member Function Documentation

alloc_buffers()

void Radiation::alloc_buffers

(

)

{
    // 1d size (m_ngas)
    const Real* mol_weight_gas_p = m_mol_weight_gas.data();
    const std::string* gas_names_p = m_gas_names.data();
    m_gas_mol_weights = real1d_k("m_gas_mol_weights", m_ngas);
    realHost1d_k m_gas_mol_weights_h("m_gas_mol_weights_h", m_ngas);
    gas_names_offset.clear(); gas_names_offset.resize(m_ngas);
    std::string* gas_names_offset_p = gas_names_offset.data();
    Kokkos::parallel_for(Kokkos::RangePolicy<Kokkos::Serial>(0, m_ngas),
                         [&] (int igas)
    {
        m_gas_mol_weights_h(igas) = mol_weight_gas_p[igas];
        gas_names_offset_p[igas]  = gas_names_p[igas];
    });
    Kokkos::deep_copy(m_gas_mol_weights, m_gas_mol_weights_h);
 
    // 1d size (1 or nlay)
    m_o3_size = m_o3vmr.size();
    AMREX_ALWAYS_ASSERT_WITH_MESSAGE(((m_o3_size==1) || (m_o3_size==m_nlay)),
                                     "O3 VMR array must be length 1 or nlay");
    Real* o3vmr_p = m_o3vmr.data();
    o3_lay = real1d_k("o3_lay", m_o3_size);
    realHost1d_k o3_lay_h("o3_lay_h", m_o3_size);
    Kokkos::parallel_for(Kokkos::RangePolicy<Kokkos::Serial>(0, m_o3_size),
                         [&] (int io3)
    {
        o3_lay_h(io3) = o3vmr_p[io3];
    });
    Kokkos::deep_copy(o3_lay, o3_lay_h);
 
    // 1d size (ncol)
    mu0              = real1d_k("mu0"             , m_ncol);
    sfc_alb_dir_vis  = real1d_k("sfc_alb_dir_vis" , m_ncol);
    sfc_alb_dir_nir  = real1d_k("sfc_alb_dir_nir" , m_ncol);
    sfc_alb_dif_vis  = real1d_k("sfc_alb_dif_vis" , m_ncol);
    sfc_alb_dif_nir  = real1d_k("sfc_alb_dif_nir" , m_ncol);
    sfc_flux_dir_vis = real1d_k("sfc_flux_dir_vis", m_ncol);
    sfc_flux_dir_nir = real1d_k("sfc_flux_dir_nir", m_ncol);
    sfc_flux_dif_vis = real1d_k("sfc_flux_dif_vis", m_ncol);
    sfc_flux_dif_nir = real1d_k("sfc_flux_dif_nir", m_ncol);
    lat              = real1d_k("lat"             , m_ncol);
    lon              = real1d_k("lon"             , m_ncol);
    sfc_emis         = real1d_k("sfc_emis"        , m_ncol);
    t_sfc            = real1d_k("t_sfc"           , m_ncol);
    lw_src           = real1d_k("lw_src"          , m_ncol);
 
    // 2d size (ncol, nlay)
    r_lay         = real2d_k("r_lay"        , m_ncol, m_nlay);
    p_lay         = real2d_k("p_lay"        , m_ncol, m_nlay);
    t_lay         = real2d_k("t_lay"        , m_ncol, m_nlay);
    z_del         = real2d_k("z_del"        , m_ncol, m_nlay);
    qv_lay        = real2d_k("qv"           , m_ncol, m_nlay);
    qc_lay        = real2d_k("qc"           , m_ncol, m_nlay);
    qi_lay        = real2d_k("qi"           , m_ncol, m_nlay);
    cldfrac_tot   = real2d_k("cldfrac_tot"  , m_ncol, m_nlay);
    eff_radius_qc = real2d_k("eff_radius_qc", m_ncol, m_nlay);
    eff_radius_qi = real2d_k("eff_radius_qi", m_ncol, m_nlay);
    lwp           = real2d_k("lwp"          , m_ncol, m_nlay);
    iwp           = real2d_k("iwp"          , m_ncol, m_nlay);
    sw_heating    = real2d_k("sw_heating"   , m_ncol, m_nlay);
    lw_heating    = real2d_k("lw_heating"   , m_ncol, m_nlay);
    if (datalog_int > 0) {
        sw_clrsky_heating = real2d_k("sw_clrsky_heating", m_ncol, m_nlay);
        lw_clrsky_heating = real2d_k("lw_clrsky_heating", m_ncol, m_nlay);
    }
 
    // 2d size (ncol, nlay+1)
    d_tint                   = real2d_k("d_tint"               , m_ncol, m_nlay+1);
    p_lev                    = real2d_k("p_lev"                , m_ncol, m_nlay+1);
    t_lev                    = real2d_k("t_lev"                , m_ncol, m_nlay+1);
 
    sw_flux_up               = real2d_k("sw_flux_up"           , m_ncol, m_nlay+1);
    sw_flux_dn               = real2d_k("sw_flux_dn"           , m_ncol, m_nlay+1);
    sw_flux_dn_dir           = real2d_k("sw_flux_dn_dir"       , m_ncol, m_nlay+1);
 
    lw_flux_up               = real2d_k("lw_flux_up"           , m_ncol, m_nlay+1);
    lw_flux_dn               = real2d_k("lw_flux_dn"           , m_ncol, m_nlay+1);
 
    // Clear-sky flux arrays are always needed
    if (datalog_int > 0) {
        sw_clrsky_flux_up        = real2d_k("sw_clrsky_flux_up"    , m_ncol, m_nlay+1);
        sw_clrsky_flux_dn        = real2d_k("sw_clrsky_flux_dn"    , m_ncol, m_nlay+1);
        sw_clrsky_flux_dn_dir    = real2d_k("sw_clrsky_flux_dn_dir", m_ncol, m_nlay+1);
        lw_clrsky_flux_up        = real2d_k("lw_clrsky_flux_up"    , m_ncol, m_nlay+1);
        lw_clrsky_flux_dn        = real2d_k("lw_clrsky_flux_dn"    , m_ncol, m_nlay+1);
    } else {
        sw_clrsky_flux_up        = real2d_k("sw_clrsky_flux_up"    , m_ncol_chunk, m_nlay+1);
        sw_clrsky_flux_dn        = real2d_k("sw_clrsky_flux_dn"    , m_ncol_chunk, m_nlay+1);
        sw_clrsky_flux_dn_dir    = real2d_k("sw_clrsky_flux_dn_dir", m_ncol_chunk, m_nlay+1);
        lw_clrsky_flux_up        = real2d_k("lw_clrsky_flux_up"    , m_ncol_chunk, m_nlay+1);
        lw_clrsky_flux_dn        = real2d_k("lw_clrsky_flux_dn"    , m_ncol_chunk, m_nlay+1);
    }
 
    // Clean-clear-sky diagnostic fluxes (only when enabled)
    if (m_extra_clnclrsky_diag) {
        sw_clnclrsky_flux_up     = real2d_k("sw_clnclrsky_flux_up"    , m_ncol, m_nlay+1);
        sw_clnclrsky_flux_dn     = real2d_k("sw_clnclrsky_flux_dn"    , m_ncol, m_nlay+1);
        sw_clnclrsky_flux_dn_dir = real2d_k("sw_clnclrsky_flux_dn_dir", m_ncol, m_nlay+1);
        lw_clnclrsky_flux_up     = real2d_k("lw_clnclrsky_flux_up"    , m_ncol, m_nlay+1);
        lw_clnclrsky_flux_dn     = real2d_k("lw_clnclrsky_flux_dn"    , m_ncol, m_nlay+1);
    } else {
        sw_clnclrsky_flux_up     = real2d_k("sw_clnclrsky_flux_up"    , 1, 1);
        sw_clnclrsky_flux_dn     = real2d_k("sw_clnclrsky_flux_dn"    , 1, 1);
        sw_clnclrsky_flux_dn_dir = real2d_k("sw_clnclrsky_flux_dn_dir", 1, 1);
        lw_clnclrsky_flux_up     = real2d_k("lw_clnclrsky_flux_up"    , 1, 1);
        lw_clnclrsky_flux_dn     = real2d_k("lw_clnclrsky_flux_dn"    , 1, 1);
    }
 
    // Clean-sky diagnostic fluxes (only when enabled)
    if (m_extra_clnsky_diag) {
        sw_clnsky_flux_up        = real2d_k("sw_clnsky_flux_up"       , m_ncol, m_nlay+1);
        sw_clnsky_flux_dn        = real2d_k("sw_clnsky_flux_dn"       , m_ncol, m_nlay+1);
        sw_clnsky_flux_dn_dir    = real2d_k("sw_clnsky_flux_dn_dir"   , m_ncol, m_nlay+1);
        lw_clnsky_flux_up        = real2d_k("lw_clnsky_flux_up"       , m_ncol, m_nlay+1);
        lw_clnsky_flux_dn        = real2d_k("lw_clnsky_flux_dn"       , m_ncol, m_nlay+1);
    } else {
        sw_clnsky_flux_up        = real2d_k("sw_clnsky_flux_up"       , 1, 1);
        sw_clnsky_flux_dn        = real2d_k("sw_clnsky_flux_dn"       , 1, 1);
        sw_clnsky_flux_dn_dir    = real2d_k("sw_clnsky_flux_dn_dir"   , 1, 1);
        lw_clnsky_flux_up        = real2d_k("lw_clnsky_flux_up"       , 1, 1);
        lw_clnsky_flux_dn        = real2d_k("lw_clnsky_flux_dn"       , 1, 1);
    }
 
    // 3d size (ncol_chunk, nlay+1, nswbands)
    sw_bnd_flux_up  = real3d_k("sw_bnd_flux_up" , m_ncol_chunk, m_nlay+1, m_nswbands);
    sw_bnd_flux_dn  = real3d_k("sw_bnd_flux_dn" , m_ncol_chunk, m_nlay+1, m_nswbands);
    sw_bnd_flux_dir = real3d_k("sw_bnd_flux_dir", m_ncol_chunk, m_nlay+1, m_nswbands);
    sw_bnd_flux_dif = real3d_k("sw_bnd_flux_dif", m_ncol_chunk, m_nlay+1, m_nswbands);
 
    // 3d size (ncol_chunk, nlay+1, nlwbands)
    lw_bnd_flux_up = real3d_k("lw_bnd_flux_up" , m_ncol_chunk, m_nlay+1, m_nlwbands);
    lw_bnd_flux_dn = real3d_k("lw_bnd_flux_dn" , m_ncol_chunk, m_nlay+1, m_nlwbands);
 
    // 2d size (ncol, nswbands)
    sfc_alb_dir = real2d_k("sfc_alb_dir", m_ncol, m_nswbands);
    sfc_alb_dif = real2d_k("sfc_alb_dif", m_ncol, m_nswbands);
 
    // Aerosol optical properties — allocated only when aerosol coupling is on.
    // The flag gates allocation so today (coupling not implemented, abort fires
    // in the constructor) these stay as empty Views and cost nothing. When a
    // future aerosol scheme populates them, hook up the plumbing into
    // rrtmgp_main as well.
    if (m_do_aerosol_rad) {
        aero_tau_sw = real3d_k("aero_tau_sw", m_ncol, m_nlay, m_nswbands);
        aero_ssa_sw = real3d_k("aero_ssa_sw", m_ncol, m_nlay, m_nswbands);
        aero_g_sw   = real3d_k("aero_g_sw",   m_ncol, m_nlay, m_nswbands);
        aero_tau_lw = real3d_k("aero_tau_lw", m_ncol, m_nlay, m_nlwbands);
    }
}

Here is the call graph for this function:

dealloc_buffers()

void Radiation::dealloc_buffers

(

)

{
    // 1d size (m_ngas)
    m_gas_mol_weights = real1d_k();
 
    // 1d size (1 or nlay)
    o3_lay = real1d_k();
 
    // 1d size (ncol)
    mu0              = real1d_k();
    sfc_alb_dir_vis  = real1d_k();
    sfc_alb_dir_nir  = real1d_k();
    sfc_alb_dif_vis  = real1d_k();
    sfc_alb_dif_nir  = real1d_k();
    sfc_flux_dir_vis = real1d_k();
    sfc_flux_dir_nir = real1d_k();
    sfc_flux_dif_vis = real1d_k();
    sfc_flux_dif_nir = real1d_k();
    lat              = real1d_k();
    lon              = real1d_k();
    sfc_emis         = real1d_k();
    t_sfc            = real1d_k();
    lw_src           = real1d_k();
 
    // 2d size (ncol, nlay)
    r_lay             = real2d_k();
    p_lay             = real2d_k();
    t_lay             = real2d_k();
    z_del             = real2d_k();
    qv_lay            = real2d_k();
    qc_lay            = real2d_k();
    qi_lay            = real2d_k();
    cldfrac_tot       = real2d_k();
    eff_radius_qc     = real2d_k();
    eff_radius_qi     = real2d_k();
    lwp               = real2d_k();
    iwp               = real2d_k();
    sw_heating        = real2d_k();
    lw_heating        = real2d_k();
    sw_clrsky_heating = real2d_k();
    lw_clrsky_heating = real2d_k();
 
    // 2d size (ncol, nlay+1)
    d_tint                   = real2d_k();
    p_lev                    = real2d_k();
    t_lev                    = real2d_k();
    sw_flux_up               = real2d_k();
    sw_flux_dn               = real2d_k();
    sw_flux_dn_dir           = real2d_k();
    lw_flux_up               = real2d_k();
    lw_flux_dn               = real2d_k();
    sw_clnclrsky_flux_up     = real2d_k();
    sw_clnclrsky_flux_dn     = real2d_k();
    sw_clnclrsky_flux_dn_dir = real2d_k();
    sw_clrsky_flux_up        = real2d_k();
    sw_clrsky_flux_dn        = real2d_k();
    sw_clrsky_flux_dn_dir    = real2d_k();
    sw_clnsky_flux_up        = real2d_k();
    sw_clnsky_flux_dn        = real2d_k();
    sw_clnsky_flux_dn_dir    = real2d_k();
    lw_clnclrsky_flux_up     = real2d_k();
    lw_clnclrsky_flux_dn     = real2d_k();
    lw_clrsky_flux_up        = real2d_k();
    lw_clrsky_flux_dn        = real2d_k();
    lw_clnsky_flux_up        = real2d_k();
    lw_clnsky_flux_dn        = real2d_k();
 
    // 3d size (ncol, nlay+1, nswbands)
    sw_bnd_flux_up  = real3d_k();
    sw_bnd_flux_dn  = real3d_k();
    sw_bnd_flux_dir = real3d_k();
    sw_bnd_flux_dif = real3d_k();
 
    // 3d size (ncol, nlay+1, nlwbands)
    lw_bnd_flux_up = real3d_k();
    lw_bnd_flux_dn = real3d_k();
 
    // 2d size (ncol, nswbands)
    sfc_alb_dir = real2d_k();
    sfc_alb_dif = real2d_k();
 
    // Aerosol scaffolding (no-op unless m_do_aerosol_rad enabled allocation above)
    aero_tau_sw = real3d_k();
    aero_ssa_sw = real3d_k();
    aero_g_sw   = real3d_k();
    aero_tau_lw = real3d_k();
}

finalize_impl()

void Radiation::finalize_impl

(

amrex::Vector< amrex::MultiFab * > &

lsm_output_ptrs

)

{
    // Reset gas concentrations (k-dist data persists across steps)
    m_gas_concs.reset();
 
    // Fill the AMReX MFs from Kokkos Views
    kokkos_buffers_to_mf(lsm_output_ptrs);
 
    // Write fluxes if requested
    if (m_rad_write_fluxes) { write_rrtmgp_fluxes(); }
 
    // Fill output data for datalog before deallocating
    if (datalog_int > 0) {
        rrtmgp::compute_heating_rate(sw_clrsky_flux_up, sw_clrsky_flux_dn, r_lay, z_del, sw_clrsky_heating);
        rrtmgp::compute_heating_rate(lw_clrsky_flux_up, lw_clrsky_flux_dn, r_lay, z_del, lw_clrsky_heating);
        Kokkos::fence();
        populateDatalogMF();
    }
 
    // Deallocate the buffer arrays
    dealloc_buffers();
}

Referenced by rad_run_impl().

Here is the call graph for this function:

Here is the caller graph for this function:

get_lsm_input_varnames()

virtual amrex::Vector<std::string> Radiation::get_lsm_input_varnames ( )

inlineoverridevirtual

Reimplemented from IRadiation.

    {
        return m_lsm_input_names;
    }

get_lsm_output_varnames()

virtual amrex::Vector<std::string> Radiation::get_lsm_output_varnames ( )

inlineoverridevirtual

Reimplemented from IRadiation.

    {
        return m_lsm_output_names;
    }

Init()

virtual void Radiation::Init ( const amrex::Geometry &  geom, const amrex::BoxArray &  ba, amrex::MultiFab *  cons_in  )

inlineoverridevirtual

Implements IRadiation.

    {
        // Ensure the boxes span klo -> khi
        int klo = geom.Domain().smallEnd(2);
        int khi = geom.Domain().bigEnd(2);
 
        // Reset vector of offsets for columnar data
        m_nlay = geom.Domain().length(2);
 
        m_ncol = 0;
        m_col_offsets.clear();
        m_col_offsets.resize(int(ba.size()));
        for (amrex::MFIter mfi(*cons_in); mfi.isValid(); ++mfi) {
            const amrex::Box& vbx = mfi.validbox();
            AMREX_ALWAYS_ASSERT_WITH_MESSAGE((klo == vbx.smallEnd(2)) &&
                                             (khi == vbx.bigEnd(2)),
                                             "Vertical decomposition with radiation is not allowed.");
            int nx = vbx.length(0);
            int ny = vbx.length(1);
            m_col_offsets[mfi.index()] = m_ncol;
            m_ncol += nx * ny;
        }
 
        // Recompute the effective chunk size from the user's request rather than
        // clamping the stored value in place. Init() is called again for every
        // MakeNewLevelFromScratch / MakeNewLevelFromCoarse / RemakeLevel on the same
        // Radiation object, so an in-place clamp could only ever shrink: a rank that
        // owned no boxes at some Init (m_ncol == 0) would pin the chunk size at zero
        // and then spin forever in the run_impl chunk loop once a regrid handed it
        // real boxes. Deriving from the immutable request instead lets the chunk size
        // grow back. It is zero only when m_ncol is zero, and a rank with no columns
        // does no radiation work.
        m_ncol_chunk = amrex::min(m_ncol, m_ncol_chunk_requested);
    };

Here is the call graph for this function:

initialize_impl()

void Radiation::initialize_impl

(

)

{
    // Initialize gas concentrations for this step
    m_gas_concs.init(gas_names_offset, m_ncol, m_nlay);
 
    // Load k-distribution and cloud optics data only once.
    // These are static lookup tables that never change.
    // Size the memory pool for the requested chunk size (not the effective one, and
    // not min with the current m_ncol) so that the pool remains valid even if m_ncol
    // grows after regridding/load balancing. The pool is created once and never
    // resized, whereas the effective chunk size is recomputed at every Init() and is
    // bounded above by the request.
    if (!rrtmgp::initialized) {
        gas_concs_t gas_concs_pool;
        gas_concs_pool.init(gas_names_offset, m_ncol_chunk_requested, m_nlay);
        rrtmgp::rrtmgp_initialize(gas_concs_pool,
                                  rrtmgp_coeffs_file_sw      , rrtmgp_coeffs_file_lw      ,
                                  rrtmgp_cloud_optics_file_sw, rrtmgp_cloud_optics_file_lw,
                                  m_rad_nvar);
        gas_concs_pool.reset();
    }
}

Referenced by rad_run_impl().

Here is the call graph for this function:

Here is the caller graph for this function:

kokkos_buffers_to_mf()

void Radiation::kokkos_buffers_to_mf

(

amrex::Vector< amrex::MultiFab * > &

lsm_output_ptrs

)

{
    // Heating rate, fluxes, zenith, lsm ptrs
 
    Table2D<Real,Order::C> p_lay_tab(p_lay.data(), {0,0}, {static_cast<int>(p_lay.extent(0)),static_cast<int>(p_lay.extent(1))});
    Table2D<Real,Order::C> sw_heating_tab(sw_heating.data(), {0,0}, {static_cast<int>(sw_heating.extent(0)),static_cast<int>(sw_heating.extent(1))});
    Table2D<Real,Order::C> lw_heating_tab(lw_heating.data(), {0,0}, {static_cast<int>(lw_heating.extent(0)),static_cast<int>(lw_heating.extent(1))});
    Table2D<Real,Order::C> sw_flux_up_tab(sw_flux_up.data(), {0,0}, {static_cast<int>(sw_flux_up.extent(0)),static_cast<int>(sw_flux_up.extent(1))});
    Table2D<Real,Order::C> sw_flux_dn_tab(sw_flux_dn.data(), {0,0}, {static_cast<int>(sw_flux_dn.extent(0)),static_cast<int>(sw_flux_dn.extent(1))});
    Table2D<Real,Order::C> lw_flux_up_tab(lw_flux_up.data(), {0,0}, {static_cast<int>(lw_flux_up.extent(0)),static_cast<int>(lw_flux_up.extent(1))});
    Table2D<Real,Order::C> lw_flux_dn_tab(lw_flux_dn.data(), {0,0}, {static_cast<int>(lw_flux_dn.extent(0)),static_cast<int>(lw_flux_dn.extent(1))});
 
    TableData<Real,1> sfc_flux_sw_dn; sfc_flux_sw_dn.resize({0}, {static_cast<int>(sw_flux_dn.extent(0))});
    TableData<Real,1> sfc_flux_lw_dn; sfc_flux_lw_dn.resize({0}, {static_cast<int>(lw_flux_dn.extent(0))});
    Table1D<Real> sfc_flux_sw_dn_tab = sfc_flux_sw_dn.table();
    Table1D<Real> sfc_flux_lw_dn_tab = sfc_flux_lw_dn.table();
    Table1D<Real> sfc_flux_sw_dir_vis_tab(sfc_flux_dir_vis.data(), {0}, {static_cast<int>(sfc_flux_dir_vis.extent(0))});
    Table1D<Real> sfc_flux_sw_dir_nir_tab(sfc_flux_dir_nir.data(), {0}, {static_cast<int>(sfc_flux_dir_nir.extent(0))});
    Table1D<Real> sfc_flux_sw_dif_vis_tab(sfc_flux_dif_vis.data(), {0}, {static_cast<int>(sfc_flux_dif_vis.extent(0))});
    Table1D<Real> sfc_flux_sw_dif_nir_tab(sfc_flux_dif_nir.data(), {0}, {static_cast<int>(sfc_flux_dif_nir.extent(0))});
    Table1D<Real>              mu0_tab(mu0.data(),              {0}, {static_cast<int>(mu0.extent(0))});
    Vector<Table1D<Real>> rrtmgp_out_vars = {mu0_tab                , sfc_flux_sw_dn_tab     ,
                                             sfc_flux_sw_dir_vis_tab, sfc_flux_sw_dir_nir_tab,
                                             sfc_flux_sw_dif_vis_tab, sfc_flux_sw_dif_nir_tab,
                                             sfc_flux_lw_dn_tab     };
 
    for (MFIter mfi(*m_cons_in); mfi.isValid(); ++mfi) {
        const auto& vbx      = mfi.validbox();
        const auto& sbx      = makeSlab(vbx,2,vbx.smallEnd(2));
        const int nx         = vbx.length(0);
        const int imin       = vbx.smallEnd(0);
        const int jmin       = vbx.smallEnd(1);
        const int offset     = m_col_offsets[mfi.index()];
        const Array4<Real>& q_arr = m_qheating_rates->array(mfi);
        const Array4<Real>& f_arr = m_rad_fluxes->array(mfi);
        ParallelFor(vbx, [=]
                    AMREX_GPU_DEVICE (int i, int j, int k)
        {
            // map [i,j,k] 0-based to [icol, ilay] 0-based
            const int icol = (j-jmin)*nx + (i-imin) + offset;
            const int ilay = k;
 
            // Temperature heating rate for SW and LW
            q_arr(i,j,k,0) = sw_heating_tab(icol,ilay);
            q_arr(i,j,k,1) = lw_heating_tab(icol,ilay);
 
            // Convert the dT/dz to dTheta/dz
            Real iexner = one/getExnergivenP(Real(p_lay_tab(icol,ilay)), RdoCp);
            q_arr(i,j,k,0) *= iexner;
            q_arr(i,j,k,1) *= iexner;
 
            // Populate the fluxes
            f_arr(i,j,k,0) = sw_flux_up_tab(icol,ilay);
            f_arr(i,j,k,1) = sw_flux_dn_tab(icol,ilay);
            f_arr(i,j,k,2) = lw_flux_up_tab(icol,ilay);
            f_arr(i,j,k,3) = lw_flux_dn_tab(icol,ilay);
 
            if (k==0) {
                sfc_flux_sw_dn_tab(icol) = sw_flux_dn_tab(icol,ilay);
                sfc_flux_lw_dn_tab(icol) = lw_flux_dn_tab(icol,ilay);
            }
        });
        for (int ivar(0); ivar<lsm_output_ptrs.size(); ivar++) {
            if (lsm_output_ptrs[ivar]) {
                auto rrtmgp_for_fill = rrtmgp_out_vars[ivar];
                const Array4<Real>& lsm_out_arr = lsm_output_ptrs[ivar]->array(mfi);
                ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
                {
                    // map [i,j,k] 0-based to [icol, ilay] 0-based
                    const int icol   = (j-jmin)*nx + (i-imin) + offset;
 
                    // export the desired variable at surface
                    lsm_out_arr(i,j,k) = rrtmgp_for_fill(icol);
                });
            } // valid ptr
        } // ivar
    }// mfi
}

Here is the call graph for this function:

mf_to_kokkos_buffers()

void Radiation::mf_to_kokkos_buffers	(	amrex::iMultiFab *	lmask,
		amrex::MultiFab *	t_surf,
		amrex::Vector< amrex::MultiFab * > &	lsm_input_ptrs
	)

{
    Table2D<Real,Order::C> r_lay_tab(r_lay.data(), {0,0}, {static_cast<int>(r_lay.extent(0)),static_cast<int>(r_lay.extent(1))});
    Table2D<Real,Order::C> p_lay_tab(p_lay.data(), {0,0}, {static_cast<int>(p_lay.extent(0)),static_cast<int>(p_lay.extent(1))});
    Table2D<Real,Order::C> t_lay_tab(t_lay.data(), {0,0}, {static_cast<int>(t_lay.extent(0)),static_cast<int>(t_lay.extent(1))});
    Table2D<Real,Order::C> z_del_tab(z_del.data(), {0,0}, {static_cast<int>(z_del.extent(0)),static_cast<int>(z_del.extent(1))});
    Table2D<Real,Order::C> qv_lay_tab(qv_lay.data(), {0,0}, {static_cast<int>(qv_lay.extent(0)),static_cast<int>(qv_lay.extent(1))});
    Table2D<Real,Order::C> qc_lay_tab(qc_lay.data(), {0,0}, {static_cast<int>(qc_lay.extent(0)),static_cast<int>(qc_lay.extent(1))});
    Table2D<Real,Order::C> qi_lay_tab(qi_lay.data(), {0,0}, {static_cast<int>(qi_lay.extent(0)),static_cast<int>(qi_lay.extent(1))});
    Table2D<Real,Order::C> cldfrac_tot_tab(cldfrac_tot.data(), {0,0}, {static_cast<int>(cldfrac_tot.extent(0)),static_cast<int>(cldfrac_tot.extent(1))});
 
    Table2D<Real,Order::C> lwp_tab(lwp.data(), {0,0}, {static_cast<int>(lwp.extent(0)),static_cast<int>(lwp.extent(1))});
    Table2D<Real,Order::C> iwp_tab(iwp.data(), {0,0}, {static_cast<int>(iwp.extent(0)),static_cast<int>(iwp.extent(1))});
    Table2D<Real,Order::C> eff_radius_qc_tab(eff_radius_qc.data(), {0,0}, {static_cast<int>(eff_radius_qc.extent(0)),static_cast<int>(eff_radius_qc.extent(1))});
    Table2D<Real,Order::C> eff_radius_qi_tab(eff_radius_qi.data(), {0,0}, {static_cast<int>(eff_radius_qi.extent(0)),static_cast<int>(eff_radius_qi.extent(1))});
 
    Table2D<Real,Order::C> p_lev_tab(p_lev.data(), {0,0}, {static_cast<int>(p_lev.extent(0)),static_cast<int>(p_lev.extent(1))});
    Table2D<Real,Order::C> t_lev_tab(t_lev.data(), {0,0}, {static_cast<int>(t_lev.extent(0)),static_cast<int>(t_lev.extent(1))});
 
    Table1D<Real> lat_tab(lat.data(), {0}, {static_cast<int>(lat.extent(0))});
    Table1D<Real> lon_tab(lon.data(), {0}, {static_cast<int>(lon.extent(0))});
    Table1D<Real> t_sfc_tab(t_sfc.data(), {0}, {static_cast<int>(t_sfc.extent(0))});
 
    bool moist = m_moist;
    bool ice   = m_ice;
    const int qi_comp = m_qi_comp;
    const bool has_lsm = m_lsm;
    const bool has_lat = m_lat;
    const bool has_lon = m_lon;
    const bool has_surflayer = (t_surf);
    int  ncol  = m_ncol;
    int  nlay  = m_nlay;
    Real dz    = m_geom.CellSize(2);
    Real cons_lat = m_lat_cons;
    Real cons_lon = m_lon_cons;
    Real rad_t_sfc = m_rad_t_sfc;
 
    for (MFIter mfi(*m_cons_in); mfi.isValid(); ++mfi) {
        const auto& vbx  = mfi.validbox();
        const int nx     = vbx.length(0);
        const int imin   = vbx.smallEnd(0);
        const int jmin   = vbx.smallEnd(1);
        const int offset = m_col_offsets[mfi.index()];
        const Array4<const Real>& cons_arr = m_cons_in->const_array(mfi);
        const Array4<const Real>& z_arr    = (m_z_phys) ? m_z_phys->const_array(mfi) :
                                                          Array4<const Real>{};
        const Array4<const Real>& lat_arr  = (m_lat)    ? m_lat->const_array(mfi) :
                                                          Array4<const Real>{};
        const Array4<const Real>& lon_arr  = (m_lon)    ? m_lon->const_array(mfi) :
                                                          Array4<const Real>{};
        ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
        {
            // map [i,j,k] 0-based to [icol, ilay] 0-based
            const int icol   = (j-jmin)*nx + (i-imin) + offset;
            const int ilay   = k;
 
            // EOS input (at CC)
            Real r  = cons_arr(i,j,k,Rho_comp);
            Real rt = cons_arr(i,j,k,RhoTheta_comp);
            Real qv = (moist) ? std::max(cons_arr(i,j,k,RhoQ1_comp)/r,Real(0.)) : Real(0.);
            Real qc = (moist) ? std::max(cons_arr(i,j,k,RhoQ2_comp)/r,Real(0.)) : Real(0.);
            Real qi = (ice && qi_comp >= 0) ? std::max(cons_arr(i,j,k,qi_comp)/r,Real(0.)) : Real(0.);
 
            // EOS avg to z-face
            Real r_lo   = cons_arr(i,j,k-1,Rho_comp);
            Real rt_lo  = cons_arr(i,j,k-1,RhoTheta_comp);
            Real qv_lo  = (moist) ? cons_arr(i,j,k-1,RhoQ1_comp)/r_lo : Real(0.);
            Real dz_k   = (z_arr) ? Real(0.125) * ( (z_arr(i  ,j  ,k+1) - z_arr(i  ,j  ,k))
                                                  + (z_arr(i+1,j  ,k+1) - z_arr(i+1,j  ,k))
                                                  + (z_arr(i  ,j+1,k+1) - z_arr(i  ,j+1,k))
                                                  + (z_arr(i+1,j+1,k+1) - z_arr(i+1,j+1,k)) ) : Real(0.5)*dz; // Dist from w-face to CC at k
            Real dz_km1 = (z_arr) ? Real(0.125) * ( (z_arr(i  ,j  ,k  ) - z_arr(i  ,j  ,k-1))
                                                  + (z_arr(i+1,j  ,k  ) - z_arr(i+1,j  ,k-1))
                                                  + (z_arr(i  ,j+1,k  ) - z_arr(i  ,j+1,k-1))
                                                  + (z_arr(i+1,j+1,k  ) - z_arr(i+1,j+1,k-1)) ) : Real(0.5)*dz; // Dist from w-face to CC at k-1
            // NOTE: Linear interpolation to the w-face weights each CC value by the
            //       distance from the face to the *opposite* CC (inverse distance)
            Real r_avg  = (dz_km1*r  + dz_k*r_lo ) / (dz_k + dz_km1);
            Real rt_avg = (dz_km1*rt + dz_k*rt_lo) / (dz_k + dz_km1);
            Real qv_avg = (dz_km1*qv + dz_k*qv_lo) / (dz_k + dz_km1);
 
            // Views at CC
            r_lay_tab(icol,ilay) = r;
 
            p_lay_tab(icol,ilay) = getPgivenRTh(rt, qv);
            t_lay_tab(icol,ilay) = getTgivenRandRTh(r, rt, qv);
            z_del_tab(icol,ilay) = (z_arr) ? Real(0.25) * ( (z_arr(i  ,j  ,k+1) - z_arr(i  ,j  ,k))
                                                          + (z_arr(i+1,j  ,k+1) - z_arr(i+1,j  ,k))
                                                          + (z_arr(i  ,j+1,k+1) - z_arr(i  ,j+1,k))
                                                          + (z_arr(i+1,j+1,k+1) - z_arr(i+1,j+1,k)) ) : dz;
            qv_lay_tab(icol,ilay) = qv;
            qc_lay_tab(icol,ilay) = qc;
            qi_lay_tab(icol,ilay) = qi;
            cldfrac_tot_tab(icol,ilay) = ((qc+qi)>Real(0.)) ? Real(1.) : Real(0.);
 
            // NOTE: These are populated in 'mixing_ratio_to_cloud_mass'
            lwp_tab(icol,ilay) = Real(0.);
            iwp_tab(icol,ilay) = Real(0.);
 
            // NOTE: These would be populated from P3 (we use the constants in p3_main_impl.hpp)
            // NOTE: These are in units of micron!
            eff_radius_qc_tab(icol,ilay) = (qc>Real(0.)) ? Real(10.0) : Real(0.);
            eff_radius_qi_tab(icol,ilay) = (qi>Real(0.)) ? Real(25.0) : Real(0.);
 
            // Buffers on z-faces (nlay+1)
            p_lev_tab(icol,ilay) = getPgivenRTh(rt_avg, qv_avg);
            t_lev_tab(icol,ilay) = getTgivenRandRTh(r_avg, rt_avg, qv_avg);
            if (ilay==(nlay-1)) {
                Real r_hi  = cons_arr(i,j,k+1,Rho_comp);
                Real rt_hi = cons_arr(i,j,k+1,RhoTheta_comp);
                Real qv_hi = (moist) ? std::max(cons_arr(i,j,k+1,RhoQ1_comp)/r_hi,Real(0.)) : Real(0.);
                Real dz_kp1 = (z_arr) ? Real(0.125) * ( (z_arr(i  ,j  ,k+2) - z_arr(i  ,j  ,k+1))
                                                      + (z_arr(i+1,j  ,k+2) - z_arr(i+1,j  ,k+1))
                                                      + (z_arr(i  ,j+1,k+2) - z_arr(i  ,j+1,k+1))
                                                      + (z_arr(i+1,j+1,k+2) - z_arr(i+1,j+1,k+1)) ) : Real(0.5)*dz; // Dist from w-face to CC at k+1
                r_avg  = (dz_kp1*r  + dz_k*r_hi ) / (dz_k + dz_kp1);
                rt_avg = (dz_kp1*rt + dz_k*rt_hi) / (dz_k + dz_kp1);
                qv_avg = (dz_kp1*qv + dz_k*qv_hi) / (dz_k + dz_kp1);
                p_lev_tab(icol,ilay+1) = getPgivenRTh(rt_avg, qv_avg);
                t_lev_tab(icol,ilay+1) = getTgivenRandRTh(r_avg, rt_avg, qv_avg);
            }
 
            // 1D data structures
            if (k==0) {
                lat_tab(icol) = (has_lat) ? lat_arr(i,j,0) : cons_lat;
                lon_tab(icol) = (has_lon) ? lon_arr(i,j,0) : cons_lon;
            }
 
        });
    } // mfi
 
    // Populate vars LSM would provide
    if (!has_lsm && !has_surflayer) {
        // Parsed surface temp
        Kokkos::deep_copy(t_sfc, rad_t_sfc);
 
        // EAMXX dummy atmos constants
        Kokkos::deep_copy(sfc_alb_dir_vis, Real(0.06));
        Kokkos::deep_copy(sfc_alb_dir_nir, Real(0.06));
        Kokkos::deep_copy(sfc_alb_dif_vis, Real(0.06));
        Kokkos::deep_copy(sfc_alb_dif_nir, Real(0.06));
 
        // AML NOTE: These are not used in current EAMXX, I've left
        //           the code to plug into these if we need it.
        //
        // Current EAMXX constants
        Kokkos::deep_copy(sfc_emis, Real(0.98));
        Kokkos::deep_copy(lw_src  , zero );
    } else {
        Vector<real1d_k> rrtmgp_in_vars = {t_sfc, sfc_emis,
                                           sfc_alb_dir_vis, sfc_alb_dir_nir,
                                           sfc_alb_dif_vis, sfc_alb_dif_nir};
        Vector<Real> rrtmgp_default_vals = {rad_t_sfc, Real(0.98),
                                            Real(0.06), Real(0.06),
                                            Real(0.06), Real(0.06)};
        for (int ivar(0); ivar<lsm_input_ptrs.size(); ivar++) {
            auto rrtmgp_default_val = rrtmgp_default_vals[ivar];
            auto rrtmgp_to_fill_k   = rrtmgp_in_vars[ivar];
            amrex::Table1D<amrex::Real> rrtmgp_to_fill(rrtmgp_to_fill_k.data(),
                                                       0, rrtmgp_to_fill_k.extent(0));
            for (MFIter mfi(*m_cons_in); mfi.isValid(); ++mfi) {
                const auto& vbx  = mfi.validbox();
                const auto& sbx  = makeSlab(vbx,2,vbx.smallEnd(2));
                const int nx     = vbx.length(0);
                const int imin   = vbx.smallEnd(0);
                const int jmin   = vbx.smallEnd(1);
                const int offset = m_col_offsets[mfi.index()];
                const Array4<const int>& lmask_arr   = (lmask)   ? lmask->const_array(mfi) :
                                                                   Array4<const int> {};
                const Array4<const Real>& tsurf_arr  = (t_surf) ? t_surf->const_array(mfi) :
                                                                  Array4<const Real> {};
                const Array4<      Real>& lsm_in_arr = (lsm_input_ptrs[ivar]) ? lsm_input_ptrs[ivar]->array(mfi) :
                                                                                Array4<      Real> {};
                ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
                {
                    // map [i,j,k] 0-based to [icol, ilay] 0-based
                    const int icol   = (j-jmin)*nx + (i-imin) + offset;
 
                    // Check if over land
                    bool is_land = (lmask_arr) ? lmask_arr(i,j,k) : 1;
 
                    // Check if valid LSM data
                    bool valid_lsm_data = (lsm_in_arr && (lsm_in_arr(i,j,k) < lsm_undefined));
 
                    // Have LSM and are over land
                    if (is_land && valid_lsm_data) {
                        rrtmgp_to_fill(icol) = lsm_in_arr(i,j,k);
                    }
                    // We have a SurfLayer (enforce consistency with temperature)
                    else if (tsurf_arr && (ivar==0)) {
                        rrtmgp_to_fill(icol) = tsurf_arr(i,j,k);
                        if (lsm_in_arr) { lsm_in_arr(i,j,k) = tsurf_arr(i,j,k); }
                    }
                    // Use the default value
                    else {
                        rrtmgp_to_fill(icol) = rrtmgp_default_val;
                        if (lsm_in_arr) { lsm_in_arr(i,j,k) = rrtmgp_default_val; }
                    }
                });
            } //mfi
        } // ivar
        Kokkos::deep_copy(lw_src, zero );
    } // have lsm
 
    // Enforce consistency between t_sfc and t_lev at bottom surface
    Kokkos::parallel_for(Kokkos::RangePolicy(0, ncol),
                         KOKKOS_LAMBDA (int icol)
    {
        t_lev_tab(icol,0) = t_sfc_tab(icol);
    });
}

Here is the call graph for this function:

populateDatalogMF()

void Radiation::populateDatalogMF

(

)

{
    Table2D<Real,Order::C> sw_flux_up_tab(sw_flux_up.data(), {0,0}, {static_cast<int>(sw_flux_up.extent(0)),static_cast<int>(sw_flux_up.extent(1))});
    Table2D<Real,Order::C> sw_flux_dn_tab(sw_flux_dn.data(), {0,0}, {static_cast<int>(sw_flux_dn.extent(0)),static_cast<int>(sw_flux_dn.extent(1))});
    Table2D<Real,Order::C> sw_flux_dn_dir_tab(sw_flux_dn_dir.data(), {0,0}, {static_cast<int>(sw_flux_dn_dir.extent(0)),static_cast<int>(sw_flux_dn_dir.extent(1))});
    Table2D<Real,Order::C> lw_flux_up_tab(lw_flux_up.data(), {0,0}, {static_cast<int>(lw_flux_up.extent(0)),static_cast<int>(lw_flux_up.extent(1))});
    Table2D<Real,Order::C> lw_flux_dn_tab(lw_flux_dn.data(), {0,0}, {static_cast<int>(lw_flux_dn.extent(0)),static_cast<int>(lw_flux_dn.extent(1))});
 
    Table2D<Real,Order::C> sw_clrsky_flux_up_tab(sw_clrsky_flux_up.data(), {0,0},
                                                 {static_cast<int>(sw_clrsky_flux_up.extent(0)),static_cast<int>(sw_clrsky_flux_up.extent(1))});
    Table2D<Real,Order::C> sw_clrsky_flux_dn_tab(sw_clrsky_flux_dn.data(), {0,0},
                                                 {static_cast<int>(sw_clrsky_flux_dn.extent(0)),static_cast<int>(sw_clrsky_flux_dn.extent(1))});
    Table2D<Real,Order::C> sw_clrsky_flux_dn_dir_tab(sw_clrsky_flux_dn_dir.data(), {0,0},
                                                 {static_cast<int>(sw_clrsky_flux_dn_dir.extent(0)),static_cast<int>(sw_clrsky_flux_dn_dir.extent(1))});
    Table2D<Real,Order::C> lw_clrsky_flux_up_tab(lw_clrsky_flux_up.data(), {0,0},
                                                 {static_cast<int>(lw_clrsky_flux_up.extent(0)),static_cast<int>(lw_clrsky_flux_up.extent(1))});
    Table2D<Real,Order::C> lw_clrsky_flux_dn_tab(lw_clrsky_flux_dn.data(), {0,0},
                                                 {static_cast<int>(lw_clrsky_flux_dn.extent(0)),static_cast<int>(lw_clrsky_flux_dn.extent(1))});
    Table2D<Real,Order::C> sw_clrsky_heating_tab(sw_clrsky_heating.data(), {0,0},
                                                 {static_cast<int>(sw_clrsky_heating.extent(0)),static_cast<int>(sw_clrsky_heating.extent(1))});
    Table2D<Real,Order::C> lw_clrsky_heating_tab(lw_clrsky_heating.data(), {0,0},
                                                 {static_cast<int>(lw_clrsky_heating.extent(0)),static_cast<int>(lw_clrsky_heating.extent(1))});
    Table2D<Real,Order::C> sw_clnsky_flux_up_tab(sw_clnsky_flux_up.data(), {0,0},
                                                 {static_cast<int>(sw_clnsky_flux_up.extent(0)),static_cast<int>(sw_clnsky_flux_up.extent(1))});
    Table2D<Real,Order::C> sw_clnsky_flux_dn_tab(sw_clnsky_flux_dn.data(), {0,0},
                                                 {static_cast<int>(sw_clnsky_flux_dn.extent(0)),static_cast<int>(sw_clnsky_flux_dn.extent(1))});
    Table2D<Real,Order::C> sw_clnsky_flux_dn_dir_tab(sw_clnsky_flux_dn_dir.data(), {0,0},
                                                 {static_cast<int>(sw_clnsky_flux_dn_dir.extent(0)),static_cast<int>(sw_clnsky_flux_dn_dir.extent(1))});
    Table2D<Real,Order::C> lw_clnsky_flux_up_tab(lw_clnsky_flux_up.data(), {0,0},
                                                 {static_cast<int>(lw_clnsky_flux_up.extent(0)),static_cast<int>(lw_clnsky_flux_up.extent(1))});
    Table2D<Real,Order::C> lw_clnsky_flux_dn_tab(lw_clnsky_flux_dn.data(), {0,0},
                                                 {static_cast<int>(lw_clnsky_flux_dn.extent(0)),static_cast<int>(lw_clnsky_flux_dn.extent(1))});
    Table2D<Real,Order::C> sw_clnclrsky_flux_up_tab(sw_clnclrsky_flux_up.data(), {0,0},
                                                 {static_cast<int>(sw_clnclrsky_flux_up.extent(0)),static_cast<int>(sw_clnclrsky_flux_up.extent(1))});
    Table2D<Real,Order::C> sw_clnclrsky_flux_dn_tab(sw_clnclrsky_flux_dn.data(), {0,0},
                                                 {static_cast<int>(sw_clnclrsky_flux_dn.extent(0)),static_cast<int>(sw_clnclrsky_flux_dn.extent(1))});
    Table2D<Real,Order::C> sw_clnclrsky_flux_dn_dir_tab(sw_clnclrsky_flux_dn_dir.data(), {0,0},
                                                 {static_cast<int>(sw_clnclrsky_flux_dn_dir.extent(0)),static_cast<int>(sw_clnclrsky_flux_dn_dir.extent(1))});
    Table2D<Real,Order::C> lw_clnclrsky_flux_up_tab(lw_clnclrsky_flux_up.data(), {0,0},
                                                 {static_cast<int>(lw_clnclrsky_flux_up.extent(0)),static_cast<int>(lw_clnclrsky_flux_up.extent(1))});
    Table2D<Real,Order::C> lw_clnclrsky_flux_dn_tab(lw_clnclrsky_flux_dn.data(), {0,0},
                                                 {static_cast<int>(lw_clnclrsky_flux_dn.extent(0)),static_cast<int>(lw_clnclrsky_flux_dn.extent(1))});
 
    Table1D<Real>              mu0_tab(mu0.data(),              {0}, {static_cast<int>(mu0.extent(0))});
 
    auto extra_clnsky_diag    = m_extra_clnsky_diag;
    auto extra_clnclrsky_diag = m_extra_clnclrsky_diag;
 
    for (MFIter mfi(datalog_mf); mfi.isValid(); ++mfi) {
        const auto& vbx      = mfi.validbox();
        const int nx         = vbx.length(0);
        const int imin       = vbx.smallEnd(0);
        const int jmin       = vbx.smallEnd(1);
        const int offset     = m_col_offsets[mfi.index()];
        const Array4<Real>& dst_arr = datalog_mf.array(mfi);
        const Array4<Real>&   q_arr = m_qheating_rates->array(mfi);
        ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
        {
            // map [i,j,k] 0-based to [icol, ilay] 0-based
            const int icol = (j-jmin)*nx + (i-imin) + offset;
            const int ilay = k;
 
            dst_arr(i,j,k,0) = q_arr(i, j, k, 0);
            dst_arr(i,j,k,1) = q_arr(i, j, k, 1);
 
            // SW and LW fluxes
            dst_arr(i,j,k,2) = sw_flux_up_tab(icol,ilay);
            dst_arr(i,j,k,3) = sw_flux_dn_tab(icol,ilay);
            dst_arr(i,j,k,4) = sw_flux_dn_dir_tab(icol,ilay);
            dst_arr(i,j,k,5) = lw_flux_up_tab(icol,ilay);
            dst_arr(i,j,k,6) = lw_flux_dn_tab(icol,ilay);
 
            // Cosine zenith angle
            dst_arr(i,j,k,7) = mu0_tab(icol);
 
            // Clear sky heating rates and fluxes:
            dst_arr(i,j,k,8) = sw_clrsky_heating_tab(icol, ilay);
            dst_arr(i,j,k,9) = lw_clrsky_heating_tab(icol, ilay);
 
            dst_arr(i,j,k,10) = sw_clrsky_flux_up_tab(icol,ilay);
            dst_arr(i,j,k,11) = sw_clrsky_flux_dn_tab(icol,ilay);
            dst_arr(i,j,k,12) = sw_clrsky_flux_dn_dir_tab(icol,ilay);
            dst_arr(i,j,k,13) = lw_clrsky_flux_up_tab(icol,ilay);
            dst_arr(i,j,k,14) = lw_clrsky_flux_dn_tab(icol,ilay);
 
            // Clean sky fluxes:
            if (extra_clnsky_diag) {
                dst_arr(i,j,k,15) = sw_clnsky_flux_up_tab(icol,ilay);
                dst_arr(i,j,k,16) = sw_clnsky_flux_dn_tab(icol,ilay);
                dst_arr(i,j,k,17) = sw_clnsky_flux_dn_dir_tab(icol,ilay);
                dst_arr(i,j,k,18) = lw_clnsky_flux_up_tab(icol,ilay);
                dst_arr(i,j,k,19) = lw_clnsky_flux_dn_tab(icol,ilay);
            }
 
            // Clean-clear sky fluxes:
            if (extra_clnclrsky_diag) {
                dst_arr(i,j,k,20) = sw_clnclrsky_flux_up_tab(icol,ilay);
                dst_arr(i,j,k,21) = sw_clnclrsky_flux_dn_tab(icol,ilay);
                dst_arr(i,j,k,22) = sw_clnclrsky_flux_dn_dir_tab(icol,ilay);
                dst_arr(i,j,k,23) = lw_clnclrsky_flux_up_tab(icol,ilay);
                dst_arr(i,j,k,24) = lw_clnclrsky_flux_dn_tab(icol,ilay);
            }
        });
   }
}

Here is the call graph for this function:

rad_run_impl()

void Radiation::rad_run_impl ( amrex::Vector< amrex::MultiFab * > &  lsm_output_ptrs )

inline

    {
        if (m_update_rad) {
            amrex::Print() << "Radiation advancing level " << m_lev << " at (YY-MM-DD SS) " << m_orbital_year << '-'
                << m_orbital_mon << '-' << m_orbital_day << ' ' << m_orbital_sec << " ...";
            this->initialize_impl();
            this->run_impl();
            this->finalize_impl(lsm_output_ptrs);
            amrex::Print() << "DONE\n";
        }
    }

Referenced by Run().

Here is the call graph for this function:

Here is the caller graph for this function:

Run()

virtual void Radiation::Run ( int &  level, int &  step, double &  time, const double &  dt, const amrex::BoxArray &  ba, amrex::Geometry &  geom, amrex::MultiFab *  cons_in, amrex::iMultiFab *  lmask, amrex::MultiFab *  t_surf, amrex::Vector< amrex::MultiFab * > &  lsm_input_ptrs, amrex::Vector< amrex::MultiFab * > &  lsm_output_ptrs, amrex::MultiFab *  qheating_rates, amrex::MultiFab *  rad_fluxes, amrex::MultiFab *  z_phys, amrex::MultiFab *  lat_ptr, amrex::MultiFab *  lon_ptr, const bool  updated_lsm  )

inlineoverridevirtual

Implements IRadiation.

    {
        set_grids(level, step, time, dt, ba, geom,
                  cons_in, lmask, t_surf,
                  lsm_input_ptrs, qheating_rates,
                  rad_fluxes, z_phys, lat_ptr, lon_ptr,
                  updated_lsm);
        rad_run_impl(lsm_output_ptrs);
    }

Here is the call graph for this function:

run_impl()

void Radiation::run_impl

(

)

{
    // A rank that owns no boxes on this level has no columns and therefore no
    // radiation work to do. Bail out before the chunk loop; there are no MPI
    // collectives in this routine, so returning early cannot deadlock.
    if (m_ncol == 0) { return; }
 
    // Local copies
    const auto ncol     = m_ncol;
    const auto nlay     = m_nlay;
    const auto nswbands = m_nswbands;
 
    // Compute orbital parameters; these are used both for computing
    // the solar zenith angle and also for computing total solar
    // irradiance scaling (tsi_scaling).
    double obliqr, lambm0, mvelpp;
    // Any of eccen/obliq/mvelp that the user set (i.e. is non-negative) is used
    // as is by orbital_params; the rest are computed from the orbital year.
    int  orbital_year = m_orbital_year;
    double eccen      = m_orbital_eccen;
    double obliq      = m_orbital_obliq;
    double mvelp      = m_orbital_mvelp;
    orbital_params(orbital_year, eccen, obliq,
                   mvelp, obliqr, lambm0, mvelpp);
 
    // Use the orbital parameters to calculate the solar declination and eccentricity factor
    double delta, eccf;
    // Want day + fraction; calday 1 == Jan 1 0Z
    static constexpr double dpy[] = {zero  ,  Real(31.0),  Real(59.0),  Real(90.0), Real(120.0), Real(151.0),
                                     Real(181.0), Real(212.0), Real(243.0), Real(273.0), Real(304.0), Real(334.0)};
    bool leap = (m_orbital_year % 4 == 0 && (!(m_orbital_year % 100 == 0) || (m_orbital_year % 400 == 0))) ? true : false;
    double calday = one + dpy[m_orbital_mon-1] + (m_orbital_day-one) + m_orbital_sec/Real(86400.0);
    // add extra day if leap year and past February
    if (leap && m_orbital_mon>2) { calday += one; }
    orbital_decl(calday, eccen, mvelpp, lambm0, obliqr, delta, eccf);
 
    // Overwrite eccf if using a fixed solar constant.
    auto fixed_total_solar_irradiance = m_fixed_total_solar_irradiance;
    if (fixed_total_solar_irradiance >= 0){
       eccf = fixed_total_solar_irradiance/Real(1360.9);
    }
 
    // Precompute volume mixing ratio (VMR) for all gases
    //
    // H2O is obtained from qv.
    // O3 may be a constant or a 1D vector
    // All other comps are set to constants for now
    Vector<real2d_k> vmr_full_vec(m_ngas);
    for (int igas(0); igas < m_ngas; ++igas) {
        auto name = m_gas_names[igas];
        vmr_full_vec[igas] = real2d_k("vmr_full_" + name, ncol, nlay);
        auto tmp2d = vmr_full_vec[igas];
        auto gas_mol_weight = m_mol_weight_gas[igas];
        if (name == "H2O") {
            auto qv_lay_d = qv_lay;
            Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
                                 KOKKOS_LAMBDA (int icol, int ilay)
            {
                tmp2d(icol,ilay) = qv_lay_d(icol,ilay) * mwdair/gas_mol_weight;
            });
        } else if (name == "CO2") {
            Kokkos::deep_copy(tmp2d, m_co2vmr);
        } else if (name == "O3")  {
            if (m_o3_size==1) {
                Kokkos::deep_copy(tmp2d, m_o3vmr[0] );
            } else {
                auto o3_lay_d = o3_lay;
                Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
                                     KOKKOS_LAMBDA (int icol, int ilay)
                {
                    tmp2d(icol,ilay) = o3_lay_d(ilay);
                });
            }
        } else if (name == "N2O") {
            Kokkos::deep_copy(tmp2d, m_n2ovmr);
        } else if (name == "CO")  {
            Kokkos::deep_copy(tmp2d, m_covmr );
        } else if (name == "CH4") {
            Kokkos::deep_copy(tmp2d, m_ch4vmr);
        } else if (name == "O2") {
            Kokkos::deep_copy(tmp2d, m_o2vmr );
        } else if (name == "N2") {
            Kokkos::deep_copy(tmp2d, m_n2vmr );
        } else {
            Abort("Radiation: Unknown gas component.");
        }
 
        // Populate GasConcs object
        m_gas_concs.set_vmr(name, tmp2d);
        Kokkos::fence();
    }
 
    // Populate mu0 1D array
    // This must be done on HOST and copied to device.
    auto h_mu0 = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), mu0);
    if (m_fixed_solar_zenith_angle > 0) {
        Kokkos::deep_copy(h_mu0, m_fixed_solar_zenith_angle);
    } else {
        auto h_lat = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), lat);
        auto h_lon = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), lon);
        double dt  = double(m_dt);
        auto rad_freq_in_steps = m_rad_freq_in_steps;
        Kokkos::parallel_for(Kokkos::RangePolicy<Kokkos::Serial>(0, ncol),
                             [&,PI_d=PI] (int icol)
        {
            // Convert lat/lon to radians
            double lat_col = h_lat(icol)*PI_d/Real(180.0);
            double lon_col = h_lon(icol)*PI_d/Real(180.0);
            double lcalday = calday;
            double ldelta  = delta;
            double dt_avg  = static_cast<double>(rad_freq_in_steps) * dt;
            h_mu0(icol)    = Real(orbital_cos_zenith(lcalday, lat_col, lon_col, ldelta, dt_avg));
        });
    }
    Kokkos::deep_copy(mu0, h_mu0);
 
    // Compute layer cloud mass per unit area (populates lwp/iwp)
    rrtmgp::mixing_ratio_to_cloud_mass(qc_lay, cldfrac_tot, r_lay, z_del, lwp);
    rrtmgp::mixing_ratio_to_cloud_mass(qi_lay, cldfrac_tot, r_lay, z_del, iwp);
 
    // Convert to g/m2 (needed by RRTMGP)
    Table2D<Real,Order::C> lwp_tab(lwp.data(), {0,0}, {static_cast<int>(lwp.extent(0)),static_cast<int>(lwp.extent(1))});
    Table2D<Real,Order::C> iwp_tab(iwp.data(), {0,0}, {static_cast<int>(iwp.extent(0)),static_cast<int>(iwp.extent(1))});
    Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
                         KOKKOS_LAMBDA (int icol, int ilay)
    {
        lwp_tab(icol,ilay) *= Real(1.e3);
        iwp_tab(icol,ilay) *= Real(1.e3);
    });
 
    // -----------------------------------------------------------------------
    // Process radiation in column chunks to limit peak GPU memory.
    // Radiation columns are independent (no horizontal coupling), so
    // chunking produces bit-identical results.
    // -----------------------------------------------------------------------
    const int ncol_chunk = std::min(m_ncol_chunk, ncol);
    const int kbot = 0;
 
    for (int col_s = 0; col_s < ncol; col_s += ncol_chunk) {
        const int ncol_c = std::min(ncol_chunk, ncol - col_s);
        const int col_e  = col_s + ncol_c;
        auto cr = std::make_pair(col_s, col_e);
 
        // --- Chunk subviews: 1D (ncol) ---
        real1d_k mu0_c              (mu0.data()              + col_s, ncol_c);
        real1d_k sfc_alb_dir_vis_c  (sfc_alb_dir_vis.data()  + col_s, ncol_c);
        real1d_k sfc_alb_dir_nir_c  (sfc_alb_dir_nir.data()  + col_s, ncol_c);
        real1d_k sfc_alb_dif_vis_c  (sfc_alb_dif_vis.data()  + col_s, ncol_c);
        real1d_k sfc_alb_dif_nir_c  (sfc_alb_dif_nir.data()  + col_s, ncol_c);
        real1d_k sfc_flux_dir_vis_c (sfc_flux_dir_vis.data() + col_s, ncol_c);
        real1d_k sfc_flux_dir_nir_c (sfc_flux_dir_nir.data() + col_s, ncol_c);
        real1d_k sfc_flux_dif_vis_c (sfc_flux_dif_vis.data() + col_s, ncol_c);
        real1d_k sfc_flux_dif_nir_c (sfc_flux_dif_nir.data() + col_s, ncol_c);
        real1d_k t_sfc_c            (t_sfc.data()            + col_s, ncol_c);
        real1d_k sfc_emis_c         (sfc_emis.data()         + col_s, ncol_c);
        real1d_k lw_src_c           (lw_src.data()           + col_s, ncol_c);
 
        // --- Chunk subviews: 2D (ncol, nlay) via LayoutRight pointer offset ---
        const int stride2_nlay   = nlay;
        const int stride2_nlayp1 = nlay + 1;
        real2d_k p_lay_c        (p_lay.data()         + col_s*stride2_nlay, ncol_c, nlay);
        real2d_k t_lay_c        (t_lay.data()         + col_s*stride2_nlay, ncol_c, nlay);
        real2d_k r_lay_c        (r_lay.data()         + col_s*stride2_nlay, ncol_c, nlay);
        real2d_k z_del_c        (z_del.data()         + col_s*stride2_nlay, ncol_c, nlay);
        real2d_k lwp_c          (lwp.data()           + col_s*stride2_nlay, ncol_c, nlay);
        real2d_k iwp_c          (iwp.data()           + col_s*stride2_nlay, ncol_c, nlay);
        real2d_k eff_radius_qc_c(eff_radius_qc.data() + col_s*stride2_nlay, ncol_c, nlay);
        real2d_k eff_radius_qi_c(eff_radius_qi.data() + col_s*stride2_nlay, ncol_c, nlay);
        real2d_k cldfrac_tot_c  (cldfrac_tot.data()   + col_s*stride2_nlay, ncol_c, nlay);
        real2d_k sw_heating_c   (sw_heating.data()    + col_s*stride2_nlay, ncol_c, nlay);
        real2d_k lw_heating_c   (lw_heating.data()    + col_s*stride2_nlay, ncol_c, nlay);
 
        // --- Chunk subviews: 2D (ncol, nlay+1) ---
        real2d_k p_lev_c                   (p_lev.data()                 + col_s*stride2_nlayp1, ncol_c, nlay+1);
        real2d_k t_lev_c                   (t_lev.data()                 + col_s*stride2_nlayp1, ncol_c, nlay+1);
        real2d_k sw_flux_up_c              (sw_flux_up.data()            + col_s*stride2_nlayp1, ncol_c, nlay+1);
        real2d_k sw_flux_dn_c              (sw_flux_dn.data()            + col_s*stride2_nlayp1, ncol_c, nlay+1);
        real2d_k sw_flux_dn_dir_c          (sw_flux_dn_dir.data()        + col_s*stride2_nlayp1, ncol_c, nlay+1);
        real2d_k lw_flux_up_c              (lw_flux_up.data()            + col_s*stride2_nlayp1, ncol_c, nlay+1);
        real2d_k lw_flux_dn_c              (lw_flux_dn.data()            + col_s*stride2_nlayp1, ncol_c, nlay+1);
        // Clear-sky flux subviews (always active)
        // NOTE: once on m_ncol_chunk if not writing a datalog
        real2d_k sw_clrsky_flux_up_c, sw_clrsky_flux_dn_c, sw_clrsky_flux_dn_dir_c;
        real2d_k lw_clrsky_flux_up_c, lw_clrsky_flux_dn_c;
        if (datalog_int > 0) {
            sw_clrsky_flux_up_c     = real2d_k(sw_clrsky_flux_up.data()     + col_s*stride2_nlayp1, ncol_c, nlay+1);
            sw_clrsky_flux_dn_c     = real2d_k(sw_clrsky_flux_dn.data()     + col_s*stride2_nlayp1, ncol_c, nlay+1);
            sw_clrsky_flux_dn_dir_c = real2d_k(sw_clrsky_flux_dn_dir.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
            lw_clrsky_flux_up_c     = real2d_k(lw_clrsky_flux_up.data()     + col_s*stride2_nlayp1, ncol_c, nlay+1);
            lw_clrsky_flux_dn_c     = real2d_k(lw_clrsky_flux_dn.data()     + col_s*stride2_nlayp1, ncol_c, nlay+1);
        } else {
            sw_clrsky_flux_up_c     = real2d_k(sw_clrsky_flux_up.data()     , ncol_c, nlay+1);
            sw_clrsky_flux_dn_c     = real2d_k(sw_clrsky_flux_dn.data()     , ncol_c, nlay+1);
            sw_clrsky_flux_dn_dir_c = real2d_k(sw_clrsky_flux_dn_dir.data() , ncol_c, nlay+1);
            lw_clrsky_flux_up_c     = real2d_k(lw_clrsky_flux_up.data()     , ncol_c, nlay+1);
            lw_clrsky_flux_dn_c     = real2d_k(lw_clrsky_flux_dn.data()     , ncol_c, nlay+1);
        }
 
        // Diagnostic flux subviews (placeholder when disabled)
        real2d_k sw_clnclrsky_flux_up_c, sw_clnclrsky_flux_dn_c, sw_clnclrsky_flux_dn_dir_c;
        real2d_k lw_clnclrsky_flux_up_c, lw_clnclrsky_flux_dn_c;
        if (m_extra_clnclrsky_diag) {
            sw_clnclrsky_flux_up_c     = real2d_k(sw_clnclrsky_flux_up.data()     + col_s*stride2_nlayp1, ncol_c, nlay+1);
            sw_clnclrsky_flux_dn_c     = real2d_k(sw_clnclrsky_flux_dn.data()     + col_s*stride2_nlayp1, ncol_c, nlay+1);
            sw_clnclrsky_flux_dn_dir_c = real2d_k(sw_clnclrsky_flux_dn_dir.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
            lw_clnclrsky_flux_up_c     = real2d_k(lw_clnclrsky_flux_up.data()     + col_s*stride2_nlayp1, ncol_c, nlay+1);
            lw_clnclrsky_flux_dn_c     = real2d_k(lw_clnclrsky_flux_dn.data()     + col_s*stride2_nlayp1, ncol_c, nlay+1);
        } else {
            sw_clnclrsky_flux_up_c     = real2d_k("sw_clnclrsky_flux_up_c"    , 1, 1);
            sw_clnclrsky_flux_dn_c     = real2d_k("sw_clnclrsky_flux_dn_c"    , 1, 1);
            sw_clnclrsky_flux_dn_dir_c = real2d_k("sw_clnclrsky_flux_dn_dir_c", 1, 1);
            lw_clnclrsky_flux_up_c     = real2d_k("lw_clnclrsky_flux_up_c"    , 1, 1);
            lw_clnclrsky_flux_dn_c     = real2d_k("lw_clnclrsky_flux_dn_c"    , 1, 1);
        }
 
        real2d_k sw_clnsky_flux_up_c, sw_clnsky_flux_dn_c, sw_clnsky_flux_dn_dir_c;
        real2d_k lw_clnsky_flux_up_c, lw_clnsky_flux_dn_c;
        if (m_extra_clnsky_diag) {
            sw_clnsky_flux_up_c     = real2d_k(sw_clnsky_flux_up.data()     + col_s*stride2_nlayp1, ncol_c, nlay+1);
            sw_clnsky_flux_dn_c     = real2d_k(sw_clnsky_flux_dn.data()     + col_s*stride2_nlayp1, ncol_c, nlay+1);
            sw_clnsky_flux_dn_dir_c = real2d_k(sw_clnsky_flux_dn_dir.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
            lw_clnsky_flux_up_c     = real2d_k(lw_clnsky_flux_up.data()     + col_s*stride2_nlayp1, ncol_c, nlay+1);
            lw_clnsky_flux_dn_c     = real2d_k(lw_clnsky_flux_dn.data()     + col_s*stride2_nlayp1, ncol_c, nlay+1);
        } else {
            sw_clnsky_flux_up_c     = real2d_k("sw_clnsky_flux_up_c"    , 1, 1);
            sw_clnsky_flux_dn_c     = real2d_k("sw_clnsky_flux_dn_c"    , 1, 1);
            sw_clnsky_flux_dn_dir_c = real2d_k("sw_clnsky_flux_dn_dir_c", 1, 1);
            lw_clnsky_flux_up_c     = real2d_k("lw_clnsky_flux_up_c"    , 1, 1);
            lw_clnsky_flux_dn_c     = real2d_k("lw_clnsky_flux_dn_c"    , 1, 1);
        }
 
        // --- Chunk subviews: 2D (ncol, nswbands) ---
        real2d_k sfc_alb_dir_c(sfc_alb_dir.data() + col_s*nswbands, ncol_c, nswbands);
        real2d_k sfc_alb_dif_c(sfc_alb_dif.data() + col_s*nswbands, ncol_c, nswbands);
 
        // --- Chunk subviews: 3D (ncol, nlay+1, nbands) ---
        // NOTE: Allocate these once on m_ncol_chunk and use what we need in the chunk loop
        real3d_k sw_bnd_flux_up_c (sw_bnd_flux_up.data() , ncol_c, nlay+1, nswbands);
        real3d_k sw_bnd_flux_dn_c (sw_bnd_flux_dn.data() , ncol_c, nlay+1, nswbands);
        real3d_k sw_bnd_flux_dir_c(sw_bnd_flux_dir.data(), ncol_c, nlay+1, nswbands);
        real3d_k sw_bnd_flux_dif_c(sw_bnd_flux_dif.data(), ncol_c, nlay+1, nswbands);
        real3d_k lw_bnd_flux_up_c (lw_bnd_flux_up.data() , ncol_c, nlay+1, m_nlwbands);
        real3d_k lw_bnd_flux_dn_c (lw_bnd_flux_dn.data() , ncol_c, nlay+1, m_nlwbands);
 
        // --- Create chunk gas concentrations by subsetting from pre-fetched VMR ---
        gas_concs_t gas_concs_c;
        gas_concs_c.init(gas_names_offset, ncol_c, nlay);
        for (int igas = 0; igas < m_ngas; ++igas) {
            real2d_k vmr_c("vmr_c", ncol_c, nlay);
            auto vmr_full = vmr_full_vec[igas];
            auto cs = col_s;
            Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol_c, nlay}),
                                 KOKKOS_LAMBDA (int i, int j) {
                vmr_c(i, j) = vmr_full(cs + i, j);
            });
            gas_concs_c.set_vmr(m_gas_names[igas], vmr_c);
        }
 
        // Expand surface albedos along nswbands for this chunk
        rrtmgp::compute_band_by_band_surface_albedos(ncol_c, nswbands,
                                                     sfc_alb_dir_vis_c, sfc_alb_dir_nir_c,
                                                     sfc_alb_dif_vis_c, sfc_alb_dif_nir_c,
                                                     sfc_alb_dir_c    , sfc_alb_dif_c);
 
        // Run RRTMGP driver for this column chunk
        rrtmgp::rrtmgp_main(ncol_c, m_nlay,
                            p_lay_c, t_lay_c,
                            p_lev_c, t_lev_c,
                            gas_concs_c,
                            sfc_alb_dir_c, sfc_alb_dif_c, mu0_c,
                            t_sfc_c, sfc_emis_c, lw_src_c,
                            lwp_c, iwp_c, eff_radius_qc_c, eff_radius_qi_c, cldfrac_tot_c,
                            sw_flux_up_c, sw_flux_dn_c, sw_flux_dn_dir_c,
                            lw_flux_up_c, lw_flux_dn_c,
                            sw_clnclrsky_flux_up_c, sw_clnclrsky_flux_dn_c, sw_clnclrsky_flux_dn_dir_c,
                            sw_clrsky_flux_up_c, sw_clrsky_flux_dn_c, sw_clrsky_flux_dn_dir_c,
                            sw_clnsky_flux_up_c, sw_clnsky_flux_dn_c, sw_clnsky_flux_dn_dir_c,
                            lw_clnclrsky_flux_up_c, lw_clnclrsky_flux_dn_c,
                            lw_clrsky_flux_up_c, lw_clrsky_flux_dn_c,
                            lw_clnsky_flux_up_c, lw_clnsky_flux_dn_c,
                            sw_bnd_flux_up_c, sw_bnd_flux_dn_c, sw_bnd_flux_dir_c,
                            lw_bnd_flux_up_c, lw_bnd_flux_dn_c,
                            eccf, m_extra_clnclrsky_diag, m_extra_clnsky_diag);
 
        // Compute heating rates for this chunk
        rrtmgp::compute_heating_rate(sw_flux_up_c, sw_flux_dn_c, r_lay_c, z_del_c, sw_heating_c);
        rrtmgp::compute_heating_rate(lw_flux_up_c, lw_flux_dn_c, r_lay_c, z_del_c, lw_heating_c);
 
        // Compute diffuse band fluxes and broadband surface fluxes for this chunk
        Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<3>>({0, 0, 0}, {ncol_c, nlay+1, nswbands}),
                             KOKKOS_LAMBDA (int icol, int ilay, int ibnd)
        {
            sw_bnd_flux_dif_c(icol,ilay,ibnd) = sw_bnd_flux_dn_c(icol,ilay,ibnd) - sw_bnd_flux_dir_c(icol,ilay,ibnd);
        });
        rrtmgp::compute_broadband_surface_fluxes(ncol_c, kbot, nswbands,
                                                 sw_bnd_flux_dir_c , sw_bnd_flux_dif_c ,
                                                 sfc_flux_dir_vis_c, sfc_flux_dir_nir_c,
                                                 sfc_flux_dif_vis_c, sfc_flux_dif_nir_c);
 
        gas_concs_c.reset();
    } // end column chunk loop
}

Referenced by rad_run_impl().

Here is the call graph for this function:

Here is the caller graph for this function:

set_grids()

void Radiation::set_grids	(	int &	level,
		int &	step,
		double &	time,
		const double &	dt,
		const amrex::BoxArray &	ba,
		amrex::Geometry &	geom,
		amrex::MultiFab *	cons_in,
		amrex::iMultiFab *	lmask,
		amrex::MultiFab *	t_surf,
		amrex::Vector< amrex::MultiFab * > &	lsm_input_ptrs,
		amrex::MultiFab *	qheating_rates,
		amrex::MultiFab *	rad_fluxes,
		amrex::MultiFab *	z_phys,
		amrex::MultiFab *	lat,
		amrex::MultiFab *	lon,
		const bool	updated_lsm
	)

{
    // Set data members that may change
    m_lev            = level;
    m_step           = step;
    m_time           = time;
    m_dt             = dt;
    m_geom           = geom;
    m_cons_in        = cons_in;
    m_qheating_rates = qheating_rates;
    m_rad_fluxes     = rad_fluxes;
    m_z_phys         = z_phys;
    m_lat            = lat;
    m_lon            = lon;
 
    // Update the day and month
    time_t timestamp = time_t(time);
    struct tm *timeinfo = gmtime(&timestamp);
    if (m_fixed_orbital_year) {
        m_orbital_mon  = timeinfo->tm_mon + 1;
        m_orbital_day  = timeinfo->tm_mday;
        m_orbital_sec  = timeinfo->tm_hour*3600 + timeinfo->tm_min*60 + timeinfo->tm_sec;
    } else {
        m_orbital_year = timeinfo->tm_year + 1900;
        m_orbital_mon  = timeinfo->tm_mon  + 1;
        m_orbital_day  = timeinfo->tm_mday;
        m_orbital_sec  = timeinfo->tm_hour*3600 + timeinfo->tm_min*60 + timeinfo->tm_sec;
    }
 
    // Only allocate and proceed if we are going to update radiation
    m_update_rad = false;
    if (m_rad_freq_in_steps > 0) { m_update_rad = ( (m_step == 0) || (m_step % m_rad_freq_in_steps == 0) || updated_lsm); }
 
    if (m_update_rad) {
        // Call to Init() has set the dimensions: ncol & nlay
 
        // Allocate the buffer arrays
        alloc_buffers();
 
        // Fill the KOKKOS Views from AMReX MFs
        mf_to_kokkos_buffers(lmask, t_surf, lsm_input_ptrs);
 
        // (Re)define the datalog MF whenever the grids change; this must always
        // match the layout of cons_in since populateDatalogMF() iterates over it
        // while indexing m_col_offsets and m_qheating_rates.
        if (datalog_int > 0) {
            bool needs_define = ( (datalog_mf.boxArray()        != cons_in->boxArray()) ||
                                  (datalog_mf.DistributionMap() != cons_in->DistributionMap()) );
            if (needs_define) {
                datalog_mf.define(cons_in->boxArray(), cons_in->DistributionMap(), 25, 0);
                datalog_mf.setVal(0.0);
            }
        }
    }
}

Referenced by Run().

Here is the caller graph for this function:

write_rrtmgp_fluxes()

void Radiation::write_rrtmgp_fluxes

(

)

{
    Table2D<Real,Order::C> sw_flux_up_tab(sw_flux_up.data(), {0,0}, {static_cast<int>(sw_flux_up.extent(0)),static_cast<int>(sw_flux_up.extent(1))});
    Table2D<Real,Order::C> sw_flux_dn_tab(sw_flux_dn.data(), {0,0}, {static_cast<int>(sw_flux_dn.extent(0)),static_cast<int>(sw_flux_dn.extent(1))});
    Table2D<Real,Order::C> sw_flux_dn_dir_tab(sw_flux_dn_dir.data(), {0,0}, {static_cast<int>(sw_flux_dn_dir.extent(0)),static_cast<int>(sw_flux_dn_dir.extent(1))});
    Table2D<Real,Order::C> lw_flux_up_tab(lw_flux_up.data(), {0,0}, {static_cast<int>(lw_flux_up.extent(0)),static_cast<int>(lw_flux_up.extent(1))});
    Table2D<Real,Order::C> lw_flux_dn_tab(lw_flux_dn.data(), {0,0}, {static_cast<int>(lw_flux_dn.extent(0)),static_cast<int>(lw_flux_dn.extent(1))});
 
    int n_fluxes = 5;
    MultiFab mf_flux(m_cons_in->boxArray(), m_cons_in->DistributionMap(), n_fluxes, 0);
 
   for (MFIter mfi(mf_flux); mfi.isValid(); ++mfi) {
        const auto& vbx      = mfi.validbox();
        const int nx         = vbx.length(0);
        const int imin       = vbx.smallEnd(0);
        const int jmin       = vbx.smallEnd(1);
        const int offset     = m_col_offsets[mfi.index()];
        const Array4<Real>& dst_arr = mf_flux.array(mfi);
        ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
        {
            // map [i,j,k] 0-based to [icol, ilay] 0-based
            const int icol = (j-jmin)*nx + (i-imin) + offset;
            const int ilay = k;
 
            // SW and LW fluxes
            dst_arr(i,j,k,0) = sw_flux_up_tab(icol,ilay);
            dst_arr(i,j,k,1) = sw_flux_dn_tab(icol,ilay);
            dst_arr(i,j,k,2) = sw_flux_dn_dir_tab(icol,ilay);
            dst_arr(i,j,k,3) = lw_flux_up_tab(icol,ilay);
            dst_arr(i,j,k,4) = lw_flux_dn_tab(icol,ilay);
        });
   }
 
 
   std::string plotfilename = amrex::Concatenate("plt_rad", m_step, 5);
   Vector<std::string> flux_names = {"sw_flux_up", "sw_flux_dn", "sw_flux_dir",
                                     "lw_flux_up", "lw_flux_dn"};
   WriteSingleLevelPlotfile(plotfilename, mf_flux, flux_names, m_geom, static_cast<Real>(m_time), m_step);
}

Here is the call graph for this function:

WriteDataLog()

void Radiation::WriteDataLog ( const double &  time )

overridevirtual

Implements IRadiation.

{
    constexpr int datwidth = 14;
    constexpr int datprecision = 9;
    constexpr int timeprecision = 13;
 
    Gpu::HostVector<Real> h_avg_radqrsw, h_avg_radqrlw, h_avg_sw_up, h_avg_sw_dn, h_avg_sw_dn_dir, h_avg_lw_up, h_avg_lw_dn, h_avg_zenith;
    // Clear sky
    Gpu::HostVector<Real> h_avg_radqrcsw, h_avg_radqrclw, h_avg_sw_clr_up, h_avg_sw_clr_dn, h_avg_sw_clr_dn_dir, h_avg_lw_clr_up, h_avg_lw_clr_dn;
    // Clean sky
    Gpu::HostVector<Real> h_avg_sw_cln_up, h_avg_sw_cln_dn, h_avg_sw_cln_dn_dir, h_avg_lw_cln_up, h_avg_lw_cln_dn;
    // Clean clear sky
    Gpu::HostVector<Real> h_avg_sw_clnclr_up, h_avg_sw_clnclr_dn, h_avg_sw_clnclr_dn_dir, h_avg_lw_clnclr_up, h_avg_lw_clnclr_dn;
 
 
    auto domain = m_geom.Domain();
    h_avg_radqrsw    = sumToLine(datalog_mf, 0, 1, domain, 2);
    h_avg_radqrlw    = sumToLine(datalog_mf, 1, 1, domain, 2);
    h_avg_sw_up      = sumToLine(datalog_mf, 2, 1, domain, 2);
    h_avg_sw_dn      = sumToLine(datalog_mf, 3, 1, domain, 2);
    h_avg_sw_dn_dir  = sumToLine(datalog_mf, 4, 1, domain, 2);
    h_avg_lw_up      = sumToLine(datalog_mf, 5, 1, domain, 2);
    h_avg_lw_dn      = sumToLine(datalog_mf, 6, 1, domain, 2);
    h_avg_zenith     = sumToLine(datalog_mf, 7, 1, domain, 2);
 
    h_avg_radqrcsw       = sumToLine(datalog_mf, 8, 1, domain, 2);
    h_avg_radqrclw       = sumToLine(datalog_mf, 9, 1, domain, 2);
    h_avg_sw_clr_up      = sumToLine(datalog_mf, 10, 1, domain, 2);
    h_avg_sw_clr_dn      = sumToLine(datalog_mf, 11, 1, domain, 2);
    h_avg_sw_clr_dn_dir  = sumToLine(datalog_mf, 12, 1, domain, 2);
    h_avg_lw_clr_up      = sumToLine(datalog_mf, 13, 1, domain, 2);
    h_avg_lw_clr_dn      = sumToLine(datalog_mf, 14, 1, domain, 2);
 
    if (m_extra_clnsky_diag) {
        h_avg_sw_cln_up      = sumToLine(datalog_mf, 15, 1, domain, 2);
        h_avg_sw_cln_dn      = sumToLine(datalog_mf, 16, 1, domain, 2);
        h_avg_sw_cln_dn_dir  = sumToLine(datalog_mf, 17, 1, domain, 2);
        h_avg_lw_cln_up      = sumToLine(datalog_mf, 18, 1, domain, 2);
        h_avg_lw_cln_dn      = sumToLine(datalog_mf, 19, 1, domain, 2);
    }
 
    if (m_extra_clnclrsky_diag) {
        h_avg_sw_clnclr_up      = sumToLine(datalog_mf, 20, 1, domain, 2);
        h_avg_sw_clnclr_dn      = sumToLine(datalog_mf, 21, 1, domain, 2);
        h_avg_sw_clnclr_dn_dir  = sumToLine(datalog_mf, 22, 1, domain, 2);
        h_avg_lw_clnclr_up      = sumToLine(datalog_mf, 23, 1, domain, 2);
        h_avg_lw_clnclr_dn      = sumToLine(datalog_mf, 24, 1, domain, 2);
    }
 
    Real area_z = static_cast<Real>(domain.length(0)*domain.length(1));
    int nz = domain.length(2);
    for (int k = 0; k < nz; k++) {
        h_avg_radqrsw[k] /= area_z;
        h_avg_radqrlw[k] /= area_z;
        h_avg_sw_up[k] /= area_z;
        h_avg_sw_dn[k] /= area_z;
        h_avg_sw_dn_dir[k] /= area_z;
        h_avg_lw_up[k] /= area_z;
        h_avg_lw_dn[k] /= area_z;
        h_avg_zenith[k] /= area_z;
 
        h_avg_radqrcsw[k] /= area_z;
        h_avg_radqrclw[k] /= area_z;
        h_avg_sw_clr_up[k] /= area_z;
        h_avg_sw_clr_dn[k] /= area_z;
        h_avg_sw_clr_dn_dir[k] /= area_z;
        h_avg_lw_clr_up[k] /= area_z;
        h_avg_lw_clr_dn[k] /= area_z;
    }
 
    if (m_extra_clnsky_diag) {
        for (int k = 0; k < nz; k++) {
            h_avg_sw_cln_up[k] /= area_z;
            h_avg_sw_cln_dn[k] /= area_z;
            h_avg_sw_cln_dn_dir[k] /= area_z;
            h_avg_lw_cln_up[k] /= area_z;
            h_avg_lw_cln_dn[k] /= area_z;
        }
    }
 
    if (m_extra_clnclrsky_diag) {
        for (int k = 0; k < nz; k++) {
            h_avg_sw_clnclr_up[k] /= area_z;
            h_avg_sw_clnclr_dn[k] /= area_z;
            h_avg_sw_clnclr_dn_dir[k] /= area_z;
            h_avg_lw_clnclr_up[k] /= area_z;
            h_avg_lw_clnclr_dn[k] /= area_z;
        }
    }
 
    if (ParallelDescriptor::IOProcessor()) {
        std::ostream& log = *datalog;
        if (log.good()) {
 
            for (int k = 0; k < nz; k++)
            {
                Real z = k * m_geom.CellSize(2);
                log << std::setw(datwidth) << std::setprecision(timeprecision) << time << " "
                    << std::setw(datwidth) << std::setprecision(datprecision) << z << " "
                    << h_avg_radqrsw[k]   << " " << h_avg_radqrlw[k]       << " " << h_avg_sw_up[k] << " "
                    << h_avg_sw_dn[k]     << " " << h_avg_sw_dn_dir[k]     << " " << h_avg_lw_up[k] << " "
                    << h_avg_lw_dn[k]     << " " << h_avg_zenith[k]        << " "
                    << h_avg_radqrcsw[k]  << " " << h_avg_radqrclw[k]      << " " << h_avg_sw_clr_up[k] << " "
                    << h_avg_sw_clr_dn[k] << " " << h_avg_sw_clr_dn_dir[k] << " " << h_avg_lw_clr_up[k] << " "
                    << h_avg_lw_clr_dn[k] << " ";
                    if (m_extra_clnsky_diag) {
                        log << h_avg_sw_cln_up[k] << " " << h_avg_sw_cln_dn[k] << " " << h_avg_sw_cln_dn_dir[k] << " "
                            << h_avg_lw_cln_up[k] << " " << h_avg_lw_cln_dn[k] << " ";
                    } else {
                        log << zero << " " << zero << " " << zero << " " << zero << " " << zero << " ";
                    }
 
                    if (m_extra_clnclrsky_diag) {
                        log << h_avg_sw_clnclr_up[k] << " " << h_avg_sw_clnclr_dn[k] << " " << h_avg_sw_clnclr_dn_dir[k] << " "
                            << h_avg_lw_clnclr_up[k] << " " << h_avg_lw_clnclr_dn[k] << std::endl;
                    } else {
                        log << zero << " " << zero << " " << zero << " " << zero << " " << zero << std::endl;
                    }
            }
            // Write top face values
            Real z = nz * m_geom.CellSize(2);
            log << std::setw(datwidth) << std::setprecision(timeprecision) << time << " "
                << std::setw(datwidth) << std::setprecision(datprecision) << z << " "
                << zero << " " << zero << " " << zero << " " << zero << " " << zero << " " << zero << " "
                << zero << " " << zero << " "
                << zero << " " << zero << " " << zero << " " << zero << " " << zero << " " << zero << " "
                << zero << " "
                << zero << " " << zero << " " << zero << " " << zero << " " << zero << " "
                << zero << " " << zero << " " << zero << " " << zero << " " << zero
                << std::endl;
        }
    }
}

Member Data Documentation

aero_g_sw

real3d_k Radiation::aero_g_sw

private

aero_ssa_sw

real3d_k Radiation::aero_ssa_sw

private

aero_tau_lw

real3d_k Radiation::aero_tau_lw

private

aero_tau_sw

real3d_k Radiation::aero_tau_sw

private

cldfrac_tot

real2d_k Radiation::cldfrac_tot

private

d_tint

real2d_k Radiation::d_tint

private

datalog_mf

amrex::MultiFab Radiation::datalog_mf

private

eff_radius_qc

real2d_k Radiation::eff_radius_qc

private

eff_radius_qi

real2d_k Radiation::eff_radius_qi

private

gas_names_offset

std::vector<std::string> Radiation::gas_names_offset

private

iwp

real2d_k Radiation::iwp

private

lat

real1d_k Radiation::lat

private

lon

real1d_k Radiation::lon

private

lw_bnd_flux_dn

real3d_k Radiation::lw_bnd_flux_dn

private

lw_bnd_flux_up

real3d_k Radiation::lw_bnd_flux_up

private

lw_clnclrsky_flux_dn

real2d_k Radiation::lw_clnclrsky_flux_dn

private

lw_clnclrsky_flux_up

real2d_k Radiation::lw_clnclrsky_flux_up

private

lw_clnsky_flux_dn

real2d_k Radiation::lw_clnsky_flux_dn

private

lw_clnsky_flux_up

real2d_k Radiation::lw_clnsky_flux_up

private

lw_clrsky_flux_dn

real2d_k Radiation::lw_clrsky_flux_dn

private

lw_clrsky_flux_up

real2d_k Radiation::lw_clrsky_flux_up

private

lw_clrsky_heating

real2d_k Radiation::lw_clrsky_heating

private

lw_flux_dn

real2d_k Radiation::lw_flux_dn

private

lw_flux_up

real2d_k Radiation::lw_flux_up

private

lw_heating

real2d_k Radiation::lw_heating

private

lw_src

real1d_k Radiation::lw_src

private

lwp

real2d_k Radiation::lwp

private

m_ba

amrex::BoxArray Radiation::m_ba

private

m_ch4vmr

amrex::Real Radiation::m_ch4vmr = amrex::Real(1807.851e-9)

private

m_co2vmr

amrex::Real Radiation::m_co2vmr = amrex::Real(388.717e-6)

private

m_col_offsets

amrex::Vector<int> Radiation::m_col_offsets

private

Referenced by Init().

m_cons_in

amrex::MultiFab* Radiation::m_cons_in = nullptr

private

m_covmr

amrex::Real Radiation::m_covmr = amrex::Real(1.0e-7)

private

m_do_aerosol_rad

bool Radiation::m_do_aerosol_rad = false

private

m_do_subcol_sampling

bool Radiation::m_do_subcol_sampling = true

private

m_dt

double Radiation::m_dt

private

m_extra_clnclrsky_diag

bool Radiation::m_extra_clnclrsky_diag = false

private

m_extra_clnsky_diag

bool Radiation::m_extra_clnsky_diag = false

private

m_fixed_orbital_year

bool Radiation::m_fixed_orbital_year = false

private

m_fixed_solar_zenith_angle

amrex::Real Radiation::m_fixed_solar_zenith_angle = -amrex::Real(9999.)

private

m_fixed_total_solar_irradiance

amrex::Real Radiation::m_fixed_total_solar_irradiance = -amrex::Real(9999.)

private

m_gas_concs

GasConcsK<amrex::Real, layout_t, KokkosDefaultDevice> Radiation::m_gas_concs

private

m_gas_mol_weights

real1d_k Radiation::m_gas_mol_weights

private

m_gas_names

const std::vector<std::string> Radiation::m_gas_names

private

Initial value:

= {"H2O", "CO2", "O3", "N2O",
                                                  "CO" , "CH4", "O2", "N2" }

m_geom

amrex::Geometry Radiation::m_geom

private

m_ice

bool Radiation::m_ice = false

private

m_lat

amrex::MultiFab* Radiation::m_lat = nullptr

private

m_lat_cons

amrex::Real Radiation::m_lat_cons = amrex::Real(39.809860)

private

m_lev

int Radiation::m_lev

private

Referenced by rad_run_impl().

m_lon

amrex::MultiFab* Radiation::m_lon = nullptr

private

m_lon_cons

amrex::Real Radiation::m_lon_cons = -amrex::Real(98.555183)

private

m_lsm

bool Radiation::m_lsm = false

private

m_lsm_input_names

amrex::Vector<std::string> Radiation::m_lsm_input_names

private

Initial value:

= {"t_sfc"          , "sfc_emis"       ,
                                                    "sfc_alb_dir_vis", "sfc_alb_dir_nir",
                                                    "sfc_alb_dif_vis", "sfc_alb_dif_nir"}

Referenced by get_lsm_input_varnames().

m_lsm_output_names

amrex::Vector<std::string> Radiation::m_lsm_output_names

private

Initial value:

= {"cos_zenith_angle"  , "sw_flux_dn"        ,
                                                     "sw_flux_dn_dir_vis", "sw_flux_dn_dir_nir",
                                                     "sw_flux_dn_dif_vis", "sw_flux_dn_dif_nir",
                                                     "lw_flux_dn"}

Referenced by get_lsm_output_varnames().

m_moist

bool Radiation::m_moist = false

private

m_mol_weight_gas

const std::vector<amrex::Real> Radiation::m_mol_weight_gas

private

Initial value:

= {amrex::Real(18.01528), amrex::Real(44.00950), amrex::Real(47.9982), amrex::Real(44.0128),
                                                       amrex::Real(28.01010), amrex::Real(16.04246), amrex::Real(31.9980), amrex::Real(28.0134)}

m_n2ovmr

amrex::Real Radiation::m_n2ovmr = amrex::Real(323.141e-9)

private

m_n2vmr

amrex::Real Radiation::m_n2vmr = amrex::Real(0.7906)

private

m_ncol

int Radiation::m_ncol

private

Referenced by Init().

m_ncol_chunk

int Radiation::m_ncol_chunk = 1024

private

Referenced by Init().

m_ncol_chunk_requested

int Radiation::m_ncol_chunk_requested = 1024

private

Referenced by Init().

m_ngas

int Radiation::m_ngas = 8

private

m_nlay

int Radiation::m_nlay

private

Referenced by Init().

m_nlwbands

int Radiation::m_nlwbands

private

m_nlwgpts

int Radiation::m_nlwgpts

private

m_nswbands

int Radiation::m_nswbands

private

m_nswgpts

int Radiation::m_nswgpts

private

m_o2vmr

amrex::Real Radiation::m_o2vmr = amrex::Real(0.209448)

private

m_o3_size

int Radiation::m_o3_size

private

m_o3vmr

amrex::Vector<amrex::Real> Radiation::m_o3vmr

private

m_orbital_day

int Radiation::m_orbital_day = -9999

private

Referenced by rad_run_impl().

m_orbital_eccen

amrex::Real Radiation::m_orbital_eccen = -amrex::Real(9999.)

private

m_orbital_mon

int Radiation::m_orbital_mon = -9999

private

Referenced by rad_run_impl().

m_orbital_mvelp

amrex::Real Radiation::m_orbital_mvelp = -amrex::Real(9999.)

private

m_orbital_obliq

amrex::Real Radiation::m_orbital_obliq = -amrex::Real(9999.)

private

m_orbital_sec

int Radiation::m_orbital_sec = -9999

private

Referenced by rad_run_impl().

m_orbital_year

int Radiation::m_orbital_year = -9999

private

Referenced by rad_run_impl().

m_qheating_rates

amrex::MultiFab* Radiation::m_qheating_rates = nullptr

private

m_qi_comp

int Radiation::m_qi_comp = -1

private

m_rad_fluxes

amrex::MultiFab* Radiation::m_rad_fluxes = nullptr

private

m_rad_freq_in_steps

int Radiation::m_rad_freq_in_steps = 1

private

m_rad_nvar

int Radiation::m_rad_nvar = 12

private

m_rad_t_sfc

amrex::Real Radiation::m_rad_t_sfc = -1

private

m_rad_write_fluxes

bool Radiation::m_rad_write_fluxes = false

private

m_step

int Radiation::m_step

private

m_time

double Radiation::m_time

private

m_update_rad

bool Radiation::m_update_rad = false

private

Referenced by rad_run_impl().

m_z_phys

amrex::MultiFab* Radiation::m_z_phys = nullptr

private

mu0

real1d_k Radiation::mu0

private

o3_lay

real1d_k Radiation::o3_lay

private

p_lay

real2d_k Radiation::p_lay

private

p_lev

real2d_k Radiation::p_lev

private

qc_lay

real2d_k Radiation::qc_lay

private

qi_lay

real2d_k Radiation::qi_lay

private

qv_lay

real2d_k Radiation::qv_lay

private

r_lay

real2d_k Radiation::r_lay

private

rrtmgp_cloud_optics_file_lw

std::string Radiation::rrtmgp_cloud_optics_file_lw

private

rrtmgp_cloud_optics_file_sw

std::string Radiation::rrtmgp_cloud_optics_file_sw

private

rrtmgp_cloud_optics_lw

std::string Radiation::rrtmgp_cloud_optics_lw = "rrtmgp-cloud-optics-coeffs-lw.nc"

private

rrtmgp_cloud_optics_sw

std::string Radiation::rrtmgp_cloud_optics_sw = "rrtmgp-cloud-optics-coeffs-sw.nc"

private

rrtmgp_coeffs_file_lw

std::string Radiation::rrtmgp_coeffs_file_lw

private

rrtmgp_coeffs_file_sw

std::string Radiation::rrtmgp_coeffs_file_sw

private

rrtmgp_coeffs_lw

std::string Radiation::rrtmgp_coeffs_lw = "rrtmgp-data-lw-g256-2018-12-04.nc"

private

rrtmgp_coeffs_sw

std::string Radiation::rrtmgp_coeffs_sw = "rrtmgp-data-sw-g224-2018-12-04.nc"

private

rrtmgp_file_path

std::string Radiation::rrtmgp_file_path = "."

private

sfc_alb_dif

real2d_k Radiation::sfc_alb_dif

private

sfc_alb_dif_nir

real1d_k Radiation::sfc_alb_dif_nir

private

sfc_alb_dif_vis

real1d_k Radiation::sfc_alb_dif_vis

private

sfc_alb_dir

real2d_k Radiation::sfc_alb_dir

private

sfc_alb_dir_nir

real1d_k Radiation::sfc_alb_dir_nir

private

sfc_alb_dir_vis

real1d_k Radiation::sfc_alb_dir_vis

private

sfc_emis

real1d_k Radiation::sfc_emis

private

sfc_flux_dif_nir

real1d_k Radiation::sfc_flux_dif_nir

private

sfc_flux_dif_vis

real1d_k Radiation::sfc_flux_dif_vis

private

sfc_flux_dir_nir

real1d_k Radiation::sfc_flux_dir_nir

private

sfc_flux_dir_vis

real1d_k Radiation::sfc_flux_dir_vis

private

sw_bnd_flux_dif

real3d_k Radiation::sw_bnd_flux_dif

private

sw_bnd_flux_dir

real3d_k Radiation::sw_bnd_flux_dir

private

sw_bnd_flux_dn

real3d_k Radiation::sw_bnd_flux_dn

private

sw_bnd_flux_up

real3d_k Radiation::sw_bnd_flux_up

private

sw_clnclrsky_flux_dn

real2d_k Radiation::sw_clnclrsky_flux_dn

private

sw_clnclrsky_flux_dn_dir

real2d_k Radiation::sw_clnclrsky_flux_dn_dir

private

sw_clnclrsky_flux_up

real2d_k Radiation::sw_clnclrsky_flux_up

private

sw_clnsky_flux_dn

real2d_k Radiation::sw_clnsky_flux_dn

private

sw_clnsky_flux_dn_dir

real2d_k Radiation::sw_clnsky_flux_dn_dir

private

sw_clnsky_flux_up

real2d_k Radiation::sw_clnsky_flux_up

private

sw_clrsky_flux_dn

real2d_k Radiation::sw_clrsky_flux_dn

private

sw_clrsky_flux_dn_dir

real2d_k Radiation::sw_clrsky_flux_dn_dir

private

sw_clrsky_flux_up

real2d_k Radiation::sw_clrsky_flux_up

private

sw_clrsky_heating

real2d_k Radiation::sw_clrsky_heating

private

sw_flux_dn

real2d_k Radiation::sw_flux_dn

private

sw_flux_dn_dir

real2d_k Radiation::sw_flux_dn_dir

private

sw_flux_up

real2d_k Radiation::sw_flux_up

private

sw_heating

real2d_k Radiation::sw_heating

private

t_lay

real2d_k Radiation::t_lay

private

t_lev

real2d_k Radiation::t_lev

private

t_sfc

real1d_k Radiation::t_sfc

private

z_del

real2d_k Radiation::z_del

private

---

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

• Source/PhysicsInterfaces/Radiation/ERF_Radiation.H
• Source/PhysicsInterfaces/Radiation/ERF_Radiation.cpp