ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
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)
 
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 = {"cos_zenith_angle", "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-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::Realm_mol_weight_gas
 
amrex::Real m_co2vmr = 388.717e-6
 
amrex::Vector< amrex::Realm_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
 
int m_nlwbands
 
int m_nswgpts
 
int m_nlwgpts
 
int m_rad_freq_in_steps = 1
 
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 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 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 
)
22 {
23  // Note that Kokkos is now initialized in main.cpp
24 
25  // Check if we have a valid moisture model
26  if (sc.moisture_type != MoistureType::None) { m_moist = true; }
27 
28  // Check if we have a moisture model with ice
29  if (sc.moisture_type == MoistureType::SAM) { m_ice = true; }
30 
31  // Check if we have a land surface model enabled
32  if (sc.lsm_type != LandSurfaceType::None) { m_lsm = true; }
33 
34  // Construct parser object for following reads
35  ParmParse pp("erf");
36 
37  // Must specify a surface temp without a LSM
38  if (!m_lsm) { pp.get("rad_t_sfc",m_rad_t_sfc); }
39 
40  // Radiation timestep, as a number of atm steps
41  pp.query("rad_freq_in_steps", m_rad_freq_in_steps);
42 
43  // Flag to write fluxes to plt file
44  pp.query("rad_write_fluxes", m_rad_write_fluxes);
45 
46  // Do MCICA subcolumn sampling
47  pp.query("rad_do_subcol_sampling", m_do_subcol_sampling);
48 
49  // Determine orbital year. If orbital_year is negative, use current year
50  // from timestamp for orbital year; if positive, use provided orbital year
51  // for duration of simulation.
52  m_fixed_orbital_year = pp.query("rad_orbital_year", m_orbital_year);
53 
54  // Get orbital parameters from inputs file
55  pp.query("rad_orbital_eccentricity", m_orbital_eccen);
56  pp.query("rad_orbital_obliquity" , m_orbital_obliq);
57  pp.query("rad_orbital_mvelp" , m_orbital_mvelp);
58 
59  // Get a constant lat/lon for idealized simulations
60  pp.query("rad_cons_lat", m_lat_cons);
61  pp.query("rad_cons_lon", m_lon_cons);
62 
63  // Value for prescribing an invariant solar constant (i.e. total solar irradiance at
64  // TOA). Used for idealized experiments such as RCE. Disabled when value is less than 0.
65  pp.query("fixed_total_solar_irradiance", m_fixed_total_solar_irradiance);
66 
67  // Determine whether or not we are using a fixed solar zenith angle (positive value)
68  pp.query("fixed_solar_zenith_angle", m_fixed_solar_zenith_angle);
69 
70  // Get prescribed surface values of greenhouse gases
71  pp.query("co2vmr", m_co2vmr);
72  pp.queryarr("o3vmr" , m_o3vmr );
73  pp.query("n2ovmr", m_n2ovmr);
74  pp.query("covmr" , m_covmr );
75  pp.query("ch4vmr", m_ch4vmr);
76  pp.query("o2vmr" , m_o2vmr );
77  pp.query("n2vmr" , m_n2vmr );
78 
79  // Required aerosol optical properties from SPA
80  pp.query("rad_do_aerosol", m_do_aerosol_rad);
81 
82  // Whether we do extra clean/clear sky calculations
83  pp.query("rad_extra_clnclrsky_diag", m_extra_clnclrsky_diag);
84  pp.query("rad_extra_clnsky_diag" , m_extra_clnsky_diag);
85 
86  // Parse path and file names
87  pp.query("rrtmgp_file_path" , rrtmgp_file_path);
88  pp.query("rrtmgp_coeffs_sw" , rrtmgp_coeffs_sw );
89  pp.query("rrtmgp_coeffs_lw" , rrtmgp_coeffs_lw );
90  pp.query("rrtmgp_cloud_optics_sw", rrtmgp_cloud_optics_sw);
91  pp.query("rrtmgp_cloud_optics_lw", rrtmgp_cloud_optics_lw);
92 
93  // Append file names to path
98 
99  // Get dimensions from lookup data
100  if (ParallelDescriptor::IOProcessor()) {
101  auto ncf_sw = ncutils::NCFile::open(rrtmgp_coeffs_file_sw, NC_CLOBBER | NC_NETCDF4);
102  m_nswbands = ncf_sw.dim("bnd").len();
103  m_nswgpts = ncf_sw.dim("gpt").len();
104  ncf_sw.close();
105 
106  auto ncf_lw = ncutils::NCFile::open(rrtmgp_coeffs_file_lw, NC_CLOBBER | NC_NETCDF4);
107  m_nlwbands = ncf_lw.dim("bnd").len();
108  m_nlwgpts = ncf_lw.dim("gpt").len();
109  ncf_lw.close();
110  }
111  int ioproc = ParallelDescriptor::IOProcessorNumber(); // I/O rank
112  ParallelDescriptor::Bcast(&m_nswbands, 1, ioproc);
113  ParallelDescriptor::Bcast(&m_nlwbands, 1, ioproc);
114  ParallelDescriptor::Bcast(&m_nswgpts, 1, ioproc);
115  ParallelDescriptor::Bcast(&m_nlwgpts, 1, ioproc);
116 
117  // Output for user
118  if (lev == 0) {
119  Print() << "Radiation interface constructed:\n";
120  Print() << "========================================================\n";
121  Print() << "Coeff SW file: " << rrtmgp_coeffs_file_sw << "\n";
122  Print() << "Coeff LW file: " << rrtmgp_coeffs_file_lw << "\n";
123  Print() << "Cloud SW file: " << rrtmgp_cloud_optics_file_sw << "\n";
124  Print() << "Cloud LW file: " << rrtmgp_cloud_optics_file_lw << "\n";
125  Print() << "Number of short/longwave bands: "
126  << m_nswbands << " " << m_nlwbands << "\n";
127  Print() << "Number of short/longwave gauss points: "
128  << m_nswgpts << " " << m_nlwgpts << "\n";
129  Print() << "========================================================\n";
130  }
131 }
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:280
std::string rrtmgp_coeffs_sw
Definition: ERF_Radiation.H:276
int m_rad_freq_in_steps
Definition: ERF_Radiation.H:358
amrex::Real m_lon_cons
Definition: ERF_Radiation.H:265
int m_nswbands
Definition: ERF_Radiation.H:352
bool m_do_aerosol_rad
Definition: ERF_Radiation.H:320
std::string rrtmgp_coeffs_lw
Definition: ERF_Radiation.H:277
bool m_moist
Definition: ERF_Radiation.H:233
std::string rrtmgp_cloud_optics_file_lw
Definition: ERF_Radiation.H:283
bool m_lsm
Definition: ERF_Radiation.H:237
bool m_do_subcol_sampling
Definition: ERF_Radiation.H:361
bool m_rad_write_fluxes
Definition: ERF_Radiation.H:228
std::string rrtmgp_cloud_optics_file_sw
Definition: ERF_Radiation.H:282
amrex::Real m_orbital_mvelp
Definition: ERF_Radiation.H:340
amrex::Vector< amrex::Real > m_o3vmr
Definition: ERF_Radiation.H:294
std::string rrtmgp_cloud_optics_sw
Definition: ERF_Radiation.H:278
int m_nlwgpts
Definition: ERF_Radiation.H:355
std::string rrtmgp_cloud_optics_lw
Definition: ERF_Radiation.H:279
amrex::Real m_co2vmr
Definition: ERF_Radiation.H:293
bool m_extra_clnsky_diag
Definition: ERF_Radiation.H:323
amrex::Real m_lat_cons
Definition: ERF_Radiation.H:264
amrex::Real m_fixed_total_solar_irradiance
Definition: ERF_Radiation.H:345
amrex::Real m_o2vmr
Definition: ERF_Radiation.H:298
amrex::Real m_n2ovmr
Definition: ERF_Radiation.H:295
amrex::Real m_rad_t_sfc
Definition: ERF_Radiation.H:248
std::string rrtmgp_file_path
Definition: ERF_Radiation.H:275
int m_nlwbands
Definition: ERF_Radiation.H:353
amrex::Real m_orbital_eccen
Definition: ERF_Radiation.H:338
amrex::Real m_covmr
Definition: ERF_Radiation.H:296
std::string rrtmgp_coeffs_file_lw
Definition: ERF_Radiation.H:281
amrex::Real m_ch4vmr
Definition: ERF_Radiation.H:297
bool m_ice
Definition: ERF_Radiation.H:234
amrex::Real m_fixed_solar_zenith_angle
Definition: ERF_Radiation.H:349
amrex::Real m_orbital_obliq
Definition: ERF_Radiation.H:339
amrex::Real m_n2vmr
Definition: ERF_Radiation.H:299
bool m_fixed_orbital_year
Definition: ERF_Radiation.H:337
int m_nswgpts
Definition: ERF_Radiation.H:354
bool m_extra_clnclrsky_diag
Definition: ERF_Radiation.H:324
int m_orbital_year
Definition: ERF_Radiation.H:329
static NCFile open(const std::string &name, const int cmode=NC_NOWRITE)
Definition: ERF_NCInterface.cpp:707
LandSurfaceType lsm_type
Definition: ERF_DataStruct.H:914
MoistureType moisture_type
Definition: ERF_DataStruct.H:911
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 ( )
202 {
203  // 1d size (m_ngas)
204  const Real* mol_weight_gas_p = m_mol_weight_gas.data();
205  const std::string* gas_names_p = m_gas_names.data();
206  m_gas_mol_weights = real1d_k("m_gas_mol_weights", m_ngas);
207  realHost1d_k m_gas_mol_weights_h("m_gas_mol_weights_h", m_ngas);
208  gas_names_offset.clear(); gas_names_offset.resize(m_ngas);
209  std::string* gas_names_offset_p = gas_names_offset.data();
210  Kokkos::parallel_for(Kokkos::RangePolicy<Kokkos::Serial>(0, m_ngas),
211  [&] (int igas)
212  {
213  m_gas_mol_weights_h(igas) = mol_weight_gas_p[igas];
214  gas_names_offset_p[igas] = gas_names_p[igas];
215  });
216  Kokkos::deep_copy(m_gas_mol_weights, m_gas_mol_weights_h);
217 
218  // 1d size (1 or nlay)
219  m_o3_size = m_o3vmr.size();
220  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(((m_o3_size==1) || (m_o3_size==m_nlay)),
221  "O3 VMR array must be length 1 or nlay");
222  Real* o3vmr_p = m_o3vmr.data();
223  o3_lay = real1d_k("o3_lay", m_o3_size);
224  realHost1d_k o3_lay_h("o3_lay_h", m_o3_size);
225  Kokkos::parallel_for(Kokkos::RangePolicy<Kokkos::Serial>(0, m_o3_size),
226  [&] (int io3)
227  {
228  o3_lay_h(io3) = o3vmr_p[io3];
229  });
230  Kokkos::deep_copy(o3_lay, o3_lay_h);
231 
232  // 1d size (ncol)
233  mu0 = real1d_k("mu0" , m_ncol);
234  sfc_alb_dir_vis = real1d_k("sfc_alb_dir_vis" , m_ncol);
235  sfc_alb_dir_nir = real1d_k("sfc_alb_dir_nir" , m_ncol);
236  sfc_alb_dif_vis = real1d_k("sfc_alb_dif_vis" , m_ncol);
237  sfc_alb_dif_nir = real1d_k("sfc_alb_dif_nir" , m_ncol);
238  sfc_flux_dir_vis = real1d_k("sfc_flux_dir_vis", m_ncol);
239  sfc_flux_dir_nir = real1d_k("sfc_flux_dir_nir", m_ncol);
240  sfc_flux_dif_vis = real1d_k("sfc_flux_dif_vis", m_ncol);
241  sfc_flux_dif_nir = real1d_k("sfc_flux_dif_nir", m_ncol);
242  lat = real1d_k("lat" , m_ncol);
243  lon = real1d_k("lon" , m_ncol);
244  sfc_emis = real1d_k("sfc_emis" , m_ncol);
245  t_sfc = real1d_k("t_sfc" , m_ncol);
246  lw_src = real1d_k("lw_src" , m_ncol);
247 
248  // 2d size (ncol, nlay)
249  r_lay = real2d_k("r_lay" , m_ncol, m_nlay);
250  p_lay = real2d_k("p_lay" , m_ncol, m_nlay);
251  t_lay = real2d_k("t_lay" , m_ncol, m_nlay);
252  z_del = real2d_k("z_del" , m_ncol, m_nlay);
253  p_del = real2d_k("p_del" , m_ncol, m_nlay);
254  qv_lay = real2d_k("qv" , m_ncol, m_nlay);
255  qc_lay = real2d_k("qc" , m_ncol, m_nlay);
256  qi_lay = real2d_k("qi" , m_ncol, m_nlay);
257  cldfrac_tot = real2d_k("cldfrac_tot" , m_ncol, m_nlay);
258  eff_radius_qc = real2d_k("eff_radius_qc", m_ncol, m_nlay);
259  eff_radius_qi = real2d_k("eff_radius_qi", m_ncol, m_nlay);
260  lwp = real2d_k("lwp" , m_ncol, m_nlay);
261  iwp = real2d_k("iwp" , m_ncol, m_nlay);
262  sw_heating = real2d_k("sw_heating" , m_ncol, m_nlay);
263  lw_heating = real2d_k("lw_heating" , m_ncol, m_nlay);
264  sw_clrsky_heating = real2d_k("sw_clrsky_heating", m_ncol, m_nlay);
265  lw_clrsky_heating = real2d_k("lw_clrsky_heating", m_ncol, m_nlay);
266 
267  // 2d size (ncol, nlay+1)
268  d_tint = real2d_k("d_tint" , m_ncol, m_nlay+1);
269  p_lev = real2d_k("p_lev" , m_ncol, m_nlay+1);
270  t_lev = real2d_k("t_lev" , m_ncol, m_nlay+1);
271  sw_flux_up = real2d_k("sw_flux_up" , m_ncol, m_nlay+1);
272  sw_flux_dn = real2d_k("sw_flux_dn" , m_ncol, m_nlay+1);
273  sw_flux_dn_dir = real2d_k("sw_flux_dn_dir" , m_ncol, m_nlay+1);
274  lw_flux_up = real2d_k("sw_flux_up" , m_ncol, m_nlay+1);
275  lw_flux_dn = real2d_k("sw_flux_dn" , m_ncol, m_nlay+1);
276  sw_clnclrsky_flux_up = real2d_k("sw_clnclrsky_flux_up" , m_ncol, m_nlay+1);
277  sw_clnclrsky_flux_dn = real2d_k("sw_clnclrsky_flux_dn" , m_ncol, m_nlay+1);
278  sw_clnclrsky_flux_dn_dir = real2d_k("sw_clnclrsky_flux_dn_dir", m_ncol, m_nlay+1);
279  sw_clrsky_flux_up = real2d_k("sw_clrsky_flux_up" , m_ncol, m_nlay+1);
280  sw_clrsky_flux_dn = real2d_k("sw_clrsky_flux_dn" , m_ncol, m_nlay+1);
281  sw_clrsky_flux_dn_dir = real2d_k("sw_clrsky_flux_dn_dir" , m_ncol, m_nlay+1);
282  sw_clnsky_flux_up = real2d_k("sw_clnsky_flux_up" , m_ncol, m_nlay+1);
283  sw_clnsky_flux_dn = real2d_k("sw_clnsky_flux_dn" , m_ncol, m_nlay+1);
284  sw_clnsky_flux_dn_dir = real2d_k("sw_clnsky_flux_dn_dir" , m_ncol, m_nlay+1);
285  lw_clnclrsky_flux_up = real2d_k("lw_clnclrsky_flux_up" , m_ncol, m_nlay+1);
286  lw_clnclrsky_flux_dn = real2d_k("lw_clnclrsky_flux_dn" , m_ncol, m_nlay+1);
287  lw_clrsky_flux_up = real2d_k("lw_clrsky_flux_up" , m_ncol, m_nlay+1);
288  lw_clrsky_flux_dn = real2d_k("lw_clrsky_flux_dn" , m_ncol, m_nlay+1);
289  lw_clnsky_flux_up = real2d_k("lw_clnsky_flux_up" , m_ncol, m_nlay+1);
290  lw_clnsky_flux_dn = real2d_k("lw_clnsky_flux_dn" , m_ncol, m_nlay+1);
291 
292  // 3d size (ncol, nlay+1, nswbands)
293  sw_bnd_flux_up = real3d_k("sw_bnd_flux_up" , m_ncol, m_nlay+1, m_nswbands);
294  sw_bnd_flux_dn = real3d_k("sw_bnd_flux_dn" , m_ncol, m_nlay+1, m_nswbands);
295  sw_bnd_flux_dir = real3d_k("sw_bnd_flux_dir", m_ncol, m_nlay+1, m_nswbands);
296  sw_bnd_flux_dif = real3d_k("sw_bnd_flux_dif", m_ncol, m_nlay+1, m_nswbands);
297 
298  // 3d size (ncol, nlay+1, nlwbands)
299  lw_bnd_flux_up = real3d_k("lw_bnd_flux_up" , m_ncol, m_nlay+1, m_nlwbands);
300  lw_bnd_flux_dn = real3d_k("lw_bnd_flux_up" , m_ncol, m_nlay+1, m_nlwbands);
301 
302  // 2d size (ncol, nswbands)
303  sfc_alb_dir = real2d_k("sfc_alb_dir", m_ncol, m_nswbands);
304  sfc_alb_dif = real2d_k("sfc_alb_dif", m_ncol, m_nswbands);
305 
306  // 2d size (ncol, nlwbands)
307  emis_sfc = real2d_k("emis_sfc", m_ncol, m_nlwbands);
308 
309  // 3d size (ncol, nlay, n[sw,lw]bands)
310  aero_tau_sw = real3d_k("aero_tau_sw", m_ncol, m_nlay, m_nswbands);
311  aero_ssa_sw = real3d_k("aero_ssa_sw", m_ncol, m_nlay, m_nswbands);
312  aero_g_sw = real3d_k("aero_g_sw" , m_ncol, m_nlay, m_nswbands);
313  aero_tau_lw = real3d_k("aero_tau_lw", m_ncol, m_nlay, m_nlwbands);
314 
315  // 3d size (ncol, nlay, n[sw,lw]bnds)
316  cld_tau_sw_bnd = real3d_k("cld_tau_sw_bnd", m_ncol, m_nlay, m_nswbands);
317  cld_tau_lw_bnd = real3d_k("cld_tau_lw_bnd", m_ncol, m_nlay, m_nlwbands);
318 
319  // 3d size (ncol, nlay, n[sw,lw]gpts)
320  cld_tau_sw_gpt = real3d_k("cld_tau_sw_gpt", m_ncol, m_nlay, m_nswgpts);
321  cld_tau_lw_gpt = real3d_k("cld_tau_lw_gpt", m_ncol, m_nlay, m_nlwgpts);
322 }
Kokkos::View< RealT *, KokkosDefaultDevice > real1d_k
Definition: ERF_Kokkos.H:18
Kokkos::View< RealT ***, layout_t, KokkosDefaultDevice > real3d_k
Definition: ERF_Kokkos.H:20
Kokkos::View< RealT **, layout_t, KokkosDefaultDevice > real2d_k
Definition: ERF_Kokkos.H:19
Kokkos::View< RealT *, KokkosHostDevice > realHost1d_k
Definition: ERF_Kokkos.H:16
amrex::Real Real
Definition: ERF_ShocInterface.H:19
real3d_k sw_bnd_flux_dn
Definition: ERF_Radiation.H:428
real2d_k lw_flux_up
Definition: ERF_Radiation.H:408
real3d_k aero_tau_sw
Definition: ERF_Radiation.H:444
int m_o3_size
Definition: ERF_Radiation.H:303
real3d_k sw_bnd_flux_dir
Definition: ERF_Radiation.H:429
real2d_k sw_clnsky_flux_dn
Definition: ERF_Radiation.H:417
real2d_k d_tint
Definition: ERF_Radiation.H:402
real2d_k lw_clnclrsky_flux_dn
Definition: ERF_Radiation.H:420
real2d_k lwp
Definition: ERF_Radiation.H:394
real3d_k cld_tau_lw_gpt
Definition: ERF_Radiation.H:455
real1d_k lw_src
Definition: ERF_Radiation.H:380
real2d_k eff_radius_qi
Definition: ERF_Radiation.H:393
real1d_k m_gas_mol_weights
Definition: ERF_Radiation.H:304
real2d_k sw_heating
Definition: ERF_Radiation.H:396
real3d_k lw_bnd_flux_dn
Definition: ERF_Radiation.H:434
real3d_k sw_bnd_flux_up
Definition: ERF_Radiation.H:427
real2d_k qv_lay
Definition: ERF_Radiation.H:388
real2d_k sw_clnclrsky_flux_dn_dir
Definition: ERF_Radiation.H:412
real1d_k sfc_flux_dif_vis
Definition: ERF_Radiation.H:374
real2d_k lw_clnclrsky_flux_up
Definition: ERF_Radiation.H:419
real3d_k cld_tau_sw_gpt
Definition: ERF_Radiation.H:454
real1d_k lat
Definition: ERF_Radiation.H:376
real2d_k qi_lay
Definition: ERF_Radiation.H:390
real2d_k sw_clrsky_flux_up
Definition: ERF_Radiation.H:413
real2d_k t_lev
Definition: ERF_Radiation.H:404
real1d_k o3_lay
Definition: ERF_Radiation.H:364
real1d_k sfc_alb_dif_vis
Definition: ERF_Radiation.H:370
real1d_k mu0
Definition: ERF_Radiation.H:367
real2d_k cldfrac_tot
Definition: ERF_Radiation.H:391
real3d_k cld_tau_lw_bnd
Definition: ERF_Radiation.H:451
real1d_k sfc_flux_dir_nir
Definition: ERF_Radiation.H:373
real2d_k sw_clnclrsky_flux_up
Definition: ERF_Radiation.H:410
real3d_k aero_g_sw
Definition: ERF_Radiation.H:446
real2d_k sw_flux_up
Definition: ERF_Radiation.H:405
real3d_k sw_bnd_flux_dif
Definition: ERF_Radiation.H:430
real2d_k sfc_alb_dif
Definition: ERF_Radiation.H:438
real2d_k emis_sfc
Definition: ERF_Radiation.H:441
real1d_k sfc_alb_dif_nir
Definition: ERF_Radiation.H:371
real2d_k lw_clnsky_flux_dn
Definition: ERF_Radiation.H:424
real2d_k p_lay
Definition: ERF_Radiation.H:384
real2d_k r_lay
Definition: ERF_Radiation.H:383
int m_ncol
Definition: ERF_Radiation.H:313
real2d_k sw_flux_dn_dir
Definition: ERF_Radiation.H:407
const std::vector< amrex::Real > m_mol_weight_gas
Definition: ERF_Radiation.H:289
std::vector< std::string > gas_names_offset
Definition: ERF_Radiation.H:305
real2d_k sw_clrsky_flux_dn
Definition: ERF_Radiation.H:414
real3d_k aero_ssa_sw
Definition: ERF_Radiation.H:445
real2d_k sw_clnsky_flux_dn_dir
Definition: ERF_Radiation.H:418
real2d_k qc_lay
Definition: ERF_Radiation.H:389
real2d_k lw_clrsky_flux_up
Definition: ERF_Radiation.H:421
real2d_k sw_flux_dn
Definition: ERF_Radiation.H:406
real2d_k lw_clrsky_flux_dn
Definition: ERF_Radiation.H:422
real2d_k z_del
Definition: ERF_Radiation.H:386
real2d_k p_del
Definition: ERF_Radiation.H:387
real3d_k aero_tau_lw
Definition: ERF_Radiation.H:447
real2d_k lw_flux_dn
Definition: ERF_Radiation.H:409
real2d_k sw_clrsky_flux_dn_dir
Definition: ERF_Radiation.H:415
real1d_k sfc_flux_dif_nir
Definition: ERF_Radiation.H:375
int m_ngas
Definition: ERF_Radiation.H:286
real2d_k lw_heating
Definition: ERF_Radiation.H:397
real1d_k sfc_alb_dir_nir
Definition: ERF_Radiation.H:369
real2d_k lw_clnsky_flux_up
Definition: ERF_Radiation.H:423
real2d_k sfc_alb_dir
Definition: ERF_Radiation.H:437
real1d_k lon
Definition: ERF_Radiation.H:377
real1d_k sfc_emis
Definition: ERF_Radiation.H:378
real3d_k cld_tau_sw_bnd
Definition: ERF_Radiation.H:450
real2d_k sw_clnsky_flux_up
Definition: ERF_Radiation.H:416
int m_nlay
Definition: ERF_Radiation.H:314
real3d_k lw_bnd_flux_up
Definition: ERF_Radiation.H:433
real1d_k sfc_alb_dir_vis
Definition: ERF_Radiation.H:368
const std::vector< std::string > m_gas_names
Definition: ERF_Radiation.H:287
real2d_k sw_clrsky_heating
Definition: ERF_Radiation.H:398
real2d_k t_lay
Definition: ERF_Radiation.H:385
real1d_k sfc_flux_dir_vis
Definition: ERF_Radiation.H:372
real2d_k iwp
Definition: ERF_Radiation.H:395
real2d_k p_lev
Definition: ERF_Radiation.H:403
real2d_k lw_clrsky_heating
Definition: ERF_Radiation.H:399
real2d_k eff_radius_qc
Definition: ERF_Radiation.H:392
real1d_k t_sfc
Definition: ERF_Radiation.H:379
real2d_k sw_clnclrsky_flux_dn
Definition: ERF_Radiation.H:411

◆ dealloc_buffers()

void Radiation::dealloc_buffers ( )
326 {
327  // 1d size (m_ngas)
329 
330  // 1d size (1 or nlay)
331  o3_lay = real1d_k();
332 
333  // 1d size (ncol)
334  mu0 = real1d_k();
343  lat = real1d_k();
344  lon = real1d_k();
345  sfc_emis = real1d_k();
346  t_sfc = real1d_k();
347  lw_src = real1d_k();
348 
349  // 2d size (ncol, nlay)
350  r_lay = real2d_k();
351  p_lay = real2d_k();
352  t_lay = real2d_k();
353  z_del = real2d_k();
354  p_del = real2d_k();
355  qv_lay = real2d_k();
356  qc_lay = real2d_k();
357  qi_lay = real2d_k();
358  cldfrac_tot = real2d_k();
361  lwp = real2d_k();
362  iwp = real2d_k();
363  sw_heating = real2d_k();
364  lw_heating = real2d_k();
367  sw_heating = real2d_k();
368  lw_heating = real2d_k();
371 
372  // 2d size (ncol, nlay+1)
373  d_tint = real2d_k();
374  p_lev = real2d_k();
375  t_lev = real2d_k();
376  sw_flux_up = real2d_k();
377  sw_flux_dn = real2d_k();
379  lw_flux_up = real2d_k();
380  lw_flux_dn = real2d_k();
396 
397  // 3d size (ncol, nlay+1, nswbands)
402 
403  // 3d size (ncol, nlay+1, nlwbands)
406 
407  // 2d size (ncol, nswbands)
408  sfc_alb_dir = real2d_k();
409  sfc_alb_dif = real2d_k();
410 
411  // 2d size (ncol, nlwbands)
412  emis_sfc = real2d_k();
413 
414  // 3d size (ncol, nlay, n[sw,lw]bands)
415  aero_tau_sw = real3d_k();
416  aero_ssa_sw = real3d_k();
417  aero_g_sw = real3d_k();
418  aero_tau_lw = real3d_k();
419 
420  // 3d size (ncol, nlay, n[sw,lw]bnds)
423 
424  // 3d size (ncol, nlay, n[sw,lw]gpts)
427 }

◆ finalize_impl()

void Radiation::finalize_impl ( amrex::Vector< amrex::MultiFab * > &  lsm_output_ptrs)
1243 {
1244  // Finish rrtmgp
1245  m_gas_concs.reset();
1247 
1248  // Fill the AMReX MFs from Kokkos Views
1249  kokkos_buffers_to_mf(lsm_output_ptrs);
1250 
1251  // Write fluxes if requested
1253 
1254  // Fill output data for datalog before deallocating
1255  if (datalog_int > 0) {
1258 
1260  }
1261 
1262  // Deallocate the buffer arrays
1263  dealloc_buffers();
1264 }
int datalog_int
Definition: ERF_RadiationInterface.H:87
void dealloc_buffers()
Definition: ERF_Radiation.cpp:325
void populateDatalogMF()
Definition: ERF_Radiation.cpp:751
void kokkos_buffers_to_mf(amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
Definition: ERF_Radiation.cpp:610
GasConcsK< amrex::Real, layout_t, KokkosDefaultDevice > m_gas_concs
Definition: ERF_Radiation.H:307
void write_rrtmgp_fluxes()
Definition: ERF_Radiation.cpp:710
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.

180  {
181  return m_lsm_input_names;
182  }
amrex::Vector< std::string > m_lsm_input_names
Definition: ERF_Radiation.H:240

◆ get_lsm_output_varnames()

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

Reimplemented from IRadiation.

188  {
189  return m_lsm_output_names;
190  }
amrex::Vector< std::string > m_lsm_output_names
Definition: ERF_Radiation.H:245

◆ 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:317

◆ initialize_impl()

void Radiation::initialize_impl ( )
976 {
977  // Call API to initialize
982 }
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)
611 {
612  // Heating rate, fluxes, zenith, lsm ptrs
613  Vector<real2d_k> rrtmgp_out_vars = {sw_flux_dn, lw_flux_dn};
614 
615  // Expose for device
616  auto sw_heating_d = sw_heating;
617  auto lw_heating_d = lw_heating;
618  auto p_lay_d = p_lay;
619  auto sfc_flux_dir_vis_d = sfc_flux_dir_vis;
620  auto sfc_flux_dir_nir_d = sfc_flux_dir_nir;
621  auto sfc_flux_dif_vis_d = sfc_flux_dif_vis;
622  auto sfc_flux_dif_nir_d = sfc_flux_dif_nir;
623  auto lw_flux_dn_d = lw_flux_dn;
624  auto mu0_d = mu0;
625 
626  for (MFIter mfi(*m_cons_in); mfi.isValid(); ++mfi) {
627  const auto& vbx = mfi.validbox();
628  const auto& sbx = makeSlab(vbx,2,vbx.smallEnd(2));
629  const int nx = vbx.length(0);
630  const int imin = vbx.smallEnd(0);
631  const int jmin = vbx.smallEnd(1);
632  const int offset = m_col_offsets[mfi.index()];
633  const Array4<Real>& q_arr = m_qheating_rates->array(mfi);
634  ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
635  {
636  // map [i,j,k] 0-based to [icol, ilay] 0-based
637  const int icol = (j-jmin)*nx + (i-imin) + offset;
638  const int ilay = k;
639 
640  // Temperature heating rate for SW and LW
641  q_arr(i,j,k,0) = sw_heating_d(icol,ilay);
642  q_arr(i,j,k,1) = lw_heating_d(icol,ilay);
643 
644  // Convert the rates for theta_d
645  Real exner = getExnergivenP(Real(p_lay_d(icol,ilay)), R_d/Cp_d);
646  q_arr(i,j,k,0) *= exner;
647  q_arr(i,j,k,1) *= exner;
648  });
649  if (m_lsm_fluxes) {
650  const Array4<Real>& lsm_arr = m_lsm_fluxes->array(mfi);
651  ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
652  {
653  // map [i,j,k] 0-based to [icol, ilay] 0-based
654  const int icol = (j-jmin)*nx + (i-imin) + offset;
655 
656  // SW fluxes for LSM
657  lsm_arr(i,j,k,0) = sfc_flux_dir_vis_d(icol);
658  lsm_arr(i,j,k,1) = sfc_flux_dir_nir_d(icol);
659  lsm_arr(i,j,k,2) = sfc_flux_dif_vis_d(icol);
660  lsm_arr(i,j,k,3) = sfc_flux_dif_nir_d(icol);
661 
662  // Net SW flux for LSM
663  lsm_arr(i,j,k,4) = sfc_flux_dir_vis_d(icol) + sfc_flux_dir_nir_d(icol)
664  + sfc_flux_dif_vis_d(icol) + sfc_flux_dif_nir_d(icol);
665 
666  // LW flux for LSM (at bottom surface)
667  lsm_arr(i,j,k,5) = lw_flux_dn_d(icol,0);
668  });
669  }
670  if (m_lsm_zenith) {
671  const Array4<Real>& lsm_zenith_arr = m_lsm_zenith->array(mfi);
672  ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
673  {
674  // map [i,j,k] 0-based to [icol, ilay] 0-based
675  const int icol = (j-jmin)*nx + (i-imin) + offset;
676 
677  // export cosine zenith angle for LSM
678  lsm_zenith_arr(i,j,k) = mu0_d(icol);
679  });
680  }
681  for (int ivar(0); ivar<lsm_output_ptrs.size(); ivar++) {
682  if (lsm_output_ptrs[ivar]) {
683  const Array4<Real>& lsm_out_arr = lsm_output_ptrs[ivar]->array(mfi);
684  if (ivar==0) {
685  ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
686  {
687  // map [i,j,k] 0-based to [icol, ilay] 0-based
688  const int icol = (j-jmin)*nx + (i-imin) + offset;
689 
690  // export the desired variable at surface
691  lsm_out_arr(i,j,k) = mu0_d(icol);
692  });
693  } else {
694  auto rrtmgp_for_fill = rrtmgp_out_vars[ivar-1];
695  ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
696  {
697  // map [i,j,k] 0-based to [icol, ilay] 0-based
698  const int icol = (j-jmin)*nx + (i-imin) + offset;
699 
700  // export the desired variable at surface
701  lsm_out_arr(i,j,k) = rrtmgp_for_fill(icol,0);
702  });
703  } // ivar
704  } // valid ptr
705  } // ivar
706  }// mfi
707 }
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:251
amrex::MultiFab * m_lsm_fluxes
Definition: ERF_Radiation.H:268
amrex::MultiFab * m_qheating_rates
Definition: ERF_Radiation.H:254
amrex::MultiFab * m_lsm_zenith
Definition: ERF_Radiation.H:269
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)
432 {
433  // Expose for device
434  auto r_lay_d = r_lay;
435  auto p_lay_d = p_lay;
436  auto t_lay_d = t_lay;
437  auto z_del_d = z_del;
438  auto qv_lay_d = qv_lay;
439  auto qc_lay_d = qc_lay;
440  auto qi_lay_d = qi_lay;
441  auto cldfrac_tot_d = cldfrac_tot;
442  auto lwp_d = lwp;
443  auto iwp_d = iwp;
444  auto eff_radius_qc_d = eff_radius_qc;
445  auto eff_radius_qi_d = eff_radius_qi;
446  auto p_lev_d = p_lev;
447  auto t_lev_d = t_lev;
448  auto lat_d = lat;
449  auto lon_d = lon;
450  auto t_sfc_d = t_sfc;
451  auto p_del_d = p_del;
452 
453  bool moist = m_moist;
454  bool ice = m_ice;
455  const bool has_lsm = m_lsm;
456  const bool has_lat = m_lat;
457  const bool has_lon = m_lon;
458  int ncol = m_ncol;
459  int nlay = m_nlay;
460  Real dz = m_geom.CellSize(2);
461  Real cons_lat = m_lat_cons;
462  Real cons_lon = m_lon_cons;
463  Real rad_t_sfc = m_rad_t_sfc;
464  for (MFIter mfi(*m_cons_in); mfi.isValid(); ++mfi) {
465  const auto& vbx = mfi.validbox();
466  const int nx = vbx.length(0);
467  const int imin = vbx.smallEnd(0);
468  const int jmin = vbx.smallEnd(1);
469  const int offset = m_col_offsets[mfi.index()];
470  const Array4<const Real>& cons_arr = m_cons_in->const_array(mfi);
471  const Array4<const Real>& z_arr = (m_z_phys) ? m_z_phys->const_array(mfi) :
472  Array4<const Real>{};
473  const Array4<const Real>& lat_arr = (m_lat) ? m_lat->const_array(mfi) :
474  Array4<const Real>{};
475  const Array4<const Real>& lon_arr = (m_lon) ? m_lon->const_array(mfi) :
476  Array4<const Real>{};
477  ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
478  {
479  // map [i,j,k] 0-based to [icol, ilay] 0-based
480  const int icol = (j-jmin)*nx + (i-imin) + offset;
481  const int ilay = k;
482 
483  // EOS input (at CC)
484  Real r = cons_arr(i,j,k,Rho_comp);
485  Real rt = cons_arr(i,j,k,RhoTheta_comp);
486  Real qv = (moist) ? std::max(cons_arr(i,j,k,RhoQ1_comp)/r,0.0) : 0.0;
487  Real qc = (moist) ? std::max(cons_arr(i,j,k,RhoQ2_comp)/r,0.0) : 0.0;
488  Real qi = (ice) ? std::max(cons_arr(i,j,k,RhoQ3_comp)/r,0.0) : 0.0;
489 
490  // EOS avg to z-face
491  Real r_lo = cons_arr(i,j,k-1,Rho_comp);
492  Real rt_lo = cons_arr(i,j,k-1,RhoTheta_comp);
493  Real qv_lo = (moist) ? cons_arr(i,j,k-1,RhoQ1_comp)/r_lo : 0.0;
494  Real r_avg = 0.5 * (r + r_lo);
495  Real rt_avg = 0.5 * (rt + rt_lo);
496  Real qv_avg = 0.5 * (qv + qv_lo);
497 
498  // Views at CC
499  r_lay_d(icol,ilay) = r;
500  p_lay_d(icol,ilay) = getPgivenRTh(rt, qv);
501  t_lay_d(icol,ilay) = getTgivenRandRTh(r, rt, qv);
502  z_del_d(icol,ilay) = (z_arr) ? 0.25 * ( (z_arr(i ,j ,k+1) - z_arr(i ,j ,k))
503  + (z_arr(i+1,j ,k+1) - z_arr(i+1,j ,k))
504  + (z_arr(i ,j+1,k+1) - z_arr(i ,j+1,k))
505  + (z_arr(i+1,j+1,k+1) - z_arr(i+1,j+1,k)) ) : dz;
506  qv_lay_d(icol,ilay) = qv;
507  qc_lay_d(icol,ilay) = qc;
508  qi_lay_d(icol,ilay) = qi;
509  cldfrac_tot_d(icol,ilay) = ((qc+qi)>0.0) ? 1. : 0.;
510 
511  // NOTE: These are populated in 'mixing_ratio_to_cloud_mass'
512  lwp_d(icol,ilay) = 0.0;
513  iwp_d(icol,ilay) = 0.0;
514 
515  // NOTE: These would be populated from P3 (we use the constants in p3_main_impl.hpp)
516  eff_radius_qc_d(icol,ilay) = (qc>0.0) ? 10.0e-6 : 0.0;
517  eff_radius_qi_d(icol,ilay) = (qi>0.0) ? 25.0e-6 : 0.0;
518 
519  // Buffers on z-faces (nlay+1)
520  p_lev_d(icol,ilay) = getPgivenRTh(rt_avg, qv_avg);
521  t_lev_d(icol,ilay) = getTgivenRandRTh(r_avg, rt_avg, qv_avg);
522  if (ilay==(nlay-1)) {
523  Real r_hi = cons_arr(i,j,k+1,Rho_comp);
524  Real rt_hi = cons_arr(i,j,k+1,RhoTheta_comp);
525  Real qv_hi = (moist) ? std::max(cons_arr(i,j,k+1,RhoQ1_comp)/r_hi,0.0) : 0.0;
526  r_avg = 0.5 * (r + r_hi);
527  rt_avg = 0.5 * (rt + rt_hi);
528  qv_avg = 0.5 * (qv + qv_hi);
529  p_lev_d(icol,ilay+1) = getPgivenRTh(rt_avg, qv_avg);
530  t_lev_d(icol,ilay+1) = getTgivenRandRTh(r_avg, rt_avg, qv_avg);
531  }
532 
533  // 1D data structures
534  if (k==0) {
535  lat_d(icol) = (has_lat) ? lat_arr(i,j,0) : cons_lat;
536  lon_d(icol) = (has_lon) ? lon_arr(i,j,0) : cons_lon;
537  }
538 
539  });
540  } // mfi
541 
542  // EAMXX delP is positive
543  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
544  KOKKOS_LAMBDA (int icol, int ilay)
545  {
546  p_del_d(icol,ilay) = std::abs(p_lev_d(icol,ilay+1) - p_lev_d(icol,ilay));
547  });
548 
549  // Populate vars LSM would provide
550  if (!has_lsm) {
551  // Parsed surface temp
552  Kokkos::deep_copy(t_sfc, rad_t_sfc);
553 
554  // EAMXX dummy atmos constants
555  Kokkos::deep_copy(sfc_alb_dir_vis, 0.06);
556  Kokkos::deep_copy(sfc_alb_dir_nir, 0.06);
557  Kokkos::deep_copy(sfc_alb_dif_vis, 0.06);
558  Kokkos::deep_copy(sfc_alb_dif_nir, 0.06);
559 
560  // AML NOTE: These are not used in current EAMXX, I've left
561  // the code to plug into these if we need it.
562  //
563  // Current EAMXX constants
564  Kokkos::deep_copy(emis_sfc, 0.98);
565  Kokkos::deep_copy(lw_src , 0.0 );
566  } else {
567  Vector<real1d_k> rrtmgp_in_vars = {t_sfc, sfc_emis,
570  Vector<Real> rrtmgp_default_vals = {rad_t_sfc, 0.98,
571  0.06, 0.06,
572  0.06, 0.06};
573  for (int ivar(0); ivar<lsm_input_ptrs.size(); ivar++) {
574  auto rrtmgp_to_fill = rrtmgp_in_vars[ivar];
575  if (!lsm_input_ptrs[ivar]) {
576  Kokkos::deep_copy(rrtmgp_to_fill, rrtmgp_default_vals[ivar]);
577  } else {
578  for (MFIter mfi(*m_cons_in); mfi.isValid(); ++mfi) {
579  const auto& vbx = mfi.validbox();
580  const auto& sbx = makeSlab(vbx,2,vbx.smallEnd(2));
581  const int nx = vbx.length(0);
582  const int imin = vbx.smallEnd(0);
583  const int jmin = vbx.smallEnd(1);
584  const int offset = m_col_offsets[mfi.index()];
585  const Array4<const Real>& lsm_in_arr = lsm_input_ptrs[ivar]->const_array(mfi);
586  ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
587  {
588  // map [i,j,k] 0-based to [icol, ilay] 0-based
589  const int icol = (j-jmin)*nx + (i-imin) + offset;
590 
591  // 2D mf and 1D view
592  rrtmgp_to_fill(icol) = lsm_in_arr(i,j,k);
593  });
594  } //mfi
595  } // valid lsm ptr
596  } // ivar
597  Kokkos::deep_copy(lw_src, 0.0 );
598  } // have lsm
599 
600  // Enforce consistency between t_sfc and t_lev at bottom surface
601  Kokkos::parallel_for(Kokkos::RangePolicy(0, ncol),
602  KOKKOS_LAMBDA (int icol)
603  {
604  t_lev_d(icol,0) = t_sfc_d(icol);
605  });
606 }
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:261
amrex::MultiFab * m_z_phys
Definition: ERF_Radiation.H:257
amrex::Geometry m_geom
Definition: ERF_Radiation.H:219
amrex::MultiFab * m_lat
Definition: ERF_Radiation.H:260
@ qv
Definition: ERF_Kessler.H:28
@ qc
Definition: ERF_SatAdj.H:36
Here is the call graph for this function:

◆ populateDatalogMF()

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

◆ rad_run_impl()

void Radiation::rad_run_impl ( amrex::Vector< amrex::MultiFab * > &  lsm_output_ptrs)
inline
165  {
166  if (m_update_rad) {
167  amrex::Print() << "Radiation advancing level " << m_lev << " at (YY-MM-DD SS) " << m_orbital_year << '-'
168  << m_orbital_mon << '-' << m_orbital_day << ' ' << m_orbital_sec << " ...";
169  this->initialize_impl();
170  this->run_impl();
171  this->finalize_impl(lsm_output_ptrs);
172  amrex::Print() << "DONE\n";
173  }
174  }
void initialize_impl()
Definition: ERF_Radiation.cpp:975
void run_impl()
Definition: ERF_Radiation.cpp:986
int m_orbital_mon
Definition: ERF_Radiation.H:330
void finalize_impl(amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
Definition: ERF_Radiation.cpp:1242
int m_orbital_sec
Definition: ERF_Radiation.H:332
int m_lev
Definition: ERF_Radiation.H:207
bool m_update_rad
Definition: ERF_Radiation.H:225
int m_orbital_day
Definition: ERF_Radiation.H:331

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:164
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:134
Here is the call graph for this function:

◆ run_impl()

void Radiation::run_impl ( )
987 {
988  // Local copies
989  const auto ncol = m_ncol;
990  const auto nlay = m_nlay;
991  const auto nswbands = m_nswbands;
992 
993  // Compute orbital parameters; these are used both for computing
994  // the solar zenith angle and also for computing total solar
995  // irradiance scaling (tsi_scaling).
996  double obliqr, lambm0, mvelpp;
997  int orbital_year = m_orbital_year;
998  double eccen = m_orbital_eccen;
999  double obliq = m_orbital_obliq;
1000  double mvelp = m_orbital_mvelp;
1001  if (eccen >= 0 && obliq >= 0 && mvelp >= 0) {
1002  // fixed orbital parameters forced with orbital_year == ORB_UNDEF_INT
1003  orbital_year = ORB_UNDEF_INT;
1004  }
1005  orbital_params(orbital_year, eccen, obliq,
1006  mvelp, obliqr, lambm0, mvelpp);
1007 
1008  // Use the orbital parameters to calculate the solar declination and eccentricity factor
1009  double delta, eccf;
1010  // Want day + fraction; calday 1 == Jan 1 0Z
1011  static constexpr double dpy[] = {0.0 , 31.0, 59.0, 90.0, 120.0, 151.0,
1012  181.0, 212.0, 243.0, 273.0, 304.0, 334.0};
1013  bool leap = (m_orbital_year % 4 == 0 && (!(m_orbital_year % 100 == 0) || (m_orbital_year % 400 == 0))) ? true : false;
1014  double calday = dpy[m_orbital_mon-1] + (m_orbital_day-1.0) + m_orbital_sec/86400.0;
1015  // add extra day if leap year
1016  if (leap) { calday += 1.0; }
1017  orbital_decl(calday, eccen, mvelpp, lambm0, obliqr, delta, eccf);
1018 
1019  // Overwrite eccf if using a fixed solar constant.
1020  auto fixed_total_solar_irradiance = m_fixed_total_solar_irradiance;
1021  if (fixed_total_solar_irradiance >= 0){
1022  eccf = fixed_total_solar_irradiance/1360.9;
1023  }
1024 
1025  // Precompute volume mixing ratio (VMR) for all gases
1026  //
1027  // H2O is obtained from qv.
1028  // O3 may be a constant or a 1D vector
1029  // All other comps are set to constants for now
1030  for (int igas(0); igas < m_ngas; ++igas) {
1031  auto tmp2d = Kokkos::View<RealT**,layout_t,KokkosDefaultMem>("tmp2d", ncol, nlay);
1032  auto name = m_gas_names[igas];
1033  auto gas_mol_weight = m_mol_weight_gas[igas];
1034  if (name == "H2O") {
1035  auto qv_lay_d = qv_lay;
1036  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
1037  KOKKOS_LAMBDA (int icol, int ilay)
1038  {
1039  tmp2d(icol,ilay) = qv_lay_d(icol,ilay) * mwdair/gas_mol_weight;
1040  });
1041  } else if (name == "CO2") {
1042  Kokkos::deep_copy(tmp2d, m_co2vmr);
1043  } else if (name == "O3") {
1044  if (m_o3_size==1) {
1045  Kokkos::deep_copy(tmp2d, m_o3vmr[0] );
1046  } else {
1047  auto o3_lay_d = o3_lay;
1048  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
1049  KOKKOS_LAMBDA (int icol, int ilay)
1050  {
1051  tmp2d(icol,ilay) = o3_lay_d(ilay);
1052  });
1053  }
1054  } else if (name == "N2O") {
1055  Kokkos::deep_copy(tmp2d, m_n2ovmr);
1056  } else if (name == "CO") {
1057  Kokkos::deep_copy(tmp2d, m_covmr );
1058  } else if (name == "CH4") {
1059  Kokkos::deep_copy(tmp2d, m_ch4vmr);
1060  } else if (name == "O2") {
1061  Kokkos::deep_copy(tmp2d, m_o2vmr );
1062  } else if (name == "N2") {
1063  Kokkos::deep_copy(tmp2d, m_n2vmr );
1064  } else {
1065  Abort("Radiation: Unknown gas component.");
1066  }
1067 
1068  // Populate GasConcs object
1069  m_gas_concs.set_vmr(name, tmp2d);
1070  Kokkos::fence();
1071  }
1072 
1073  // Populate mu0 1D array
1074  // This must be done on HOST and copied to device.
1075  auto h_mu0 = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), mu0);
1076  if (m_fixed_solar_zenith_angle > 0) {
1077  Kokkos::deep_copy(h_mu0, m_fixed_solar_zenith_angle);
1078  } else {
1079  auto h_lat = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), lat);
1080  auto h_lon = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), lon);
1081  double dt = double(m_dt);
1082  auto rad_freq_in_steps = m_rad_freq_in_steps;
1083  Kokkos::parallel_for(Kokkos::RangePolicy<Kokkos::Serial>(0, ncol),
1084  [&] (int icol)
1085  {
1086  // Convert lat/lon to radians
1087  double lat_col = h_lat(icol)*PI/180.0;
1088  double lon_col = h_lon(icol)*PI/180.0;
1089  double lcalday = calday;
1090  double ldelta = delta;
1091  h_mu0(icol) = Real(orbital_cos_zenith(lcalday, lat_col, lon_col, ldelta, rad_freq_in_steps * dt));
1092  });
1093  }
1094  Kokkos::deep_copy(mu0, h_mu0);
1095 
1096  // Compute layer cloud mass per unit area (populates lwp/iwp)
1099 
1100  // Convert to g/m2 (needed by RRTMGP)
1101  auto lwp_d = lwp;
1102  auto iwp_d = iwp;
1103  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
1104  KOKKOS_LAMBDA (int icol, int ilay)
1105  {
1106  lwp_d(icol,ilay) *= 1.e3;
1107  iwp_d(icol,ilay) *= 1.e3;
1108  });
1109 
1110  // Expand surface_albedos along nswbands.
1111  // This is needed since rrtmgp require band-by-band.
1116  // Run RRTMGP driver
1118  p_lay, t_lay,
1119  p_lev, t_lev,
1120  m_gas_concs,
1122  t_sfc, emis_sfc, lw_src,
1138 
1139 #if 0
1140  // UNIT TEST
1141  //================================================================================
1142  Kokkos::deep_copy(mu0, 0.86);
1143  Kokkos::deep_copy(sfc_alb_dir_vis, 0.06);
1144  Kokkos::deep_copy(sfc_alb_dir_nir, 0.06);
1145  Kokkos::deep_copy(sfc_alb_dif_vis, 0.06);
1146  Kokkos::deep_copy(sfc_alb_dif_nir, 0.06);
1147 
1148  Kokkos::deep_copy(aero_tau_sw, 0.0);
1149  Kokkos::deep_copy(aero_ssa_sw, 0.0);
1150  Kokkos::deep_copy(aero_g_sw , 0.0);
1151  Kokkos::deep_copy(aero_tau_lw, 0.0);
1152 
1153  // Generate some fake liquid and ice water data. We pick values to be midway between
1154  // the min and max of the valid lookup table values for effective radii
1155  real rel_val = 0.5 * (rrtmgp::cloud_optics_sw_k->get_min_radius_liq()
1156  + rrtmgp::cloud_optics_sw_k->get_max_radius_liq());
1157  real rei_val = 0.5 * (rrtmgp::cloud_optics_sw_k->get_min_radius_ice()
1158  + rrtmgp::cloud_optics_sw_k->get_max_radius_ice());
1159 
1160  // Restrict clouds to troposphere (> 100 hPa = 100*100 Pa) and not very close to the ground (< 900 hPa), and
1161  // put them in 2/3 of the columns since that's roughly the total cloudiness of earth.
1162  // Set sane values for liquid and ice water path.
1163  // NOTE: these "sane" values are in g/m2!
1164  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
1165  KOKKOS_LAMBDA (int icol, int ilay)
1166  {
1167  cldfrac_tot(icol,ilay) = (p_lay(icol,ilay) > 100. * 100.) &&
1168  (p_lay(icol,ilay) < 900. * 100.) &&
1169  (icol%3 != 0);
1170  // Ice and liquid will overlap in a few layers
1171  lwp(icol,ilay) = (cldfrac_tot(icol,ilay) && t_lay(icol,ilay) > 263.) ? 10. : 0.;
1172  iwp(icol,ilay) = (cldfrac_tot(icol,ilay) && t_lay(icol,ilay) < 273.) ? 10. : 0.;
1173  eff_radius_qc(icol,ilay) = (lwp(icol,ilay) > 0.) ? rel_val : 0.;
1174  eff_radius_qi(icol,ilay) = (iwp(icol,ilay) > 0.) ? rei_val : 0.;
1175  });
1176 
1181 
1183  p_lay, t_lay,
1184  p_lev, t_lev,
1185  m_gas_concs,
1187  t_sfc, emis_sfc, lw_src,
1202  1.0, false, false);
1203  //================================================================================
1204 #endif
1205 
1206  // Update heating tendency
1209 
1210  /*
1211  // AML DEBUG
1212  Kokkos::parallel_for(nlay+1, KOKKOS_LAMBDA (int ilay)
1213  {
1214  printf("Fluxes: %i %e %e %e %e %e\n",ilay,
1215  sw_flux_up(5,ilay), sw_flux_dn(5,ilay), sw_flux_dn_dir(5,ilay),
1216  lw_flux_up(5,ilay), lw_flux_dn(5,ilay));
1217  });
1218  Kokkos::parallel_for(nlay, KOKKOS_LAMBDA (int ilay)
1219  {
1220  printf("Heating Rate: %i %e %e\n",ilay,sw_heating(5,ilay),lw_heating(5,ilay));
1221  });
1222  */
1223 
1224  // Compute surface fluxes
1225  const int kbot = 0;
1226  auto sw_bnd_flux_dif_d = sw_bnd_flux_dif;
1227  auto sw_bnd_flux_dn_d = sw_bnd_flux_dn;
1228  auto sw_bnd_flux_dir_d = sw_bnd_flux_dir;
1229  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<3>>({0, 0, 0}, {ncol, nlay+1, nswbands}),
1230  KOKKOS_LAMBDA (int icol, int ilay, int ibnd)
1231  {
1232  sw_bnd_flux_dif_d(icol,ilay,ibnd) = sw_bnd_flux_dn_d(icol,ilay,ibnd) - sw_bnd_flux_dir_d(icol,ilay,ibnd);
1233  });
1234  rrtmgp::compute_broadband_surface_fluxes(ncol, kbot, nswbands,
1238 }
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:216
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 
)
149 {
150  // Set data members that may change
151  m_lev = level;
152  m_step = step;
153  m_time = time;
154  m_dt = dt;
155  m_geom = geom;
156  m_cons_in = cons_in;
157  m_lsm_fluxes = lsm_fluxes;
158  m_lsm_zenith = lsm_zenith;
159  m_qheating_rates = qheating_rates;
160  m_z_phys = z_phys;
161  m_lat = lat;
162  m_lon = lon;
163 
164  // Update the day and month
165  time_t timestamp = time_t(time);
166  struct tm *timeinfo = gmtime(&timestamp);
167  if (m_fixed_orbital_year) {
168  m_orbital_mon = timeinfo->tm_mon + 1;
169  m_orbital_day = timeinfo->tm_mday;
170  m_orbital_sec = timeinfo->tm_hour*3600 + timeinfo->tm_min*60 + timeinfo->tm_sec;
171  } else {
172  m_orbital_year = timeinfo->tm_year + 1900;
173  m_orbital_mon = timeinfo->tm_mon + 1;
174  m_orbital_day = timeinfo->tm_mday;
175  m_orbital_sec = timeinfo->tm_hour*3600 + timeinfo->tm_min*60 + timeinfo->tm_sec;
176  }
177 
178  // Only allocate and proceed if we are going to update radiation
179  m_update_rad = false;
180  if (m_rad_freq_in_steps > 0) { m_update_rad = ( (m_step == 0) || (m_step % m_rad_freq_in_steps == 0) ); }
181 
182  if (m_update_rad) {
183  // Call to Init() has set the dimensions: ncol & nlay
184 
185  // Allocate the buffer arrays
186  alloc_buffers();
187 
188  // Fill the KOKKOS Views from AMReX MFs
189  mf_to_kokkos_buffers(lsm_input_ptrs);
190 
191  // Initialize datalog MF on first step
192  if (m_first_step) {
193  m_first_step = false;
194  datalog_mf.define(cons_in->boxArray(), cons_in->DistributionMap(), 25, 0);
195  datalog_mf.setVal(0.0);
196  }
197  }
198 }
int m_step
Definition: ERF_Radiation.H:210
void mf_to_kokkos_buffers(amrex::Vector< amrex::MultiFab * > &lsm_input_ptrs)
Definition: ERF_Radiation.cpp:431
bool m_first_step
Definition: ERF_Radiation.H:230
void alloc_buffers()
Definition: ERF_Radiation.cpp:201
amrex::Real m_time
Definition: ERF_Radiation.H:213

Referenced by Run().

Here is the caller graph for this function:

◆ write_rrtmgp_fluxes()

void Radiation::write_rrtmgp_fluxes ( )
711 {
712  // Expose for device
713  auto sw_flux_up_d = sw_flux_up;
714  auto sw_flux_dn_d = sw_flux_dn;
715  auto sw_flux_dn_dir_d = sw_flux_dn_dir;
716  auto lw_flux_up_d = lw_flux_up;
717  auto lw_flux_dn_d = lw_flux_dn;
718 
719  int n_fluxes = 5;
720  MultiFab mf_flux(m_cons_in->boxArray(), m_cons_in->DistributionMap(), n_fluxes, 0);
721 
722  for (MFIter mfi(mf_flux); mfi.isValid(); ++mfi) {
723  const auto& vbx = mfi.validbox();
724  const int nx = vbx.length(0);
725  const int imin = vbx.smallEnd(0);
726  const int jmin = vbx.smallEnd(1);
727  const int offset = m_col_offsets[mfi.index()];
728  const Array4<Real>& dst_arr = mf_flux.array(mfi);
729  ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
730  {
731  // map [i,j,k] 0-based to [icol, ilay] 0-based
732  const int icol = (j-jmin)*nx + (i-imin) + offset;
733  const int ilay = k;
734 
735  // SW and LW fluxes
736  dst_arr(i,j,k,0) = sw_flux_up_d(icol,ilay);
737  dst_arr(i,j,k,1) = sw_flux_dn_d(icol,ilay);
738  dst_arr(i,j,k,2) = sw_flux_dn_dir_d(icol,ilay);
739  dst_arr(i,j,k,3) = lw_flux_up_d(icol,ilay);
740  dst_arr(i,j,k,4) = lw_flux_dn_d(icol,ilay);
741  });
742  }
743 
744 
745  std::string plotfilename = amrex::Concatenate("plt_rad", m_step, 5);
746  Vector<std::string> flux_names = {"sw_flux_up", "sw_flux_dn", "sw_flux_dir",
747  "lw_flux_up", "lw_flux_dn"};
748  WriteSingleLevelPlotfile(plotfilename, mf_flux, flux_names, m_geom, m_time, m_step);
749 }
Here is the call graph for this function:

◆ WriteDataLog()

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

Implements IRadiation.

841 {
842  constexpr int datwidth = 14;
843  constexpr int datprecision = 9;
844  constexpr int timeprecision = 13;
845 
846  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;
847  // Clear sky
848  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;
849  // Clean sky
850  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;
851  // Clean clear sky
852  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;
853 
854 
855  auto domain = m_geom.Domain();
856  h_avg_radqrsw = sumToLine(datalog_mf, 0, 1, domain, 2);
857  h_avg_radqrlw = sumToLine(datalog_mf, 1, 1, domain, 2);
858  h_avg_sw_up = sumToLine(datalog_mf, 2, 1, domain, 2);
859  h_avg_sw_dn = sumToLine(datalog_mf, 3, 1, domain, 2);
860  h_avg_sw_dn_dir = sumToLine(datalog_mf, 4, 1, domain, 2);
861  h_avg_lw_up = sumToLine(datalog_mf, 5, 1, domain, 2);
862  h_avg_lw_dn = sumToLine(datalog_mf, 6, 1, domain, 2);
863  h_avg_zenith = sumToLine(datalog_mf, 7, 1, domain, 2);
864 
865  h_avg_radqrcsw = sumToLine(datalog_mf, 8, 1, domain, 2);
866  h_avg_radqrclw = sumToLine(datalog_mf, 9, 1, domain, 2);
867  h_avg_sw_clr_up = sumToLine(datalog_mf, 10, 1, domain, 2);
868  h_avg_sw_clr_dn = sumToLine(datalog_mf, 11, 1, domain, 2);
869  h_avg_sw_clr_dn_dir = sumToLine(datalog_mf, 12, 1, domain, 2);
870  h_avg_lw_clr_up = sumToLine(datalog_mf, 13, 1, domain, 2);
871  h_avg_lw_clr_dn = sumToLine(datalog_mf, 14, 1, domain, 2);
872 
873  if (m_extra_clnsky_diag) {
874  h_avg_sw_cln_up = sumToLine(datalog_mf, 15, 1, domain, 2);
875  h_avg_sw_cln_dn = sumToLine(datalog_mf, 16, 1, domain, 2);
876  h_avg_sw_cln_dn_dir = sumToLine(datalog_mf, 17, 1, domain, 2);
877  h_avg_lw_cln_up = sumToLine(datalog_mf, 18, 1, domain, 2);
878  h_avg_lw_cln_dn = sumToLine(datalog_mf, 19, 1, domain, 2);
879  }
880 
882  h_avg_sw_clnclr_up = sumToLine(datalog_mf, 20, 1, domain, 2);
883  h_avg_sw_clnclr_dn = sumToLine(datalog_mf, 21, 1, domain, 2);
884  h_avg_sw_clnclr_dn_dir = sumToLine(datalog_mf, 22, 1, domain, 2);
885  h_avg_lw_clnclr_up = sumToLine(datalog_mf, 23, 1, domain, 2);
886  h_avg_lw_clnclr_dn = sumToLine(datalog_mf, 24, 1, domain, 2);
887  }
888 
889  Real area_z = static_cast<Real>(domain.length(0)*domain.length(1));
890  int nz = domain.length(2);
891  for (int k = 0; k < nz; k++) {
892  h_avg_radqrsw[k] /= area_z;
893  h_avg_radqrlw[k] /= area_z;
894  h_avg_sw_up[k] /= area_z;
895  h_avg_sw_dn[k] /= area_z;
896  h_avg_sw_dn_dir[k] /= area_z;
897  h_avg_lw_up[k] /= area_z;
898  h_avg_lw_dn[k] /= area_z;
899  h_avg_zenith[k] /= area_z;
900 
901  h_avg_radqrcsw[k] /= area_z;
902  h_avg_radqrclw[k] /= area_z;
903  h_avg_sw_clr_up[k] /= area_z;
904  h_avg_sw_clr_dn[k] /= area_z;
905  h_avg_sw_clr_dn_dir[k] /= area_z;
906  h_avg_lw_clr_up[k] /= area_z;
907  h_avg_lw_clr_dn[k] /= area_z;
908  }
909 
910  if (m_extra_clnsky_diag) {
911  for (int k = 0; k < nz; k++) {
912  h_avg_sw_cln_up[k] /= area_z;
913  h_avg_sw_cln_dn[k] /= area_z;
914  h_avg_sw_cln_dn_dir[k] /= area_z;
915  h_avg_lw_cln_up[k] /= area_z;
916  h_avg_lw_cln_dn[k] /= area_z;
917  }
918  }
919 
921  for (int k = 0; k < nz; k++) {
922  h_avg_sw_clnclr_up[k] /= area_z;
923  h_avg_sw_clnclr_dn[k] /= area_z;
924  h_avg_sw_clnclr_dn_dir[k] /= area_z;
925  h_avg_lw_clnclr_up[k] /= area_z;
926  h_avg_lw_clnclr_dn[k] /= area_z;
927  }
928  }
929 
930  if (ParallelDescriptor::IOProcessor()) {
931  std::ostream& log = *datalog;
932  if (log.good()) {
933 
934  for (int k = 0; k < nz; k++)
935  {
936  Real z = k * m_geom.CellSize(2);
937  log << std::setw(datwidth) << std::setprecision(timeprecision) << time << " "
938  << std::setw(datwidth) << std::setprecision(datprecision) << z << " "
939  << h_avg_radqrsw[k] << " " << h_avg_radqrlw[k] << " " << h_avg_sw_up[k] << " "
940  << h_avg_sw_dn[k] << " " << h_avg_sw_dn_dir[k] << " " << h_avg_lw_up[k] << " "
941  << h_avg_lw_dn[k] << " " << h_avg_zenith[k] << " "
942  << h_avg_radqrcsw[k] << " " << h_avg_radqrclw[k] << " " << h_avg_sw_clr_up[k] << " "
943  << h_avg_sw_clr_dn[k] << " " << h_avg_sw_clr_dn_dir[k] << " " << h_avg_lw_clr_up[k] << " "
944  << h_avg_lw_clr_dn[k] << " ";
945  if (m_extra_clnsky_diag) {
946  log << h_avg_sw_cln_up[k] << " " << h_avg_sw_cln_dn[k] << " " << h_avg_sw_cln_dn_dir[k] << " "
947  << h_avg_lw_cln_up[k] << " " << h_avg_lw_cln_dn[k] << " ";
948  } else {
949  log << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " ";
950  }
951 
953  log << h_avg_sw_clnclr_up[k] << " " << h_avg_sw_clnclr_dn[k] << " " << h_avg_sw_clnclr_dn_dir[k] << " "
954  << h_avg_lw_clnclr_up[k] << " " << h_avg_lw_clnclr_dn[k] << std::endl;
955  } else {
956  log << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << std::endl;
957  }
958  }
959  // Write top face values
960  Real z = nz * m_geom.CellSize(2);
961  log << std::setw(datwidth) << std::setprecision(timeprecision) << time << " "
962  << std::setw(datwidth) << std::setprecision(datprecision) << z << " "
963  << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " "
964  << 0.0 << " " << 0.0 << " "
965  << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " "
966  << 0.0 << " "
967  << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " "
968  << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0 << " " << 0.0
969  << std::endl;
970  }
971  }
972 }
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

◆ 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

◆ 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

◆ 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 = {"cos_zenith_angle", "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
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 = 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 = -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 = -9999.
private

◆ m_orbital_obliq

amrex::Real Radiation::m_orbital_obliq = -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_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().

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