ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
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, amrex::Real &time, const amrex::Real &dt, const amrex::BoxArray &ba, amrex::Geometry &geom, amrex::MultiFab *cons_in, amrex::MultiFab *lsm_fluxes, amrex::MultiFab *lsm_zenith, amrex::Vector< amrex::MultiFab * > &lsm_input_ptrs, amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs, amrex::MultiFab *qheating_rates, amrex::MultiFab *z_phys, amrex::MultiFab *lat, amrex::MultiFab *lon) override
 
void set_grids (int &level, int &step, amrex::Real &time, const amrex::Real &dt, const amrex::BoxArray &ba, amrex::Geometry &geom, amrex::MultiFab *cons_in, amrex::MultiFab *lsm_fluxes, amrex::MultiFab *lsm_zenith, amrex::Vector< amrex::MultiFab * > &lsm_input_ptrs, amrex::MultiFab *qheating_rates, amrex::MultiFab *z_phys, amrex::MultiFab *lat, amrex::MultiFab *lon)
 
template<typename LandSurfaceModelType >
void set_lsm_inputs (LandSurfaceModelType *model)
 
void alloc_buffers ()
 
void dealloc_buffers ()
 
void mf_to_kokkos_buffers (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 amrex::Real &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
 
amrex::Real m_time
 
amrex::Real m_dt
 
amrex::Geometry m_geom
 
amrex::BoxArray m_ba
 
bool m_update_rad = false
 
bool m_rad_write_fluxes = false
 
bool m_first_step = true
 
bool m_moist = false
 
bool m_ice = false
 
bool m_lsm = false
 
amrex::Vector< std::string > m_lsm_input_names
 
amrex::Vector< std::string > m_lsm_output_names = {"sw_flux_dn", "lw_flux_dn"}
 
amrex::Real m_rad_t_sfc = -1
 
amrex::MultiFab * m_cons_in = nullptr
 
amrex::MultiFab * m_qheating_rates = nullptr
 
amrex::MultiFab * m_z_phys = nullptr
 
amrex::MultiFab * m_lat = nullptr
 
amrex::MultiFab * m_lon = nullptr
 
amrex::Real m_lat_cons = 39.809860
 
amrex::Real m_lon_cons = -98.555183
 
amrex::MultiFab * m_lsm_fluxes = nullptr
 
amrex::MultiFab * m_lsm_zenith = nullptr
 
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-g224-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 = 388.717e-6
 
amrex::Vector< amrex::Real > m_o3vmr
 
amrex::Real m_n2ovmr = 323.141e-9
 
amrex::Real m_covmr = 1.0e-7
 
amrex::Real m_ch4vmr = 1807.851e-9
 
amrex::Real m_o2vmr = 0.209448
 
amrex::Real m_n2vmr = 0.7906
 
int m_o3_size
 
real1d_k m_gas_mol_weights
 
std::vector< std::string > gas_names_offset
 
GasConcsK< amrex::Real, layout_t, KokkosDefaultDevicem_gas_concs
 
int m_ncol
 
int m_nlay
 
amrex::Vector< int > m_col_offsets
 
bool m_do_aerosol_rad = true
 
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 = -9999.
 
amrex::Real m_orbital_obliq = -9999.
 
amrex::Real m_orbital_mvelp = -9999.
 
amrex::Real m_fixed_total_solar_irradiance = -9999.
 
amrex::Real m_fixed_solar_zenith_angle = -9999.
 
int m_nswbands = 14
 
int m_nlwbands = 16
 
int m_nswgpts = 112
 
int m_nlwgpts = 128
 
int m_rad_freq_in_steps = 1
 
bool m_do_subcol_sampling = true
 
real1d_k o3_lay
 
real1d_k cosine_zenith
 
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 p_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 tmp2d
 
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
 
real2d_k emis_sfc
 
real3d_k aero_tau_sw
 
real3d_k aero_ssa_sw
 
real3d_k aero_g_sw
 
real3d_k aero_tau_lw
 
real3d_k cld_tau_sw_bnd
 
real3d_k cld_tau_lw_bnd
 
real3d_k cld_tau_sw_gpt
 
real3d_k cld_tau_lw_gpt
 

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 
)
21 {
22  // Note that Kokkos is now initialized in main.cpp
23 
24  // Check if we have a valid moisture model
25  if (sc.moisture_type != MoistureType::None) { m_moist = true; }
26 
27  // Check if we have a moisture model with ice
28  if (sc.moisture_type == MoistureType::SAM) { m_ice = true; }
29 
30  // Check if we have a land surface model enabled
31  if (sc.lsm_type != LandSurfaceType::None) { m_lsm = true; }
32 
33  // Construct parser object for following reads
34  ParmParse pp("erf");
35 
36  // Must specify a surface temp without a LSM
37  if (!m_lsm) { pp.get("rad_t_sfc",m_rad_t_sfc); }
38 
39  // Radiation timestep, as a number of atm steps
40  pp.query("rad_freq_in_steps", m_rad_freq_in_steps);
41 
42  // Flag to write fluxes to plt file
43  pp.query("rad_write_fluxes", m_rad_write_fluxes);
44 
45  // Do MCICA subcolumn sampling
46  pp.query("rad_do_subcol_sampling", m_do_subcol_sampling);
47 
48  // Determine orbital year. If orbital_year is negative, use current year
49  // from timestamp for orbital year; if positive, use provided orbital year
50  // for duration of simulation.
51  m_fixed_orbital_year = pp.query("rad_orbital_year", m_orbital_year);
52 
53  // Get orbital parameters from inputs file
54  pp.query("rad_orbital_eccentricity", m_orbital_eccen);
55  pp.query("rad_orbital_obliquity" , m_orbital_obliq);
56  pp.query("rad_orbital_mvelp" , m_orbital_mvelp);
57 
58  // Get a constant lat/lon for idealized simulations
59  pp.query("rad_cons_lat", m_lat_cons);
60  pp.query("rad_cons_lon", m_lon_cons);
61 
62  // Value for prescribing an invariant solar constant (i.e. total solar irradiance at
63  // TOA). Used for idealized experiments such as RCE. Disabled when value is less than 0.
64  pp.query("fixed_total_solar_irradiance", m_fixed_total_solar_irradiance);
65 
66  // Determine whether or not we are using a fixed solar zenith angle (positive value)
67  pp.query("fixed_solar_zenith_angle", m_fixed_solar_zenith_angle);
68 
69  // Get prescribed surface values of greenhouse gases
70  pp.query("co2vmr", m_co2vmr);
71  pp.queryarr("o3vmr" , m_o3vmr );
72  pp.query("n2ovmr", m_n2ovmr);
73  pp.query("covmr" , m_covmr );
74  pp.query("ch4vmr", m_ch4vmr);
75  pp.query("o2vmr" , m_o2vmr );
76  pp.query("n2vmr" , m_n2vmr );
77 
78  // Required aerosol optical properties from SPA
79  pp.query("rad_do_aerosol", m_do_aerosol_rad);
80 
81  // Whether we do extra clean/clear sky calculations
82  pp.query("rad_extra_clnclrsky_diag", m_extra_clnclrsky_diag);
83  pp.query("rad_extra_clnsky_diag" , m_extra_clnsky_diag);
84 
85  // Parse the band and gauss pt sizes
86  pp.query("nswbands", m_nswbands);
87  pp.query("nlwbands", m_nlwbands);
88  pp.query("nswgpts" , m_nswgpts );
89  pp.query("nlwgpts" , m_nlwgpts );
90 
91  // Parse path and file names
92  pp.query("rrtmgp_file_path" , rrtmgp_file_path);
93  pp.query("rrtmgp_coeffs_sw" , rrtmgp_coeffs_sw );
94  pp.query("rrtmgp_coeffs_lw" , rrtmgp_coeffs_lw );
95  pp.query("rrtmgp_cloud_optics_sw", rrtmgp_cloud_optics_sw);
96  pp.query("rrtmgp_cloud_optics_lw", rrtmgp_cloud_optics_lw);
97 
98  // Append file names to path
103 
104  // Output for user
105  if (lev == 0) {
106  Print() << "Radiation interface constructed:\n";
107  Print() << "========================================================\n";
108  Print() << "Coeff SW file: " << rrtmgp_coeffs_file_sw << "\n";
109  Print() << "Coeff LW file: " << rrtmgp_coeffs_file_lw << "\n";
110  Print() << "Cloud SW file: " << rrtmgp_cloud_optics_file_sw << "\n";
111  Print() << "Cloud LW file: " << rrtmgp_cloud_optics_file_lw << "\n";
112  Print() << "========================================================\n";
113  }
114 }
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real pp(amrex::Real y)
Definition: ERF_MicrophysicsUtils.H:230
std::string rrtmgp_coeffs_file_sw
Definition: ERF_Radiation.H:342
std::string rrtmgp_coeffs_sw
Definition: ERF_Radiation.H:338
int m_rad_freq_in_steps
Definition: ERF_Radiation.H:421
amrex::Real m_lon_cons
Definition: ERF_Radiation.H:327
int m_nswbands
Definition: ERF_Radiation.H:415
bool m_do_aerosol_rad
Definition: ERF_Radiation.H:382
std::string rrtmgp_coeffs_lw
Definition: ERF_Radiation.H:339
bool m_moist
Definition: ERF_Radiation.H:295
std::string rrtmgp_cloud_optics_file_lw
Definition: ERF_Radiation.H:345
bool m_lsm
Definition: ERF_Radiation.H:299
bool m_do_subcol_sampling
Definition: ERF_Radiation.H:424
bool m_rad_write_fluxes
Definition: ERF_Radiation.H:290
std::string rrtmgp_cloud_optics_file_sw
Definition: ERF_Radiation.H:344
amrex::Real m_orbital_mvelp
Definition: ERF_Radiation.H:402
amrex::Vector< amrex::Real > m_o3vmr
Definition: ERF_Radiation.H:356
std::string rrtmgp_cloud_optics_sw
Definition: ERF_Radiation.H:340
int m_nlwgpts
Definition: ERF_Radiation.H:418
std::string rrtmgp_cloud_optics_lw
Definition: ERF_Radiation.H:341
amrex::Real m_co2vmr
Definition: ERF_Radiation.H:355
bool m_extra_clnsky_diag
Definition: ERF_Radiation.H:385
amrex::Real m_lat_cons
Definition: ERF_Radiation.H:326
amrex::Real m_fixed_total_solar_irradiance
Definition: ERF_Radiation.H:407
amrex::Real m_o2vmr
Definition: ERF_Radiation.H:360
amrex::Real m_n2ovmr
Definition: ERF_Radiation.H:357
amrex::Real m_rad_t_sfc
Definition: ERF_Radiation.H:310
std::string rrtmgp_file_path
Definition: ERF_Radiation.H:337
int m_nlwbands
Definition: ERF_Radiation.H:416
amrex::Real m_orbital_eccen
Definition: ERF_Radiation.H:400
amrex::Real m_covmr
Definition: ERF_Radiation.H:358
std::string rrtmgp_coeffs_file_lw
Definition: ERF_Radiation.H:343
amrex::Real m_ch4vmr
Definition: ERF_Radiation.H:359
bool m_ice
Definition: ERF_Radiation.H:296
amrex::Real m_fixed_solar_zenith_angle
Definition: ERF_Radiation.H:411
amrex::Real m_orbital_obliq
Definition: ERF_Radiation.H:401
amrex::Real m_n2vmr
Definition: ERF_Radiation.H:361
bool m_fixed_orbital_year
Definition: ERF_Radiation.H:399
int m_nswgpts
Definition: ERF_Radiation.H:417
bool m_extra_clnclrsky_diag
Definition: ERF_Radiation.H:386
int m_orbital_year
Definition: ERF_Radiation.H:391
LandSurfaceType lsm_type
Definition: ERF_DataStruct.H:866
MoistureType moisture_type
Definition: ERF_DataStruct.H:863
Here is the call graph for this function:

◆ ~Radiation()

Radiation::~Radiation ( )
inline
56  {
57  // Note that Kokkos is now finalized in main.cpp
58  }

Member Function Documentation

◆ alloc_buffers()

void Radiation::alloc_buffers ( )
185 {
186  // 1d size (m_ngas)
187  const Real* mol_weight_gas_p = m_mol_weight_gas.data();
188  const std::string* gas_names_p = m_gas_names.data();
189  m_gas_mol_weights = real1d_k("m_gas_mol_weights", m_ngas);
190  realHost1d_k m_gas_mol_weights_h("m_gas_mol_weights_h", m_ngas);
191  gas_names_offset.clear(); gas_names_offset.resize(m_ngas);
192  std::string* gas_names_offset_p = gas_names_offset.data();
193  Kokkos::parallel_for(Kokkos::RangePolicy<Kokkos::Serial>(0, m_ngas),
194  KOKKOS_LAMBDA (int igas)
195  {
196  m_gas_mol_weights_h(igas) = mol_weight_gas_p[igas];
197  gas_names_offset_p[igas] = gas_names_p[igas];
198  });
199  Kokkos::deep_copy(m_gas_mol_weights, m_gas_mol_weights_h);
200 
201  // 1d size (1 or nlay)
202  m_o3_size = m_o3vmr.size();
203  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(((m_o3_size==1) || (m_o3_size==m_nlay)),
204  "O3 VMR array must be length 1 or nlay");
205  Real* o3vmr_p = m_o3vmr.data();
206  o3_lay = real1d_k("o3_lay", m_o3_size);
207  realHost1d_k o3_lay_h("o3_lay_h", m_o3_size);
208  Kokkos::parallel_for(Kokkos::RangePolicy<Kokkos::Serial>(0, m_o3_size),
209  KOKKOS_LAMBDA (int io3)
210  {
211  o3_lay_h(io3) = o3vmr_p[io3];
212  });
213  Kokkos::deep_copy(o3_lay, o3_lay_h);
214 
215  // 1d size (ncol)
216  cosine_zenith = real1d_k("cosine_zenith" , m_ncol);
217  mu0 = real1d_k("mu0" , m_ncol);
218  sfc_alb_dir_vis = real1d_k("sfc_alb_dir_vis" , m_ncol);
219  sfc_alb_dir_nir = real1d_k("sfc_alb_dir_nir" , m_ncol);
220  sfc_alb_dif_vis = real1d_k("sfc_alb_dif_vis" , m_ncol);
221  sfc_alb_dif_nir = real1d_k("sfc_alb_dif_nir" , m_ncol);
222  sfc_flux_dir_vis = real1d_k("sfc_flux_dir_vis", m_ncol);
223  sfc_flux_dir_nir = real1d_k("sfc_flux_dir_nir", m_ncol);
224  sfc_flux_dif_vis = real1d_k("sfc_flux_dif_vis", m_ncol);
225  sfc_flux_dif_nir = real1d_k("sfc_flux_dif_nir", m_ncol);
226  lat = real1d_k("lat" , m_ncol);
227  lon = real1d_k("lon" , m_ncol);
228  sfc_emis = real1d_k("sfc_emis" , m_ncol);
229  t_sfc = real1d_k("t_sfc" , m_ncol);
230  lw_src = real1d_k("lw_src" , m_ncol);
231 
232  // 2d size (ncol, nlay)
233  r_lay = real2d_k("r_lay" , m_ncol, m_nlay);
234  p_lay = real2d_k("p_lay" , m_ncol, m_nlay);
235  t_lay = real2d_k("t_lay" , m_ncol, m_nlay);
236  z_del = real2d_k("z_del" , m_ncol, m_nlay);
237  p_del = real2d_k("p_del" , m_ncol, m_nlay);
238  qv_lay = real2d_k("qv" , m_ncol, m_nlay);
239  qc_lay = real2d_k("qc" , m_ncol, m_nlay);
240  qi_lay = real2d_k("qi" , m_ncol, m_nlay);
241  cldfrac_tot = real2d_k("cldfrac_tot" , m_ncol, m_nlay);
242  eff_radius_qc = real2d_k("eff_radius_qc", m_ncol, m_nlay);
243  eff_radius_qi = real2d_k("eff_radius_qi", m_ncol, m_nlay);
244  tmp2d = real2d_k("tmp2d" , m_ncol, m_nlay);
245  lwp = real2d_k("lwp" , m_ncol, m_nlay);
246  iwp = real2d_k("iwp" , m_ncol, m_nlay);
247  sw_heating = real2d_k("sw_heating" , m_ncol, m_nlay);
248  lw_heating = real2d_k("lw_heating" , m_ncol, m_nlay);
249  sw_clrsky_heating = real2d_k("sw_clrsky_heating", m_ncol, m_nlay);
250  lw_clrsky_heating = real2d_k("lw_clrsky_heating", m_ncol, m_nlay);
251 
252  // 2d size (ncol, nlay+1)
253  d_tint = real2d_k("d_tint" , m_ncol, m_nlay+1);
254  p_lev = real2d_k("p_lev" , m_ncol, m_nlay+1);
255  t_lev = real2d_k("t_lev" , m_ncol, m_nlay+1);
256  sw_flux_up = real2d_k("sw_flux_up" , m_ncol, m_nlay+1);
257  sw_flux_dn = real2d_k("sw_flux_dn" , m_ncol, m_nlay+1);
258  sw_flux_dn_dir = real2d_k("sw_flux_dn_dir" , m_ncol, m_nlay+1);
259  lw_flux_up = real2d_k("sw_flux_up" , m_ncol, m_nlay+1);
260  lw_flux_dn = real2d_k("sw_flux_dn" , m_ncol, m_nlay+1);
261  sw_clnclrsky_flux_up = real2d_k("sw_clnclrsky_flux_up" , m_ncol, m_nlay+1);
262  sw_clnclrsky_flux_dn = real2d_k("sw_clnclrsky_flux_dn" , m_ncol, m_nlay+1);
263  sw_clnclrsky_flux_dn_dir = real2d_k("sw_clnclrsky_flux_dn_dir", m_ncol, m_nlay+1);
264  sw_clrsky_flux_up = real2d_k("sw_clrsky_flux_up" , m_ncol, m_nlay+1);
265  sw_clrsky_flux_dn = real2d_k("sw_clrsky_flux_dn" , m_ncol, m_nlay+1);
266  sw_clrsky_flux_dn_dir = real2d_k("sw_clrsky_flux_dn_dir" , m_ncol, m_nlay+1);
267  sw_clnsky_flux_up = real2d_k("sw_clnsky_flux_up" , m_ncol, m_nlay+1);
268  sw_clnsky_flux_dn = real2d_k("sw_clnsky_flux_dn" , m_ncol, m_nlay+1);
269  sw_clnsky_flux_dn_dir = real2d_k("sw_clnsky_flux_dn_dir" , m_ncol, m_nlay+1);
270  lw_clnclrsky_flux_up = real2d_k("lw_clnclrsky_flux_up" , m_ncol, m_nlay+1);
271  lw_clnclrsky_flux_dn = real2d_k("lw_clnclrsky_flux_dn" , m_ncol, m_nlay+1);
272  lw_clrsky_flux_up = real2d_k("lw_clrsky_flux_up" , m_ncol, m_nlay+1);
273  lw_clrsky_flux_dn = real2d_k("lw_clrsky_flux_dn" , m_ncol, m_nlay+1);
274  lw_clnsky_flux_up = real2d_k("lw_clnsky_flux_up" , m_ncol, m_nlay+1);
275  lw_clnsky_flux_dn = real2d_k("lw_clnsky_flux_dn" , m_ncol, m_nlay+1);
276 
277  // 3d size (ncol, nlay+1, nswbands)
278  sw_bnd_flux_up = real3d_k("sw_bnd_flux_up" , m_ncol, m_nlay+1, m_nswbands);
279  sw_bnd_flux_dn = real3d_k("sw_bnd_flux_dn" , m_ncol, m_nlay+1, m_nswbands);
280  sw_bnd_flux_dir = real3d_k("sw_bnd_flux_dir", m_ncol, m_nlay+1, m_nswbands);
281  sw_bnd_flux_dif = real3d_k("sw_bnd_flux_dif", m_ncol, m_nlay+1, m_nswbands);
282 
283  // 3d size (ncol, nlay+1, nlwbands)
284  lw_bnd_flux_up = real3d_k("lw_bnd_flux_up" , m_ncol, m_nlay+1, m_nlwbands);
285  lw_bnd_flux_dn = real3d_k("lw_bnd_flux_up" , m_ncol, m_nlay+1, m_nlwbands);
286 
287  // 2d size (ncol, nswbands)
288  sfc_alb_dir = real2d_k("sfc_alb_dir", m_ncol, m_nswbands);
289  sfc_alb_dif = real2d_k("sfc_alb_dif", m_ncol, m_nswbands);
290 
291  // 2d size (ncol, nlwbands)
292  emis_sfc = real2d_k("emis_sfc", m_ncol, m_nlwbands);
293 
294  // 3d size (ncol, nlay, n[sw,lw]bands)
295  aero_tau_sw = real3d_k("aero_tau_sw", m_ncol, m_nlay, m_nswbands);
296  aero_ssa_sw = real3d_k("aero_ssa_sw", m_ncol, m_nlay, m_nswbands);
297  aero_g_sw = real3d_k("aero_g_sw" , m_ncol, m_nlay, m_nswbands);
298  aero_tau_lw = real3d_k("aero_tau_lw", m_ncol, m_nlay, m_nlwbands);
299 
300  // 3d size (ncol, nlay, n[sw,lw]bnds)
301  cld_tau_sw_bnd = real3d_k("cld_tau_sw_bnd", m_ncol, m_nlay, m_nswbands);
302  cld_tau_lw_bnd = real3d_k("cld_tau_lw_bnd", m_ncol, m_nlay, m_nlwbands);
303 
304  // 3d size (ncol, nlay, n[sw,lw]gpts)
305  cld_tau_sw_gpt = real3d_k("cld_tau_sw_gpt", m_ncol, m_nlay, m_nswgpts);
306  cld_tau_lw_gpt = real3d_k("cld_tau_lw_gpt", m_ncol, m_nlay, m_nlwgpts);
307 }
Kokkos::View< RealT *, KokkosDefaultDevice > real1d_k
Definition: ERF_Kokkos.H:17
Kokkos::View< RealT ***, layout_t, KokkosDefaultDevice > real3d_k
Definition: ERF_Kokkos.H:19
Kokkos::View< RealT **, layout_t, KokkosDefaultDevice > real2d_k
Definition: ERF_Kokkos.H:18
Kokkos::View< RealT *, KokkosHostDevice > realHost1d_k
Definition: ERF_Kokkos.H:15
real3d_k sw_bnd_flux_dn
Definition: ERF_Radiation.H:493
real2d_k lw_flux_up
Definition: ERF_Radiation.H:473
real3d_k aero_tau_sw
Definition: ERF_Radiation.H:509
int m_o3_size
Definition: ERF_Radiation.H:365
real3d_k sw_bnd_flux_dir
Definition: ERF_Radiation.H:494
real2d_k sw_clnsky_flux_dn
Definition: ERF_Radiation.H:482
real2d_k d_tint
Definition: ERF_Radiation.H:467
real2d_k lw_clnclrsky_flux_dn
Definition: ERF_Radiation.H:485
real2d_k lwp
Definition: ERF_Radiation.H:459
real3d_k cld_tau_lw_gpt
Definition: ERF_Radiation.H:520
real1d_k lw_src
Definition: ERF_Radiation.H:444
real2d_k eff_radius_qi
Definition: ERF_Radiation.H:457
real1d_k m_gas_mol_weights
Definition: ERF_Radiation.H:366
real2d_k sw_heating
Definition: ERF_Radiation.H:461
real3d_k lw_bnd_flux_dn
Definition: ERF_Radiation.H:499
real3d_k sw_bnd_flux_up
Definition: ERF_Radiation.H:492
real2d_k qv_lay
Definition: ERF_Radiation.H:452
real2d_k sw_clnclrsky_flux_dn_dir
Definition: ERF_Radiation.H:477
real1d_k sfc_flux_dif_vis
Definition: ERF_Radiation.H:438
real2d_k lw_clnclrsky_flux_up
Definition: ERF_Radiation.H:484
real3d_k cld_tau_sw_gpt
Definition: ERF_Radiation.H:519
real1d_k lat
Definition: ERF_Radiation.H:440
real2d_k qi_lay
Definition: ERF_Radiation.H:454
real2d_k sw_clrsky_flux_up
Definition: ERF_Radiation.H:478
real2d_k t_lev
Definition: ERF_Radiation.H:469
real1d_k o3_lay
Definition: ERF_Radiation.H:427
real1d_k sfc_alb_dif_vis
Definition: ERF_Radiation.H:434
real1d_k mu0
Definition: ERF_Radiation.H:431
real2d_k cldfrac_tot
Definition: ERF_Radiation.H:455
real3d_k cld_tau_lw_bnd
Definition: ERF_Radiation.H:516
real1d_k sfc_flux_dir_nir
Definition: ERF_Radiation.H:437
real2d_k sw_clnclrsky_flux_up
Definition: ERF_Radiation.H:475
real3d_k aero_g_sw
Definition: ERF_Radiation.H:511
real2d_k sw_flux_up
Definition: ERF_Radiation.H:470
real3d_k sw_bnd_flux_dif
Definition: ERF_Radiation.H:495
real2d_k sfc_alb_dif
Definition: ERF_Radiation.H:503
real2d_k emis_sfc
Definition: ERF_Radiation.H:506
real1d_k sfc_alb_dif_nir
Definition: ERF_Radiation.H:435
real2d_k lw_clnsky_flux_dn
Definition: ERF_Radiation.H:489
real2d_k p_lay
Definition: ERF_Radiation.H:448
real2d_k r_lay
Definition: ERF_Radiation.H:447
int m_ncol
Definition: ERF_Radiation.H:375
real2d_k sw_flux_dn_dir
Definition: ERF_Radiation.H:472
const std::vector< amrex::Real > m_mol_weight_gas
Definition: ERF_Radiation.H:351
std::vector< std::string > gas_names_offset
Definition: ERF_Radiation.H:367
real2d_k sw_clrsky_flux_dn
Definition: ERF_Radiation.H:479
real3d_k aero_ssa_sw
Definition: ERF_Radiation.H:510
real2d_k sw_clnsky_flux_dn_dir
Definition: ERF_Radiation.H:483
real2d_k qc_lay
Definition: ERF_Radiation.H:453
real2d_k lw_clrsky_flux_up
Definition: ERF_Radiation.H:486
real2d_k sw_flux_dn
Definition: ERF_Radiation.H:471
real2d_k tmp2d
Definition: ERF_Radiation.H:458
real2d_k lw_clrsky_flux_dn
Definition: ERF_Radiation.H:487
real2d_k z_del
Definition: ERF_Radiation.H:450
real2d_k p_del
Definition: ERF_Radiation.H:451
real3d_k aero_tau_lw
Definition: ERF_Radiation.H:512
real2d_k lw_flux_dn
Definition: ERF_Radiation.H:474
real2d_k sw_clrsky_flux_dn_dir
Definition: ERF_Radiation.H:480
real1d_k sfc_flux_dif_nir
Definition: ERF_Radiation.H:439
int m_ngas
Definition: ERF_Radiation.H:348
real2d_k lw_heating
Definition: ERF_Radiation.H:462
real1d_k sfc_alb_dir_nir
Definition: ERF_Radiation.H:433
real2d_k lw_clnsky_flux_up
Definition: ERF_Radiation.H:488
real2d_k sfc_alb_dir
Definition: ERF_Radiation.H:502
real1d_k lon
Definition: ERF_Radiation.H:441
real1d_k sfc_emis
Definition: ERF_Radiation.H:442
real3d_k cld_tau_sw_bnd
Definition: ERF_Radiation.H:515
real2d_k sw_clnsky_flux_up
Definition: ERF_Radiation.H:481
int m_nlay
Definition: ERF_Radiation.H:376
real3d_k lw_bnd_flux_up
Definition: ERF_Radiation.H:498
real1d_k sfc_alb_dir_vis
Definition: ERF_Radiation.H:432
const std::vector< std::string > m_gas_names
Definition: ERF_Radiation.H:349
real2d_k sw_clrsky_heating
Definition: ERF_Radiation.H:463
real2d_k t_lay
Definition: ERF_Radiation.H:449
real1d_k sfc_flux_dir_vis
Definition: ERF_Radiation.H:436
real2d_k iwp
Definition: ERF_Radiation.H:460
real2d_k p_lev
Definition: ERF_Radiation.H:468
real2d_k lw_clrsky_heating
Definition: ERF_Radiation.H:464
real2d_k eff_radius_qc
Definition: ERF_Radiation.H:456
real1d_k cosine_zenith
Definition: ERF_Radiation.H:430
real1d_k t_sfc
Definition: ERF_Radiation.H:443
real2d_k sw_clnclrsky_flux_dn
Definition: ERF_Radiation.H:476

◆ dealloc_buffers()

void Radiation::dealloc_buffers ( )
311 {
312  // 1d size (m_ngas)
314 
315  // 1d size (1 or nlay)
316  o3_lay = real1d_k();
317 
318  // 1d size (ncol)
320  mu0 = real1d_k();
329  lat = real1d_k();
330  lon = real1d_k();
331  sfc_emis = real1d_k();
332  t_sfc = real1d_k();
333  lw_src = real1d_k();
334 
335  // 2d size (ncol, nlay)
336  r_lay = real2d_k();
337  p_lay = real2d_k();
338  t_lay = real2d_k();
339  z_del = real2d_k();
340  p_del = real2d_k();
341  qv_lay = real2d_k();
342  qc_lay = real2d_k();
343  qi_lay = real2d_k();
344  cldfrac_tot = real2d_k();
347  tmp2d = real2d_k();
348  lwp = real2d_k();
349  iwp = real2d_k();
350  sw_heating = real2d_k();
351  lw_heating = real2d_k();
354  sw_heating = real2d_k();
355  lw_heating = real2d_k();
358 
359  // 2d size (ncol, nlay+1)
360  d_tint = real2d_k();
361  p_lev = real2d_k();
362  t_lev = real2d_k();
363  sw_flux_up = real2d_k();
364  sw_flux_dn = real2d_k();
366  lw_flux_up = real2d_k();
367  lw_flux_dn = real2d_k();
383 
384  // 3d size (ncol, nlay+1, nswbands)
389 
390  // 3d size (ncol, nlay+1, nlwbands)
393 
394  // 2d size (ncol, nswbands)
395  sfc_alb_dir = real2d_k();
396  sfc_alb_dif = real2d_k();
397 
398  // 2d size (ncol, nlwbands)
399  emis_sfc = real2d_k();
400 
401  // 3d size (ncol, nlay, n[sw,lw]bands)
402  aero_tau_sw = real3d_k();
403  aero_ssa_sw = real3d_k();
404  aero_g_sw = real3d_k();
405  aero_tau_lw = real3d_k();
406 
407  // 3d size (ncol, nlay, n[sw,lw]bnds)
410 
411  // 3d size (ncol, nlay, n[sw,lw]gpts)
414 }

◆ finalize_impl()

void Radiation::finalize_impl ( amrex::Vector< amrex::MultiFab * > &  lsm_output_ptrs)
1219 {
1220  // Finish rrtmgp
1221  m_gas_concs.reset();
1223 
1224  // Fill the AMReX MFs from Kokkos Views
1225  kokkos_buffers_to_mf(lsm_output_ptrs);
1226 
1227  // Write fluxes if requested
1229 
1230  // Fill output data for datalog before deallocating
1231  if (datalog_int > 0 && m_step % datalog_int == 0) {
1234 
1236  }
1237 
1238  // Deallocate the buffer arrays
1239  dealloc_buffers();
1240 }
int datalog_int
Definition: ERF_RadiationInterface.H:87
void dealloc_buffers()
Definition: ERF_Radiation.cpp:310
int m_step
Definition: ERF_Radiation.H:272
void populateDatalogMF()
Definition: ERF_Radiation.cpp:726
void kokkos_buffers_to_mf(amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
Definition: ERF_Radiation.cpp:596
GasConcsK< amrex::Real, layout_t, KokkosDefaultDevice > m_gas_concs
Definition: ERF_Radiation.H:369
void write_rrtmgp_fluxes()
Definition: ERF_Radiation.cpp:685
void rrtmgp_finalize()
Definition: ERF_RRTMGP_Interface.cpp:266
void compute_heating_rate(View1 const &flux_up, View2 const &flux_dn, View3 const &rho, View4 const &dz, View5 &heating_rate)
Definition: ERF_RRTMGP_Utils.H:81

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.

242  {
243  return m_lsm_input_names;
244  }
amrex::Vector< std::string > m_lsm_input_names
Definition: ERF_Radiation.H:302

◆ get_lsm_output_varnames()

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

Reimplemented from IRadiation.

250  {
251  return m_lsm_output_names;
252  }
amrex::Vector< std::string > m_lsm_output_names
Definition: ERF_Radiation.H:307

◆ Init()

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

Implements IRadiation.

65  {
66  // Ensure the boxes span klo -> khi
67  int klo = geom.Domain().smallEnd(2);
68  int khi = geom.Domain().bigEnd(2);
69 
70  // Reset vector of offsets for columnar data
71  m_nlay = geom.Domain().length(2);
72 
73  m_ncol = 0;
74  m_col_offsets.clear();
75  m_col_offsets.resize(int(ba.size()));
76  for (amrex::MFIter mfi(*cons_in); mfi.isValid(); ++mfi) {
77  const amrex::Box& vbx = mfi.validbox();
78  AMREX_ALWAYS_ASSERT_WITH_MESSAGE((klo == vbx.smallEnd(2)) &&
79  (khi == vbx.bigEnd(2)),
80  "Vertical decomposition with radiation is not allowed.");
81  int nx = vbx.length(0);
82  int ny = vbx.length(1);
83  m_col_offsets[mfi.index()] = m_ncol;
84  m_ncol += nx * ny;
85  }
86  };
amrex::Vector< int > m_col_offsets
Definition: ERF_Radiation.H:379

◆ initialize_impl()

void Radiation::initialize_impl ( )
951 {
952  // Call API to initialize
957 }
void rrtmgp_initialize(gas_concs_t &gas_concs_k, const std::string &coefficients_file_sw, const std::string &coefficients_file_lw, const std::string &cloud_optics_file_sw, const std::string &cloud_optics_file_lw)
Definition: ERF_RRTMGP_Interface.cpp:229

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)
597 {
598  // Heating rate, fluxes, zenith, lsm ptrs
599  Vector<real2d_k> rrtmgp_out_vars = {sw_flux_dn, lw_flux_dn};
600 
601  // Expose for device
602  auto sw_heating_d = sw_heating;
603  auto lw_heating_d = lw_heating;
604  auto p_lay_d = p_lay;
605  auto sfc_flux_dir_vis_d = sfc_flux_dir_vis;
606  auto sfc_flux_dir_nir_d = sfc_flux_dir_nir;
607  auto sfc_flux_dif_vis_d = sfc_flux_dif_vis;
608  auto sfc_flux_dif_nir_d = sfc_flux_dif_nir;
609  auto lw_flux_dn_d = lw_flux_dn;
610  auto mu0_d = mu0;
611 
612  for (MFIter mfi(*m_cons_in); mfi.isValid(); ++mfi) {
613  const auto& vbx = mfi.validbox();
614  const auto& sbx = makeSlab(vbx,2,vbx.smallEnd(2));
615  const int nx = vbx.length(0);
616  const int imin = vbx.smallEnd(0);
617  const int jmin = vbx.smallEnd(1);
618  const int offset = m_col_offsets[mfi.index()];
619  const Array4<Real>& q_arr = m_qheating_rates->array(mfi);
620  ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
621  {
622  // map [i,j,k] 0-based to [icol, ilay] 0-based
623  const int icol = (j-jmin)*nx + (i-imin) + offset;
624  const int ilay = k;
625 
626  // Temperature heating rate for SW and LW
627  q_arr(i,j,k,0) = sw_heating_d(icol,ilay);
628  q_arr(i,j,k,1) = lw_heating_d(icol,ilay);
629 
630  // Convert the rates for theta_d
631  Real exner = getExnergivenP(Real(p_lay_d(icol,ilay)), R_d/Cp_d);
632  q_arr(i,j,k,0) *= exner;
633  q_arr(i,j,k,1) *= exner;
634  });
635  if (m_lsm_fluxes) {
636  const Array4<Real>& lsm_arr = m_lsm_fluxes->array(mfi);
637  ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
638  {
639  // map [i,j,k] 0-based to [icol, ilay] 0-based
640  const int icol = (j-jmin)*nx + (i-imin) + offset;
641 
642  // SW fluxes for LSM
643  lsm_arr(i,j,k,0) = sfc_flux_dir_vis_d(icol);
644  lsm_arr(i,j,k,1) = sfc_flux_dir_nir_d(icol);
645  lsm_arr(i,j,k,2) = sfc_flux_dif_vis_d(icol);
646  lsm_arr(i,j,k,3) = sfc_flux_dif_nir_d(icol);
647 
648  // Net SW flux for LSM
649  lsm_arr(i,j,k,4) = sfc_flux_dir_vis_d(icol) + sfc_flux_dir_nir_d(icol)
650  + sfc_flux_dif_vis_d(icol) + sfc_flux_dif_nir_d(icol);
651 
652  // LW flux for LSM (at bottom surface)
653  lsm_arr(i,j,k,5) = lw_flux_dn_d(icol,1);
654  });
655  }
656  if (m_lsm_zenith) {
657  const Array4<Real>& lsm_zenith_arr = m_lsm_zenith->array(mfi);
658  ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
659  {
660  // map [i,j,k] 0-based to [icol, ilay] 0-based
661  const int icol = (j-jmin)*nx + (i-imin) + offset;
662 
663  // export cosine zenith angle for LSM
664  lsm_zenith_arr(i,j,k) = mu0_d(icol);
665  });
666  }
667  for (int ivar(0); ivar<lsm_output_ptrs.size(); ivar++) {
668  auto rrtmgp_for_fill = rrtmgp_out_vars[ivar];
669  if (lsm_output_ptrs[ivar]) {
670  const Array4<Real>& lsm_out_arr = lsm_output_ptrs[ivar]->array(mfi);
671  ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
672  {
673  // map [i,j,k] 0-based to [icol, ilay] 0-based
674  const int icol = (j-jmin)*nx + (i-imin) + offset;
675 
676  // export the desired variable at surface
677  lsm_out_arr(i,j,k) = rrtmgp_for_fill(icol,0);
678  });
679  } // valid ptr
680  } // ivar
681  }// mfi
682 }
constexpr amrex::Real Cp_d
Definition: ERF_Constants.H:12
constexpr amrex::Real R_d
Definition: ERF_Constants.H:10
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getExnergivenP(const amrex::Real P, const amrex::Real rdOcp)
Definition: ERF_EOS.H:141
AMREX_FORCE_INLINE IntVect offset(const int face_dir, const int normal)
Definition: ERF_ReadBndryPlanes.cpp:28
amrex::MultiFab * m_cons_in
Definition: ERF_Radiation.H:313
amrex::MultiFab * m_lsm_fluxes
Definition: ERF_Radiation.H:330
amrex::MultiFab * m_qheating_rates
Definition: ERF_Radiation.H:316
amrex::MultiFab * m_lsm_zenith
Definition: ERF_Radiation.H:331
Here is the call graph for this function:

◆ mf_to_kokkos_buffers()

void Radiation::mf_to_kokkos_buffers ( amrex::Vector< amrex::MultiFab * > &  lsm_input_ptrs)
419 {
420  // Expose for device
421  auto r_lay_d = r_lay;
422  auto p_lay_d = p_lay;
423  auto t_lay_d = t_lay;
424  auto z_del_d = z_del;
425  auto qv_lay_d = qv_lay;
426  auto qc_lay_d = qc_lay;
427  auto qi_lay_d = qi_lay;
428  auto cldfrac_tot_d = cldfrac_tot;
429  auto lwp_d = lwp;
430  auto iwp_d = iwp;
431  auto eff_radius_qc_d = eff_radius_qc;
432  auto eff_radius_qi_d = eff_radius_qi;
433  auto p_lev_d = p_lev;
434  auto t_lev_d = t_lev;
435  auto lat_d = lat;
436  auto lon_d = lon;
437  auto t_sfc_d = t_sfc;
438  auto p_del_d = p_del;
439 
440  bool moist = m_moist;
441  bool ice = m_ice;
442  const bool has_lsm = m_lsm;
443  const bool has_lat = m_lat;
444  const bool has_lon = m_lon;
445  int ncol = m_ncol;
446  int nlay = m_nlay;
447  Real dz = m_geom.CellSize(2);
448  Real cons_lat = m_lat_cons;
449  Real cons_lon = m_lon_cons;
450  Real rad_t_sfc = m_rad_t_sfc;
451  for (MFIter mfi(*m_cons_in); mfi.isValid(); ++mfi) {
452  const auto& vbx = mfi.validbox();
453  const int nx = vbx.length(0);
454  const int imin = vbx.smallEnd(0);
455  const int jmin = vbx.smallEnd(1);
456  const int offset = m_col_offsets[mfi.index()];
457  const Array4<const Real>& cons_arr = m_cons_in->const_array(mfi);
458  const Array4<const Real>& z_arr = (m_z_phys) ? m_z_phys->const_array(mfi) :
459  Array4<const Real>{};
460  const Array4<const Real>& lat_arr = (m_lat) ? m_lat->const_array(mfi) :
461  Array4<const Real>{};
462  const Array4<const Real>& lon_arr = (m_lon) ? m_lon->const_array(mfi) :
463  Array4<const Real>{};
464  ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
465  {
466  // map [i,j,k] 0-based to [icol, ilay] 0-based
467  const int icol = (j-jmin)*nx + (i-imin) + offset;
468  const int ilay = k;
469 
470  // EOS input (at CC)
471  Real r = cons_arr(i,j,k,Rho_comp);
472  Real rt = cons_arr(i,j,k,RhoTheta_comp);
473  Real qv = (moist) ? cons_arr(i,j,k,RhoQ1_comp)/r : 0.0;
474  Real qc = (moist) ? cons_arr(i,j,k,RhoQ2_comp)/r : 0.0;
475  Real qi = (ice) ? cons_arr(i,j,k,RhoQ3_comp)/r : 0.0;
476 
477  // EOS avg to z-face
478  Real r_lo = cons_arr(i,j,k-1,Rho_comp);
479  Real rt_lo = cons_arr(i,j,k-1,RhoTheta_comp);
480  Real qv_lo = (moist) ? cons_arr(i,j,k-1,RhoQ1_comp)/r_lo : 0.0;
481 
482  // make sure qv is >= 0 for H2O VMR
483  qv_lo = std::max(0.0, qv_lo);
484  qv = std::max(0.0, qv);
485 
486  Real r_avg = 0.5 * (r + r_lo);
487  Real rt_avg = 0.5 * (rt + rt_lo);
488  Real qv_avg = 0.5 * (qv + qv_lo);
489 
490  // Views at CC
491  r_lay_d(icol,ilay) = r;
492  p_lay_d(icol,ilay) = getPgivenRTh(rt, qv);
493  t_lay_d(icol,ilay) = getTgivenRandRTh(r, rt, qv);
494  z_del_d(icol,ilay) = (z_arr) ? 0.25 * ( (z_arr(i ,j ,k+1) - z_arr(i ,j ,k))
495  + (z_arr(i+1,j ,k+1) - z_arr(i+1,j ,k))
496  + (z_arr(i ,j+1,k+1) - z_arr(i ,j+1,k))
497  + (z_arr(i+1,j ,k+1) - z_arr(i+1,j ,k)) ) : dz;
498  qv_lay_d(icol,ilay) = qv;
499  qc_lay_d(icol,ilay) = qc;
500  qi_lay_d(icol,ilay) = qi;
501  cldfrac_tot_d(icol,ilay) = ((qc+qi)>0.0) ? 1. : 0.;
502 
503  // NOTE: These are populated in 'mixing_ratio_to_cloud_mass'
504  lwp_d(icol,ilay) = 0.0;
505  iwp_d(icol,ilay) = 0.0;
506 
507  // NOTE: These would be populated from P3 (we use the constants in p3_main_impl.hpp)
508  eff_radius_qc_d(icol,ilay) = (qc>0.0) ? 10.0e-6 : 0.0;
509  eff_radius_qi_d(icol,ilay) = (qi>0.0) ? 25.0e-6 : 0.0;
510 
511  // Buffers on z-faces (nlay+1)
512  p_lev_d(icol,ilay) = getPgivenRTh(rt_avg, qv_avg);
513  t_lev_d(icol,ilay) = getTgivenRandRTh(r_avg, rt_avg, qv_avg);
514  if (ilay==(nlay-1)) {
515  Real r_hi = cons_arr(i,j,k+1,Rho_comp);
516  Real rt_hi = cons_arr(i,j,k+1,RhoTheta_comp);
517  Real qv_hi = (moist) ? cons_arr(i,j,k+1,RhoQ1_comp)/r_hi : 0.0;
518  qv_hi = std::max(0.0, qv_hi);
519  r_avg = 0.5 * (r + r_hi);
520  rt_avg = 0.5 * (rt + rt_hi);
521  qv_avg = 0.5 * (qv + qv_hi);
522  p_lev_d(icol,ilay+1) = getPgivenRTh(rt_avg, qv_avg);
523  t_lev_d(icol,ilay+1) = getTgivenRandRTh(r_avg, rt_avg, qv_avg);
524  }
525 
526  // 1D data structures
527  if (k==0) {
528  lat_d(icol) = (has_lat) ? lat_arr(i,j,0) : cons_lat;
529  lon_d(icol) = (has_lon) ? lon_arr(i,j,0) : cons_lon;
530  }
531 
532  });
533  } // mfi
534 
535  // EAMXX delP is positive
536  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
537  KOKKOS_LAMBDA (int icol, int ilay)
538  {
539  p_del_d(icol,ilay) = std::abs(p_lev_d(icol,ilay+1) - p_lev_d(icol,ilay));
540  });
541 
542  // Populate vars LSM would provide
543  if (!has_lsm) {
544  // Parsed surface temp
545  Kokkos::deep_copy(t_sfc, rad_t_sfc);
546 
547  // EAMXX dummy atmos constants
548  Kokkos::deep_copy(sfc_alb_dir_vis, 0.06);
549  Kokkos::deep_copy(sfc_alb_dir_nir, 0.06);
550  Kokkos::deep_copy(sfc_alb_dif_vis, 0.06);
551  Kokkos::deep_copy(sfc_alb_dif_nir, 0.06);
552 
553  // AML NOTE: These are not used in current EAMXX, I've left
554  // the code to plug into these if we need it.
555  //
556  // Current EAMXX constants
557  Kokkos::deep_copy(emis_sfc, 0.98);
558  Kokkos::deep_copy(lw_src , 0.0 );
559  } else {
560  Vector<real1d_k> rrtmgp_in_vars = {t_sfc, sfc_emis,
563  Vector<Real> rrtmgp_default_vals = {rad_t_sfc, 0.98,
564  0.06, 0.06,
565  0.06, 0.06};
566  for (int ivar(0); ivar<lsm_input_ptrs.size(); ivar++) {
567  auto rrtmgp_to_fill = rrtmgp_in_vars[ivar];
568  if (!lsm_input_ptrs[ivar]) {
569  Kokkos::deep_copy(rrtmgp_to_fill, rrtmgp_default_vals[ivar]);
570  } else {
571  for (MFIter mfi(*m_cons_in); mfi.isValid(); ++mfi) {
572  const auto& vbx = mfi.validbox();
573  const auto& sbx = makeSlab(vbx,2,vbx.smallEnd(2));
574  const int nx = vbx.length(0);
575  const int imin = vbx.smallEnd(0);
576  const int jmin = vbx.smallEnd(1);
577  const int offset = m_col_offsets[mfi.index()];
578  const Array4<const Real>& lsm_in_arr = lsm_input_ptrs[ivar]->const_array(mfi);
579  ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
580  {
581  // map [i,j,k] 0-based to [icol, ilay] 0-based
582  const int icol = (j-jmin)*nx + (i-imin) + offset;
583 
584  // 2D mf and 1D view
585  rrtmgp_to_fill(icol) = lsm_in_arr(i,j,k);
586  });
587  } //mfi
588  } // valid lsm ptr
589  } // ivar
590  Kokkos::deep_copy(lw_src, 0.0 );
591  } // have lsm
592 }
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getPgivenRTh(const amrex::Real rhotheta, const amrex::Real qv=0.)
Definition: ERF_EOS.H:81
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getTgivenRandRTh(const amrex::Real rho, const amrex::Real rhotheta, const amrex::Real qv=0.0)
Definition: ERF_EOS.H:46
#define Rho_comp
Definition: ERF_IndexDefines.H:36
#define RhoTheta_comp
Definition: ERF_IndexDefines.H:37
#define RhoQ2_comp
Definition: ERF_IndexDefines.H:43
#define RhoQ3_comp
Definition: ERF_IndexDefines.H:44
#define RhoQ1_comp
Definition: ERF_IndexDefines.H:42
amrex::MultiFab * m_lon
Definition: ERF_Radiation.H:323
amrex::MultiFab * m_z_phys
Definition: ERF_Radiation.H:319
amrex::Geometry m_geom
Definition: ERF_Radiation.H:281
amrex::MultiFab * m_lat
Definition: ERF_Radiation.H:322
@ qv
Definition: ERF_Kessler.H:28
@ qc
Definition: ERF_SatAdj.H:36
Here is the call graph for this function:

◆ populateDatalogMF()

void Radiation::populateDatalogMF ( )
727 {
728  // Expose for device
729  auto sw_flux_up_d = sw_flux_up;
730  auto sw_flux_dn_d = sw_flux_dn;
731  auto sw_flux_dn_dir_d = sw_flux_dn_dir;
732  auto lw_flux_up_d = lw_flux_up;
733  auto lw_flux_dn_d = lw_flux_dn;
734  auto mu0_d = mu0;
735  auto sw_clrsky_heating_d = sw_clrsky_heating;
736  auto lw_clrsky_heating_d = lw_clrsky_heating;
737 
738  auto sw_clrsky_flux_up_d = sw_clrsky_flux_up;
739  auto sw_clrsky_flux_dn_d = sw_clrsky_flux_dn;
740  auto sw_clrsky_flux_dn_dir_d = sw_clrsky_flux_dn_dir;
741  auto lw_clrsky_flux_up_d = lw_clrsky_flux_up;
742  auto lw_clrsky_flux_dn_d = lw_clrsky_flux_dn;
743  auto sw_clnsky_flux_up_d = sw_clnsky_flux_up;
744  auto sw_clnsky_flux_dn_d = sw_clnsky_flux_dn;
745  auto sw_clnsky_flux_dn_dir_d = sw_clnsky_flux_dn_dir;
746  auto lw_clnsky_flux_up_d = lw_clnsky_flux_up;
747  auto lw_clnsky_flux_dn_d = lw_clnsky_flux_dn;
748  auto sw_clnclrsky_flux_up_d = sw_clnclrsky_flux_up;
749  auto sw_clnclrsky_flux_dn_d = sw_clnclrsky_flux_dn;
750  auto sw_clnclrsky_flux_dn_dir_d = sw_clnclrsky_flux_dn_dir;
751  auto lw_clnclrsky_flux_up_d = lw_clnclrsky_flux_up;
752  auto lw_clnclrsky_flux_dn_d = lw_clnclrsky_flux_dn;
753 
754  auto extra_clnsky_diag = m_extra_clnsky_diag;
755  auto extra_clnclrsky_diag = m_extra_clnclrsky_diag;
756 
757  for (MFIter mfi(datalog_mf); mfi.isValid(); ++mfi) {
758  const auto& vbx = mfi.validbox();
759  const int nx = vbx.length(0);
760  const int imin = vbx.smallEnd(0);
761  const int jmin = vbx.smallEnd(1);
762  const int offset = m_col_offsets[mfi.index()];
763  const Array4<Real>& dst_arr = datalog_mf.array(mfi);
764  const Array4<Real>& q_arr = m_qheating_rates->array(mfi);
765  ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
766  {
767  // map [i,j,k] 0-based to [icol, ilay] 0-based
768  const int icol = (j-jmin)*nx + (i-imin) + offset;
769  const int ilay = k;
770 
771  dst_arr(i,j,k,0) = q_arr(i, j, k, 0);
772  dst_arr(i,j,k,1) = q_arr(i, j, k, 1);
773 
774  // SW and LW fluxes
775  dst_arr(i,j,k,2) = sw_flux_up_d(icol,ilay);
776  dst_arr(i,j,k,3) = sw_flux_dn_d(icol,ilay);
777  dst_arr(i,j,k,4) = sw_flux_dn_dir_d(icol,ilay);
778  dst_arr(i,j,k,5) = lw_flux_up_d(icol,ilay);
779  dst_arr(i,j,k,6) = lw_flux_dn_d(icol,ilay);
780 
781  // Cosine zenith angle
782  dst_arr(i,j,k,7) = mu0_d(icol);
783 
784  // Clear sky heating rates and fluxes:
785  dst_arr(i,j,k,8) = sw_clrsky_heating_d(icol, ilay);
786  dst_arr(i,j,k,9) = lw_clrsky_heating_d(icol, ilay);
787 
788  dst_arr(i,j,k,10) = sw_clrsky_flux_up_d(icol,ilay);
789  dst_arr(i,j,k,11) = sw_clrsky_flux_dn_d(icol,ilay);
790  dst_arr(i,j,k,12) = sw_clrsky_flux_dn_dir_d(icol,ilay);
791  dst_arr(i,j,k,13) = lw_clrsky_flux_up_d(icol,ilay);
792  dst_arr(i,j,k,14) = lw_clrsky_flux_dn_d(icol,ilay);
793 
794  // Clean sky fluxes:
795  if (extra_clnsky_diag) {
796  dst_arr(i,j,k,15) = sw_clnsky_flux_up_d(icol,ilay);
797  dst_arr(i,j,k,16) = sw_clnsky_flux_dn_d(icol,ilay);
798  dst_arr(i,j,k,17) = sw_clnsky_flux_dn_dir_d(icol,ilay);
799  dst_arr(i,j,k,18) = lw_clnsky_flux_up_d(icol,ilay);
800  dst_arr(i,j,k,19) = lw_clnsky_flux_dn_d(icol,ilay);
801  }
802 
803  // Clean-clear sky fluxes:
804  if (extra_clnclrsky_diag) {
805  dst_arr(i,j,k,20) = sw_clnclrsky_flux_up_d(icol,ilay);
806  dst_arr(i,j,k,21) = sw_clnclrsky_flux_dn_d(icol,ilay);
807  dst_arr(i,j,k,22) = sw_clnclrsky_flux_dn_dir_d(icol,ilay);
808  dst_arr(i,j,k,23) = lw_clnclrsky_flux_up_d(icol,ilay);
809  dst_arr(i,j,k,24) = lw_clnclrsky_flux_dn_d(icol,ilay);
810  }
811  });
812  }
813 }
amrex::MultiFab datalog_mf
Definition: ERF_Radiation.H:334
Here is the call graph for this function:

◆ rad_run_impl()

void Radiation::rad_run_impl ( amrex::Vector< amrex::MultiFab * > &  lsm_output_ptrs)
inline
228  {
229  if (m_update_rad) {
230  amrex::Print() << "Advancing radiation at level: " << m_lev << " ...";
231  this->initialize_impl();
232  this->run_impl();
233  this->finalize_impl(lsm_output_ptrs);
234  amrex::Print() << "DONE\n";
235  }
236  }
void initialize_impl()
Definition: ERF_Radiation.cpp:950
void run_impl()
Definition: ERF_Radiation.cpp:961
void finalize_impl(amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
Definition: ERF_Radiation.cpp:1218
int m_lev
Definition: ERF_Radiation.H:269
bool m_update_rad
Definition: ERF_Radiation.H:287

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,
amrex::Real &  time,
const amrex::Real &  dt,
const amrex::BoxArray &  ba,
amrex::Geometry &  geom,
amrex::MultiFab *  cons_in,
amrex::MultiFab *  lsm_fluxes,
amrex::MultiFab *  lsm_zenith,
amrex::Vector< amrex::MultiFab * > &  lsm_input_ptrs,
amrex::Vector< amrex::MultiFab * > &  lsm_output_ptrs,
amrex::MultiFab *  qheating_rates,
amrex::MultiFab *  z_phys,
amrex::MultiFab *  lat,
amrex::MultiFab *  lon 
)
inlineoverridevirtual

Implements IRadiation.

105  {
106  set_grids(level, step, time, dt, ba, geom,
107  cons_in, lsm_fluxes, lsm_zenith,
108  lsm_input_ptrs, qheating_rates,
109  z_phys, lat, lon);
110  rad_run_impl(lsm_output_ptrs);
111  }
void rad_run_impl(amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
Definition: ERF_Radiation.H:227
void set_grids(int &level, int &step, amrex::Real &time, const amrex::Real &dt, const amrex::BoxArray &ba, amrex::Geometry &geom, amrex::MultiFab *cons_in, amrex::MultiFab *lsm_fluxes, amrex::MultiFab *lsm_zenith, amrex::Vector< amrex::MultiFab * > &lsm_input_ptrs, amrex::MultiFab *qheating_rates, amrex::MultiFab *z_phys, amrex::MultiFab *lat, amrex::MultiFab *lon)
Definition: ERF_Radiation.cpp:117
Here is the call graph for this function:

◆ run_impl()

void Radiation::run_impl ( )
962 {
963  // Local copies
964  const auto ncol = m_ncol;
965  const auto nlay = m_nlay;
966  const auto nswbands = m_nswbands;
967 
968  // Compute orbital parameters; these are used both for computing
969  // the solar zenith angle and also for computing total solar
970  // irradiance scaling (tsi_scaling).
971  double obliqr, lambm0, mvelpp;
972  int orbital_year = m_orbital_year;
973  double eccen = m_orbital_eccen;
974  double obliq = m_orbital_obliq;
975  double mvelp = m_orbital_mvelp;
976  if (eccen >= 0 && obliq >= 0 && mvelp >= 0) {
977  // fixed orbital parameters forced with orbital_year == ORB_UNDEF_INT
978  orbital_year = ORB_UNDEF_INT;
979  }
980  orbital_params(orbital_year, eccen, obliq,
981  mvelp, obliqr, lambm0, mvelpp);
982 
983 
984  // Use the orbital parameters to calculate the solar declination and eccentricity factor
985  double delta, eccf;
986  // Want day + fraction; calday 1 == Jan 1 0Z
987  static constexpr double dpy[] = {0.0 , 31.0, 59.0, 90.0, 120.0, 151.0,
988  181.0, 212.0, 243.0, 273.0, 304.0, 334.0};
989  bool leap = (m_orbital_year % 4 == 0 && (!(m_orbital_year % 100 == 0) || (m_orbital_year % 400 == 0))) ? true : false;
990  double calday = dpy[m_orbital_mon-1] + (m_orbital_day-1.0) + m_orbital_sec/86400.0;
991  // add extra day if leap year
992  if (leap) { calday += 1.0; }
993  orbital_decl(calday, eccen, mvelpp, lambm0, obliqr, delta, eccf);
994 
995  // Overwrite eccf if using a fixed solar constant.
996  auto fixed_total_solar_irradiance = m_fixed_total_solar_irradiance;
997  if (fixed_total_solar_irradiance >= 0){
998  eccf = fixed_total_solar_irradiance/1360.9;
999  }
1000 
1001  // Precompute volume mixing ratio (VMR) for all gases
1002  //
1003  // H2O is obtained from qv.
1004  // O3 may be a constant or a 1D vector
1005  // All other comps are set to constants for now
1006  for (int igas(0); igas < m_ngas; ++igas) {
1007  auto tmp2d_d = tmp2d;
1008  auto name = m_gas_names[igas];
1009  auto gas_mol_weight = m_mol_weight_gas[igas];
1010  if (name == "H2O") {
1011  auto qv_lay_d = qv_lay;
1012  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
1013  KOKKOS_LAMBDA (int icol, int ilay)
1014  {
1015  tmp2d_d(icol,ilay) = qv_lay_d(icol,ilay) * mwdair/ gas_mol_weight;
1016  });
1017  } else if (name == "CO2") {
1018  Kokkos::deep_copy(tmp2d, m_co2vmr);
1019  } else if (name == "O3") {
1020  if (m_o3_size==1) {
1021  Kokkos::deep_copy(tmp2d, m_o3vmr[0] );
1022  } else {
1023  auto o3_lay_d = o3_lay;
1024  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
1025  KOKKOS_LAMBDA (int icol, int ilay)
1026  {
1027  tmp2d_d(icol,ilay) = o3_lay_d(ilay);
1028  });
1029  }
1030  } else if (name == "N2O") {
1031  Kokkos::deep_copy(tmp2d, m_n2ovmr);
1032  } else if (name == "CO") {
1033  Kokkos::deep_copy(tmp2d, m_covmr );
1034  } else if (name == "CH4") {
1035  Kokkos::deep_copy(tmp2d, m_ch4vmr);
1036  } else if (name == "O2") {
1037  Kokkos::deep_copy(tmp2d, m_o2vmr );
1038  } else if (name == "N2") {
1039  Kokkos::deep_copy(tmp2d, m_n2vmr );
1040  } else {
1041  Abort("Radiation: Unknown gas component.");
1042  }
1043 
1044  // Populate GasConcs object
1045  m_gas_concs.set_vmr(name, tmp2d);
1046  }
1047 
1048  // Populate mu0 1D array
1049  // This must be done on HOST and copied to device.
1050  auto h_mu0 = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), mu0);
1051  if (m_fixed_solar_zenith_angle > 0) {
1052  Kokkos::deep_copy(h_mu0, m_fixed_solar_zenith_angle);
1053  } else {
1054  auto h_lat = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), lat);
1055  auto h_lon = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), lon);
1056  double dt = double(m_dt);
1057  auto rad_freq_in_steps = m_rad_freq_in_steps;
1058  Kokkos::parallel_for(Kokkos::RangePolicy<Kokkos::Serial>(0, ncol),
1059  KOKKOS_LAMBDA (int icol)
1060  {
1061  // Convert lat/lon to radians
1062  double lat_col = h_lat(icol)*PI/180.0;
1063  double lon_col = h_lon(icol)*PI/180.0;
1064  double lcalday = calday;
1065  double ldelta = delta;
1066  h_mu0(icol) = Real(orbital_cos_zenith(lcalday, lat_col, lon_col, ldelta, rad_freq_in_steps * dt));
1067  });
1068  }
1069  Kokkos::deep_copy(mu0, h_mu0);
1070 
1071  // Compute layer cloud mass per unit area (populates lwp/iwp)
1074 
1075  // Convert to g/m2 (needed by RRTMGP)
1076  auto lwp_d = lwp;
1077  auto iwp_d = iwp;
1078  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
1079  KOKKOS_LAMBDA (int icol, int ilay)
1080  {
1081  lwp_d(icol,ilay) *= 1.e3;
1082  iwp_d(icol,ilay) *= 1.e3;
1083  });
1084 
1085  // Expand surface_albedos along nswbands.
1086  // This is needed since rrtmgp require band-by-band.
1091 
1092  // Run RRTMGP driver
1094  p_lay, t_lay,
1095  p_lev, t_lev,
1096  m_gas_concs,
1098  t_sfc, emis_sfc, lw_src,
1114 
1115 #if 0
1116  // UNIT TEST
1117  //================================================================================
1118  Kokkos::deep_copy(mu0, 0.86);
1119  Kokkos::deep_copy(sfc_alb_dir_vis, 0.06);
1120  Kokkos::deep_copy(sfc_alb_dir_nir, 0.06);
1121  Kokkos::deep_copy(sfc_alb_dif_vis, 0.06);
1122  Kokkos::deep_copy(sfc_alb_dif_nir, 0.06);
1123 
1124  Kokkos::deep_copy(aero_tau_sw, 0.0);
1125  Kokkos::deep_copy(aero_ssa_sw, 0.0);
1126  Kokkos::deep_copy(aero_g_sw , 0.0);
1127  Kokkos::deep_copy(aero_tau_lw, 0.0);
1128 
1129  // Generate some fake liquid and ice water data. We pick values to be midway between
1130  // the min and max of the valid lookup table values for effective radii
1131  real rel_val = 0.5 * (rrtmgp::cloud_optics_sw_k->get_min_radius_liq()
1132  + rrtmgp::cloud_optics_sw_k->get_max_radius_liq());
1133  real rei_val = 0.5 * (rrtmgp::cloud_optics_sw_k->get_min_radius_ice()
1134  + rrtmgp::cloud_optics_sw_k->get_max_radius_ice());
1135 
1136  // Restrict clouds to troposphere (> 100 hPa = 100*100 Pa) and not very close to the ground (< 900 hPa), and
1137  // put them in 2/3 of the columns since that's roughly the total cloudiness of earth.
1138  // Set sane values for liquid and ice water path.
1139  // NOTE: these "sane" values are in g/m2!
1140  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
1141  KOKKOS_LAMBDA (int icol, int ilay)
1142  {
1143  cldfrac_tot(icol,ilay) = (p_lay(icol,ilay) > 100. * 100.) &&
1144  (p_lay(icol,ilay) < 900. * 100.) &&
1145  (icol%3 != 0);
1146  // Ice and liquid will overlap in a few layers
1147  lwp(icol,ilay) = (cldfrac_tot(icol,ilay) && t_lay(icol,ilay) > 263.) ? 10. : 0.;
1148  iwp(icol,ilay) = (cldfrac_tot(icol,ilay) && t_lay(icol,ilay) < 273.) ? 10. : 0.;
1149  eff_radius_qc(icol,ilay) = (lwp(icol,ilay) > 0.) ? rel_val : 0.;
1150  eff_radius_qi(icol,ilay) = (iwp(icol,ilay) > 0.) ? rei_val : 0.;
1151  });
1152 
1157 
1159  p_lay, t_lay,
1160  p_lev, t_lev,
1161  m_gas_concs,
1163  t_sfc, emis_sfc, lw_src,
1178  1.0, false, false);
1179  //================================================================================
1180 #endif
1181 
1182  // Update heating tendency
1185 
1186  /*
1187  // AML DEBUG
1188  Kokkos::parallel_for(nlay+1, KOKKOS_LAMBDA (int ilay)
1189  {
1190  printf("Fluxes: %i %e %e %e %e %e\n",ilay,
1191  sw_flux_up(5,ilay), sw_flux_dn(5,ilay), sw_flux_dn_dir(5,ilay),
1192  lw_flux_up(5,ilay), lw_flux_dn(5,ilay));
1193  });
1194  Kokkos::parallel_for(nlay, KOKKOS_LAMBDA (int ilay)
1195  {
1196  printf("Heating Rate: %i %e %e\n",ilay,sw_heating(5,ilay),lw_heating(5,ilay));
1197  });
1198  */
1199 
1200  // Compute surface fluxes
1201  const int kbot = 0;
1202  auto sw_bnd_flux_dif_d = sw_bnd_flux_dif;
1203  auto sw_bnd_flux_dn_d = sw_bnd_flux_dn;
1204  auto sw_bnd_flux_dir_d = sw_bnd_flux_dir;
1205  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<3>>({0, 0, 0}, {ncol, nlay+1, nswbands}),
1206  KOKKOS_LAMBDA (int icol, int ilay, int ibnd)
1207  {
1208  sw_bnd_flux_dif_d(icol,ilay,ibnd) = sw_bnd_flux_dn_d(icol,ilay,ibnd) - sw_bnd_flux_dir_d(icol,ilay,ibnd);
1209  });
1210  rrtmgp::compute_broadband_surface_fluxes(ncol, kbot, nswbands,
1214 }
static constexpr int ORB_UNDEF_INT
Definition: ERF_Constants.H:96
constexpr amrex::Real mwdair
Definition: ERF_Constants.H:64
constexpr amrex::Real PI
Definition: ERF_Constants.H:6
double real
Definition: ERF_OrbCosZenith.H:9
AMREX_GPU_HOST AMREX_FORCE_INLINE real orbital_cos_zenith(real &jday, real &lat, real &lon, real &declin, real dt_avg=-1., real uniform_angle=-1., real constant_zenith_angle_deg=-1.)
Definition: ERF_OrbCosZenith.H:559
AMREX_GPU_HOST AMREX_FORCE_INLINE void orbital_decl(real &calday, real &eccen, real &mvelpp, real &lambm0, real &obliqr, real &delta, real &eccf)
Definition: ERF_OrbCosZenith.H:15
AMREX_GPU_HOST AMREX_FORCE_INLINE void orbital_params(int &iyear_AD, real &eccen, real &obliq, real &mvelp, real &obliqr, real &lambm0, real &mvelpp)
Definition: ERF_OrbCosZenith.H:81
amrex::Real m_dt
Definition: ERF_Radiation.H:278
int m_orbital_mon
Definition: ERF_Radiation.H:392
int m_orbital_sec
Definition: ERF_Radiation.H:394
int m_orbital_day
Definition: ERF_Radiation.H:393
void compute_band_by_band_surface_albedos(const int ncol, const int nswbands, real1d_k &sfc_alb_dir_vis, real1d_k &sfc_alb_dir_nir, real1d_k &sfc_alb_dif_vis, real1d_k &sfc_alb_dif_nir, real2d_k &sfc_alb_dir, real2d_k &sfc_alb_dif)
Definition: ERF_RRTMGP_Interface.cpp:282
void compute_broadband_surface_fluxes(const int ncol, const int kbot, const int nswbands, real3d_k &sw_bnd_flux_dir, real3d_k &sw_bnd_flux_dif, real1d_k &sfc_flux_dir_vis, real1d_k &sfc_flux_dir_nir, real1d_k &sfc_flux_dif_vis, real1d_k &sfc_flux_dif_nir)
Definition: ERF_RRTMGP_Interface.cpp:326
void mixing_ratio_to_cloud_mass(View1 const &mixing_ratio, View2 const &cloud_fraction, View3 const &rho, View4 const &dz, View5 const &cloud_mass)
Definition: ERF_RRTMGP_Utils.H:12
std::unique_ptr< cloud_optics_t > cloud_optics_sw_k
Definition: ERF_RRTMGP_Interface.cpp:18
void rrtmgp_main(const int ncol, const int nlay, real2d_k &p_lay, real2d_k &t_lay, real2d_k &p_lev, real2d_k &t_lev, gas_concs_t &gas_concs, real2d_k &sfc_alb_dir, real2d_k &sfc_alb_dif, real1d_k &mu0, real1d_k &t_sfc, real2d_k &emis_sfc, real1d_k &lw_src, real2d_k &lwp, real2d_k &iwp, real2d_k &rel, real2d_k &rei, real2d_k &cldfrac, real3d_k &aer_tau_sw, real3d_k &aer_ssa_sw, real3d_k &aer_asm_sw, real3d_k &aer_tau_lw, real3d_k &cld_tau_sw_bnd, real3d_k &cld_tau_lw_bnd, real3d_k &cld_tau_sw_gpt, real3d_k &cld_tau_lw_gpt, 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_dn_dir, real3d_k &lw_bnd_flux_up, real3d_k &lw_bnd_flux_dn, const RealT tsi_scaling, const bool extra_clnclrsky_diag, const bool extra_clnsky_diag)
Definition: ERF_RRTMGP_Interface.cpp:384

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,
amrex::Real &  time,
const amrex::Real &  dt,
const amrex::BoxArray &  ba,
amrex::Geometry &  geom,
amrex::MultiFab *  cons_in,
amrex::MultiFab *  lsm_fluxes,
amrex::MultiFab *  lsm_zenith,
amrex::Vector< amrex::MultiFab * > &  lsm_input_ptrs,
amrex::MultiFab *  qheating_rates,
amrex::MultiFab *  z_phys,
amrex::MultiFab *  lat,
amrex::MultiFab *  lon 
)
132 {
133  // Set data members that may change
134  m_lev = level;
135  m_step = step;
136  m_time = time;
137  m_dt = dt;
138  m_geom = geom;
139  m_cons_in = cons_in;
140  m_lsm_fluxes = lsm_fluxes;
141  m_lsm_zenith = lsm_zenith;
142  m_qheating_rates = qheating_rates;
143  m_z_phys = z_phys;
144  m_lat = lat;
145  m_lon = lon;
146 
147  // Update the day and month
148  time_t timestamp = time_t(time);
149  struct tm *timeinfo = gmtime(&timestamp);
150  if (m_fixed_orbital_year) {
151  m_orbital_mon = timeinfo->tm_mon + 1;
152  m_orbital_day = timeinfo->tm_mday;
153  m_orbital_sec = timeinfo->tm_hour*3600 + timeinfo->tm_min*60 + timeinfo->tm_sec;
154  } else {
155  m_orbital_year = timeinfo->tm_year + 1900;
156  m_orbital_mon = timeinfo->tm_mon + 1;
157  m_orbital_day = timeinfo->tm_mday;
158  m_orbital_sec = timeinfo->tm_hour*3600 + timeinfo->tm_min*60 + timeinfo->tm_sec;
159  }
160 
161  // Only allocate and proceed if we are going to update radiation
162  m_update_rad = false;
163  if (m_rad_freq_in_steps > 0) { m_update_rad = ( (m_step == 0) || (m_step % m_rad_freq_in_steps == 0) ); }
164 
165  if (m_update_rad) {
166  // Call to Init() has set the dimensions: ncol & nlay
167 
168  // Allocate the buffer arrays
169  alloc_buffers();
170 
171  // Fill the KOKKOS Views from AMReX MFs
172  mf_to_kokkos_buffers(lsm_input_ptrs);
173 
174  // Initialize datalog MF on first step
175  if (m_first_step) {
176  m_first_step = false;
177  datalog_mf.define(cons_in->boxArray(), cons_in->DistributionMap(), 25, 0);
178  datalog_mf.setVal(0.0);
179  }
180  }
181 }
void mf_to_kokkos_buffers(amrex::Vector< amrex::MultiFab * > &lsm_input_ptrs)
Definition: ERF_Radiation.cpp:418
bool m_first_step
Definition: ERF_Radiation.H:292
void alloc_buffers()
Definition: ERF_Radiation.cpp:184
amrex::Real m_time
Definition: ERF_Radiation.H:275

Referenced by Run().

Here is the caller graph for this function:

◆ set_lsm_inputs()

template<typename LandSurfaceModelType >
void Radiation::set_lsm_inputs ( LandSurfaceModelType *  model)
inline
132  {
133  if constexpr(std::is_same_v<LandSurfaceModelType, SLM>)
134  {
135  if (!m_update_rad)
136  {
137  // skip getting LSM inputs if radiation not updating this timestep
138  return;
139  }
140 
141  // get surface temperature, LW flux, and emissivities from SLM
142  /*
143  auto slm_to_rad_vars = model->export_to_RRTMGP();
144 
145  for (amrex::MFIter mfi(*slm_to_rad_vars[0]); mfi.isValid(); ++mfi) {
146  const auto& vbx = mfi.validbox();
147  const auto& sbx = amrex::makeSlab(vbx, 2, 0);
148 
149  const int nx = vbx.length(0);
150  const int imin = vbx.smallEnd(0);
151  const int jmin = vbx.smallEnd(1);
152  const int offset = m_col_offsets[mfi.index()];
153 
154  const auto &slm_tsurf = slm_to_rad_vars[0]->const_array(mfi);
155  const auto &slm_alb_dir_vis_arr = slm_to_rad_vars[2]->const_array(mfi);
156  const auto &slm_alb_dir_nir_arr = slm_to_rad_vars[4]->const_array(mfi);
157  const auto &slm_IR_emis_arr = slm_to_rad_vars[5]->const_array(mfi);
158  const auto &slm_net_rad_arr = slm_to_rad_vars[6]->const_array(mfi);
159  amrex::ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
160  {
161  // map [i,j,k] 0-based to [icol, ilay] 1-based
162  const int icol = (j-jmin)*nx + (i-imin) + 1 + offset;
163  const int ilay = 0;
164 
165  t_sfc(icol) = slm_tsurf(i, j, k);
166 
167  sfc_alb_dir_vis(icol) = slm_alb_dir_vis_arr(i, j, k);
168  sfc_alb_dir_nir(icol) = slm_alb_dir_nir_arr(i, j, k);
169  sfc_alb_dif_vis(icol) = slm_alb_dir_vis_arr(i, j, k);
170  sfc_alb_dif_nir(icol) = slm_alb_dir_nir_arr(i, j, k);
171 
172  sfc_emis(icol) = slm_IR_emis_arr(i, j, k);
173 
174  //lw_src(icol) = slm_net_rad_arr(i, j, k, SLM_NetRad::net_lwup1);
175  lw_src(icol) = 0.0; // TODO
176  });
177  }
178 
179  // expand sfc_emis and lw_src over nlwbands
180  // TODO: check if this is correct - should emissivity vary based on band?
181  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {m_nlwbands, m_ncol}),
182  KOKKOS_LAMBDA (int ibnd, int icol)
183  {
184  emis_sfc(icol, ibnd) = sfc_emis(icol);
185  });
186  */
187 
188  Kokkos::deep_copy(emis_sfc, 0.98);
189  Kokkos::deep_copy(lw_src, 0.0);
190  }
191  }

◆ write_rrtmgp_fluxes()

void Radiation::write_rrtmgp_fluxes ( )
686 {
687  // Expose for device
688  auto sw_flux_up_d = sw_flux_up;
689  auto sw_flux_dn_d = sw_flux_dn;
690  auto sw_flux_dn_dir_d = sw_flux_dn_dir;
691  auto lw_flux_up_d = lw_flux_up;
692  auto lw_flux_dn_d = lw_flux_dn;
693 
694  int n_fluxes = 5;
695  MultiFab mf_flux(m_cons_in->boxArray(), m_cons_in->DistributionMap(), n_fluxes, 0);
696 
697  for (MFIter mfi(mf_flux); mfi.isValid(); ++mfi) {
698  const auto& vbx = mfi.validbox();
699  const int nx = vbx.length(0);
700  const int imin = vbx.smallEnd(0);
701  const int jmin = vbx.smallEnd(1);
702  const int offset = m_col_offsets[mfi.index()];
703  const Array4<Real>& dst_arr = mf_flux.array(mfi);
704  ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
705  {
706  // map [i,j,k] 0-based to [icol, ilay] 0-based
707  const int icol = (j-jmin)*nx + (i-imin) + offset;
708  const int ilay = k;
709 
710  // SW and LW fluxes
711  dst_arr(i,j,k,0) = sw_flux_up_d(icol,ilay);
712  dst_arr(i,j,k,1) = sw_flux_dn_d(icol,ilay);
713  dst_arr(i,j,k,2) = sw_flux_dn_dir_d(icol,ilay);
714  dst_arr(i,j,k,3) = lw_flux_up_d(icol,ilay);
715  dst_arr(i,j,k,4) = lw_flux_dn_d(icol,ilay);
716  });
717  }
718 
719 
720  std::string plotfilename = amrex::Concatenate("plt_rad", m_step, 5);
721  Vector<std::string> flux_names = {"sw_flux_up", "sw_flux_dn", "sw_flux_dir",
722  "lw_flux_up", "lw_flux_dn"};
723  WriteSingleLevelPlotfile(plotfilename, mf_flux, flux_names, m_geom, m_time, m_step);
724 }
Here is the call graph for this function:

◆ WriteDataLog()

void Radiation::WriteDataLog ( const amrex::Real &  time)
overridevirtual

Implements IRadiation.

816 {
817  constexpr int datwidth = 14;
818  constexpr int datprecision = 9;
819  constexpr int timeprecision = 13;
820 
821  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;
822  // Clear sky
823  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;
824  // Clean sky
825  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;
826  // Clean clear sky
827  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;
828 
829 
830  auto domain = m_geom.Domain();
831  h_avg_radqrsw = sumToLine(datalog_mf, 0, 1, domain, 2);
832  h_avg_radqrlw = sumToLine(datalog_mf, 1, 1, domain, 2);
833  h_avg_sw_up = sumToLine(datalog_mf, 2, 1, domain, 2);
834  h_avg_sw_dn = sumToLine(datalog_mf, 3, 1, domain, 2);
835  h_avg_sw_dn_dir = sumToLine(datalog_mf, 4, 1, domain, 2);
836  h_avg_lw_up = sumToLine(datalog_mf, 5, 1, domain, 2);
837  h_avg_lw_dn = sumToLine(datalog_mf, 6, 1, domain, 2);
838  h_avg_zenith = sumToLine(datalog_mf, 7, 1, domain, 2);
839 
840  h_avg_radqrcsw = sumToLine(datalog_mf, 8, 1, domain, 2);
841  h_avg_radqrclw = sumToLine(datalog_mf, 9, 1, domain, 2);
842  h_avg_sw_clr_up = sumToLine(datalog_mf, 10, 1, domain, 2);
843  h_avg_sw_clr_dn = sumToLine(datalog_mf, 11, 1, domain, 2);
844  h_avg_sw_clr_dn_dir = sumToLine(datalog_mf, 12, 1, domain, 2);
845  h_avg_lw_clr_up = sumToLine(datalog_mf, 13, 1, domain, 2);
846  h_avg_lw_clr_dn = sumToLine(datalog_mf, 14, 1, domain, 2);
847 
848  if (m_extra_clnsky_diag) {
849  h_avg_sw_cln_up = sumToLine(datalog_mf, 15, 1, domain, 2);
850  h_avg_sw_cln_dn = sumToLine(datalog_mf, 16, 1, domain, 2);
851  h_avg_sw_cln_dn_dir = sumToLine(datalog_mf, 17, 1, domain, 2);
852  h_avg_lw_cln_up = sumToLine(datalog_mf, 18, 1, domain, 2);
853  h_avg_lw_cln_dn = sumToLine(datalog_mf, 19, 1, domain, 2);
854  }
855 
857  h_avg_sw_clnclr_up = sumToLine(datalog_mf, 20, 1, domain, 2);
858  h_avg_sw_clnclr_dn = sumToLine(datalog_mf, 21, 1, domain, 2);
859  h_avg_sw_clnclr_dn_dir = sumToLine(datalog_mf, 22, 1, domain, 2);
860  h_avg_lw_clnclr_up = sumToLine(datalog_mf, 23, 1, domain, 2);
861  h_avg_lw_clnclr_dn = sumToLine(datalog_mf, 24, 1, domain, 2);
862  }
863 
864  Real area_z = static_cast<Real>(domain.length(0)*domain.length(1));
865  int nz = domain.length(2);
866  for (int k = 0; k < nz; k++) {
867  h_avg_radqrsw[k] /= area_z;
868  h_avg_radqrlw[k] /= area_z;
869  h_avg_sw_up[k] /= area_z;
870  h_avg_sw_dn[k] /= area_z;
871  h_avg_sw_dn_dir[k] /= area_z;
872  h_avg_lw_up[k] /= area_z;
873  h_avg_lw_dn[k] /= area_z;
874  h_avg_zenith[k] /= area_z;
875 
876  h_avg_radqrcsw[k] /= area_z;
877  h_avg_radqrclw[k] /= area_z;
878  h_avg_sw_clr_up[k] /= area_z;
879  h_avg_sw_clr_dn[k] /= area_z;
880  h_avg_sw_clr_dn_dir[k] /= area_z;
881  h_avg_lw_clr_up[k] /= area_z;
882  h_avg_lw_clr_dn[k] /= area_z;
883  }
884 
885  if (m_extra_clnsky_diag) {
886  for (int k = 0; k < nz; k++) {
887  h_avg_sw_cln_up[k] /= area_z;
888  h_avg_sw_cln_dn[k] /= area_z;
889  h_avg_sw_cln_dn_dir[k] /= area_z;
890  h_avg_lw_cln_up[k] /= area_z;
891  h_avg_lw_cln_dn[k] /= area_z;
892  }
893  }
894 
896  for (int k = 0; k < nz; k++) {
897  h_avg_sw_clnclr_up[k] /= area_z;
898  h_avg_sw_clnclr_dn[k] /= area_z;
899  h_avg_sw_clnclr_dn_dir[k] /= area_z;
900  h_avg_lw_clnclr_up[k] /= area_z;
901  h_avg_lw_clnclr_dn[k] /= area_z;
902  }
903  }
904 
905  if (ParallelDescriptor::IOProcessor()) {
906  std::ostream& log = *datalog;
907  if (log.good()) {
908 
909  for (int k = 0; k < nz; k++)
910  {
911  Real z = k * m_geom.CellSize(2);
912  log << std::setw(datwidth) << std::setprecision(timeprecision) << time << " "
913  << std::setw(datwidth) << std::setprecision(datprecision) << z << " "
914  << h_avg_radqrsw[k] << " " << h_avg_radqrlw[k] << " " << h_avg_sw_up[k] << " "
915  << h_avg_sw_dn[k] << " " << h_avg_sw_dn_dir[k] << " " << h_avg_lw_up[k] << " "
916  << h_avg_lw_dn[k] << " " << h_avg_zenith[k] << " "
917  << h_avg_radqrcsw[k] << " " << h_avg_radqrclw[k] << " " << h_avg_sw_clr_up[k] << " "
918  << h_avg_sw_clr_dn[k] << " " << h_avg_sw_clr_dn_dir[k] << " " << h_avg_lw_clr_up[k] << " "
919  << h_avg_lw_clr_dn[k] << " ";
920  if (m_extra_clnsky_diag) {
921  log << h_avg_sw_cln_up[k] << " " << h_avg_sw_cln_dn[k] << " " << h_avg_sw_cln_dn_dir[k] << " "
922  << h_avg_lw_cln_up[k] << " " << h_avg_lw_cln_dn[k] << " ";
923  } else {
924  log << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " ";
925  }
926 
928  log << h_avg_sw_clnclr_up[k] << " " << h_avg_sw_clnclr_dn[k] << " " << h_avg_sw_clnclr_dn_dir[k] << " "
929  << h_avg_lw_clnclr_up[k] << " " << h_avg_lw_clnclr_dn[k] << std::endl;
930  } else {
931  log << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << std::endl;
932  }
933  }
934  // Write top face values
935  Real z = nz * m_geom.CellSize(2);
936  log << std::setw(datwidth) << std::setprecision(timeprecision) << time << " "
937  << std::setw(datwidth) << std::setprecision(datprecision) << z << " "
938  << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " "
939  << 0.0 << " " << 0.0 << " "
940  << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " "
941  << 0.0 << " "
942  << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " "
943  << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0
944  << std::endl;
945  }
946  }
947 }
std::unique_ptr< std::fstream > datalog
Definition: ERF_RadiationInterface.H:84

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

◆ cld_tau_lw_bnd

real3d_k Radiation::cld_tau_lw_bnd
private

◆ cld_tau_lw_gpt

real3d_k Radiation::cld_tau_lw_gpt
private

◆ cld_tau_sw_bnd

real3d_k Radiation::cld_tau_sw_bnd
private

◆ cld_tau_sw_gpt

real3d_k Radiation::cld_tau_sw_gpt
private

◆ cldfrac_tot

real2d_k Radiation::cldfrac_tot
private

◆ cosine_zenith

real1d_k Radiation::cosine_zenith
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

◆ emis_sfc

real2d_k Radiation::emis_sfc
private

Referenced by set_lsm_inputs().

◆ gas_names_offset

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

◆ iwp

real2d_k Radiation::iwp
private

◆ lat

real1d_k Radiation::lat
private

Referenced by Run().

◆ lon

real1d_k Radiation::lon
private

Referenced by Run().

◆ 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

Referenced by set_lsm_inputs().

◆ lwp

real2d_k Radiation::lwp
private

◆ m_ba

amrex::BoxArray Radiation::m_ba
private

◆ m_ch4vmr

amrex::Real Radiation::m_ch4vmr = 1807.851e-9
private

◆ m_co2vmr

amrex::Real Radiation::m_co2vmr = 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 = 1.0e-7
private

◆ m_do_aerosol_rad

bool Radiation::m_do_aerosol_rad = true
private

◆ m_do_subcol_sampling

bool Radiation::m_do_subcol_sampling = true
private

◆ m_dt

amrex::Real 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_first_step

bool Radiation::m_first_step = true
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 = -9999.
private

◆ m_fixed_total_solar_irradiance

amrex::Real Radiation::m_fixed_total_solar_irradiance = -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 = 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 = -98.555183
private

◆ m_lsm

bool Radiation::m_lsm = false
private

◆ m_lsm_fluxes

amrex::MultiFab* Radiation::m_lsm_fluxes = nullptr
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 = {"sw_flux_dn", "lw_flux_dn"}
private

Referenced by get_lsm_output_varnames().

◆ m_lsm_zenith

amrex::MultiFab* Radiation::m_lsm_zenith = nullptr
private

◆ 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:
= {18.01528, 44.00950, 47.9982, 44.0128,
28.01010, 16.04246, 31.9980, 28.0134}

◆ m_n2ovmr

amrex::Real Radiation::m_n2ovmr = 323.141e-9
private

◆ m_n2vmr

amrex::Real Radiation::m_n2vmr = 0.7906
private

◆ m_ncol

int Radiation::m_ncol
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 = 16
private

◆ m_nlwgpts

int Radiation::m_nlwgpts = 128
private

◆ m_nswbands

int Radiation::m_nswbands = 14
private

◆ m_nswgpts

int Radiation::m_nswgpts = 112
private

◆ m_o2vmr

amrex::Real Radiation::m_o2vmr = 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

◆ m_orbital_eccen

amrex::Real Radiation::m_orbital_eccen = -9999.
private

◆ m_orbital_mon

int Radiation::m_orbital_mon = -9999
private

◆ m_orbital_mvelp

amrex::Real Radiation::m_orbital_mvelp = -9999.
private

◆ m_orbital_obliq

amrex::Real Radiation::m_orbital_obliq = -9999.
private

◆ m_orbital_sec

int Radiation::m_orbital_sec = -9999
private

◆ m_orbital_year

int Radiation::m_orbital_year = -9999
private

◆ m_qheating_rates

amrex::MultiFab* Radiation::m_qheating_rates = nullptr
private

◆ m_rad_freq_in_steps

int Radiation::m_rad_freq_in_steps = 1
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

amrex::Real Radiation::m_time
private

◆ m_update_rad

bool Radiation::m_update_rad = false
private

Referenced by rad_run_impl(), and set_lsm_inputs().

◆ 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_del

real2d_k Radiation::p_del
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-g224-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

◆ tmp2d

real2d_k Radiation::tmp2d
private

◆ z_del

real2d_k Radiation::z_del
private

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