ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
ERF_Radiation.H
Go to the documentation of this file.
1 #ifndef ERF_RADIATION_H
2 #define ERF_RADIATION_H
3 
4 /*
5  * RTE-RRTMGP radiation model interface to ERF
6  * The original code is developed by RobertPincus, and the code is open source available at:
7  * https://github.com/earth-system-radiation/rte-rrtmgp
8  * Please reference to the following paper,
9  * https://agupubs.onlinelibrary.wiley.com/doi/10.1029/2019MS001621
10  * NOTE: we use the C++ version of RTE-RRTMGP, which is the implementation of the original Fortran
11  * code using C++ KOKKOS for CUDA, HiP and SYCL application by E3SM ECP team, the C++ version
12  * of the rte-rrtmgp code is located at:
13  * https://github.com/E3SM-Project/rte-rrtmgp
14  * The RTE-RRTMGP uses BSD-3-Clause Open Source License, if you want to make changes,
15  * and modifications to the code, please refer to BSD-3-Clause Open Source License.
16  */
17 
18 #include <ctime>
19 #include <string>
20 #include <vector>
21 #include <memory>
22 
23 #include <mo_gas_concentrations.h>
24 
25 #include <Kokkos_Core.hpp>
26 
27 #include <AMReX_ParmParse.H>
28 #include <AMReX_FArrayBox.H>
29 #include <AMReX_Geometry.H>
30 #include <AMReX_TableData.H>
31 #include <AMReX_MultiFabUtil.H>
32 #include <AMReX_PlotFileUtil.H>
33 
34 #include <ERF_RadiationInterface.H>
35 #include <ERF_RRTMGP_Interface.H>
36 #include <ERF_RRTMGP_Utils.H>
37 #include <ERF_OrbCosZenith.H>
38 #include <ERF_Constants.H>
39 #include <ERF_IndexDefines.H>
40 #include <ERF_DataStruct.H>
41 #include <ERF_EOS.H>
42 
43 #include <ERF_LandSurface.H>
44 
45 class Radiation : public IRadiation {
46 public:
47 
48  #include <ERF_Kokkos_Radiation.H>
49 
50  // Constructor
51  Radiation (const int& lev,
52  SolverChoice& sc);
53 
54  // Destructor
55  ~Radiation () = default;
56 
57  virtual void Init (const amrex::Geometry& geom,
58  const amrex::BoxArray& ba,
59  amrex::MultiFab* cons_in) override {};
60 
61  virtual void Run (int& level,
62  int& step,
63  amrex::Real& time,
64  const amrex::Real& dt,
65  const amrex::BoxArray& ba,
66  amrex::Geometry& geom,
67  amrex::MultiFab* cons_in,
68  amrex::MultiFab* lsm_fluxes,
69  amrex::MultiFab* lsm_zenith,
70  amrex::MultiFab* qheating_rates,
71  amrex::MultiFab* z_phys,
72  amrex::MultiFab* lat,
73  amrex::MultiFab* lon) override
74  {
75  set_grids(level, step, time, dt, ba, geom,
76  cons_in, lsm_fluxes, lsm_zenith,
77  qheating_rates, z_phys, lat, lon);
78  rad_run_impl();
79  }
80 
81  // Set the grid info for columnar data in KOKKOS
82  void
83  set_grids (int& level,
84  int& step,
85  amrex::Real& time,
86  const amrex::Real& dt,
87  const amrex::BoxArray& ba,
88  amrex::Geometry& geom,
89  amrex::MultiFab* cons_in,
90  amrex::MultiFab* lsm_fluxes,
91  amrex::MultiFab* lsm_zenith,
92  amrex::MultiFab* qheating_rates,
93  amrex::MultiFab* z_phys,
94  amrex::MultiFab* lat,
95  amrex::MultiFab* lon);
96 
97  template<typename LandSurfaceModelType>
98  void set_lsm_inputs (LandSurfaceModelType *model)
99  {
100  if constexpr(std::is_same_v<LandSurfaceModelType, SLM>)
101  {
102  if (!m_update_rad)
103  {
104  // skip getting LSM inputs if radiation not updating this timestep
105  return;
106  }
107 
108  // get surface temperature, LW flux, and emissivities from SLM
109  /*
110  auto slm_to_rad_vars = model->export_to_RRTMGP();
111 
112  for (amrex::MFIter mfi(*slm_to_rad_vars[0]); mfi.isValid(); ++mfi) {
113  const auto& vbx = mfi.validbox();
114  const auto& sbx = amrex::makeSlab(vbx, 2, 0);
115 
116  const int nx = vbx.length(0);
117  const int imin = vbx.smallEnd(0);
118  const int jmin = vbx.smallEnd(1);
119  const int offset = m_col_offsets[mfi.index()];
120 
121  const auto &slm_tsurf = slm_to_rad_vars[0]->const_array(mfi);
122  const auto &slm_alb_dir_vis_arr = slm_to_rad_vars[2]->const_array(mfi);
123  const auto &slm_alb_dir_nir_arr = slm_to_rad_vars[4]->const_array(mfi);
124  const auto &slm_IR_emis_arr = slm_to_rad_vars[5]->const_array(mfi);
125  const auto &slm_net_rad_arr = slm_to_rad_vars[6]->const_array(mfi);
126  amrex::ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
127  {
128  // map [i,j,k] 0-based to [icol, ilay] 1-based
129  const int icol = (j-jmin)*nx + (i-imin) + 1 + offset;
130  const int ilay = 0;
131 
132  t_sfc(icol) = slm_tsurf(i, j, k);
133 
134  sfc_alb_dir_vis(icol) = slm_alb_dir_vis_arr(i, j, k);
135  sfc_alb_dir_nir(icol) = slm_alb_dir_nir_arr(i, j, k);
136  sfc_alb_dif_vis(icol) = slm_alb_dir_vis_arr(i, j, k);
137  sfc_alb_dif_nir(icol) = slm_alb_dir_nir_arr(i, j, k);
138 
139  sfc_emis(icol) = slm_IR_emis_arr(i, j, k);
140 
141  //lw_src(icol) = slm_net_rad_arr(i, j, k, SLM_NetRad::net_lwup1);
142  lw_src(icol) = 0.0; // TODO
143  });
144  }
145 
146  // expand sfc_emis and lw_src over nlwbands
147  // TODO: check if this is correct - should emissivity vary based on band?
148  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {m_nlwbands, m_ncol}),
149  KOKKOS_LAMBDA (int ibnd, int icol)
150  {
151  emis_sfc(icol, ibnd) = sfc_emis(icol);
152  });
153  */
154 
155  Kokkos::deep_copy(emis_sfc, 0.98);
156  Kokkos::deep_copy(lw_src, 0.0);
157  }
158  }
159 
160  // Initialize the temporary variables
161  void
162  alloc_buffers ();
163 
164  // Clear the temporary variables
165  void
166  dealloc_buffers ();
167 
168  // Fill KOKKOS Views from AMReX MultiFabs
169  void
171 
172  // Fill AMReX MultiFabs from KOKKOS Views
173  void
175 
176  // Write the rrtmgp fluxes
177  void
179 
180  // Initialize the implementation
181  void
182  initialize_impl ();
183 
184  // Run the implementation
185  void
186  run_impl ();
187 
188  // Finalize the implementation
189  void
190  finalize_impl ();
191 
192  // Wrapper for implementation steps
193  void
195  {
196  if (m_update_rad) {
197  amrex::Print() << "Advancing radiation at level: " << m_lev << " ...";
198  this->initialize_impl();
199  this->run_impl();
200  this->finalize_impl();
201  amrex::Print() << "DONE\n";
202  }
203  }
204 
205  void populateDatalogMF ();
206 
207  virtual void WriteDataLog (const amrex::Real &time) override;
208 
209 private:
210 
211  // Process interface vars from ERF/AMReX
212  //===================================================================================
213 
214  // Grid level
215  int m_lev;
216 
217  // Step number
218  int m_step;
219 
220  // Current time
221  amrex::Real m_time;
222 
223  // Timestep at given level
224  amrex::Real m_dt;
225 
226  // Geometry at given level
227  amrex::Geometry m_geom;
228 
229  // Boxarray at given level
230  amrex::BoxArray m_ba;
231 
232  // Are we updating radiation?
233  bool m_update_rad = false;
234 
235  // Are we writing out the fluxes?
236  bool m_rad_write_fluxes = false;
237 
238  bool m_first_step = true;
239 
240  // Do we have moisture and cold comps?
241  bool m_moist = false;
242  bool m_ice = false;
243 
244  // Do we have a land surface model?
245  bool m_lsm = false;
246 
247  // Pointer to the CC conserved vars
248  amrex::MultiFab* m_cons_in = nullptr;
249 
250  // Pointer to the radiation source terms
251  amrex::MultiFab* m_qheating_rates = nullptr;
252 
253  // Pointer to the terrain heights
254  amrex::MultiFab* m_z_phys = nullptr;
255 
256  // Pointer to latitude and longitude
257  amrex::MultiFab* m_lat = nullptr;
258  amrex::MultiFab* m_lon = nullptr;
259 
260  // Constant lat/lon if the above MFs are not valid
261  amrex::Real m_lat_cons = 39.809860;
262  amrex::Real m_lon_cons = -98.555183;
263 
264  // Pointer to output data for LSM
265  amrex::MultiFab* m_lsm_fluxes = nullptr;
266  amrex::MultiFab* m_lsm_zenith = nullptr;
267 
268  // Holds output from KOKKOS views used in the datalog
269  amrex::MultiFab datalog_mf;
270 
271  // Path, data file, and coefficient file for K-distribution
272  std::string rrtmgp_file_path = ".";
273  std::string rrtmgp_coeffs_sw = "rrtmgp-data-sw-g224-2018-12-04.nc";
274  std::string rrtmgp_coeffs_lw = "rrtmgp-data-lw-g224-2018-12-04.nc";
275  std::string rrtmgp_cloud_optics_sw = "rrtmgp-cloud-optics-coeffs-sw.nc";
276  std::string rrtmgp_cloud_optics_lw = "rrtmgp-cloud-optics-coeffs-lw.nc";
281 
282  // Active gases
283  int m_ngas = 8;
284  const std::vector<std::string> m_gas_names = {"H2O", "CO2", "O3", "N2O",
285  "CO" , "CH4", "O2", "N2" };
286  const std::vector<amrex::Real> m_mol_weight_gas = {18.01528, 44.00950, 47.9982, 44.0128,
287  28.01010, 16.04246, 31.9980, 28.0134}; // g/mol
288 
289  // Prescribed greenhouse gas surface concentrations in moles / moles air
290  amrex::Real m_co2vmr = 388.717e-6;
291  amrex::Vector<amrex::Real> m_o3vmr;
292  amrex::Real m_n2ovmr = 323.141e-9;
293  amrex::Real m_covmr = 1.0e-7;
294  amrex::Real m_ch4vmr = 1807.851e-9;
295  amrex::Real m_o2vmr = 0.209448;
296  amrex::Real m_n2vmr = 0.7906;
297  //amrex::Real m_f11vmr = 768.7644e-12;
298  //amrex::Real m_f12vmr = 531.2820e-12;
299 
302  std::vector<std::string> gas_names_offset;
303 
304  GasConcsK<amrex::Real, layout_t, RadDefaultDevice> m_gas_concs;
305 
306  // Process interface vars modeled after EAMXX
307  //===================================================================================
308 
309  // Keep track of number of columns and levels
310  int m_ncol;
311  int m_nlay;
312 
313  // Offsets for MultiFab <-> KOKKOS transfer
314  amrex::Vector<int> m_col_offsets;
315 
316  // Whether we use aerosol forcing in radiation
317  bool m_do_aerosol_rad = true;
318 
319  // Whether we do extra aerosol forcing calls
320  bool m_extra_clnsky_diag = false;
322 
323  // The orbital year, used for zenith angle calculations:
324  // If > 0, use constant orbital year for duration of simulation
325  // If < 0, use year from timestamp for orbital parameters
326  int m_orbital_year = -9999;
327  int m_orbital_mon = -9999;
328  int m_orbital_day = -9999;
329  int m_orbital_sec = -9999;
330 
331  // Orbital parameters, used for zenith angle calculations.
332  // If >= 0, bypass computation based on orbital year and use fixed parameters
333  // If < 0, compute based on orbital year, specified above
334  bool m_fixed_orbital_year = false;
335  amrex::Real m_orbital_eccen = -9999.; // Eccentricity
336  amrex::Real m_orbital_obliq = -9999.; // Obliquity
337  amrex::Real m_orbital_mvelp = -9999.; // Vernal Equinox Mean Longitude of Perihelion
338 
339  // Value for prescribing an invariant solar constant (i.e. total solar irradiance
340  // at TOA). Used for idealized experiments such as RCE. This is only used when a
341  // positive value is supplied.
342  amrex::Real m_fixed_total_solar_irradiance = -9999.;
343 
344  // Fixed solar zenith angle to use for shortwave calculations
345  // This is only used if a positive value is supplied
346  amrex::Real m_fixed_solar_zenith_angle = -9999.;
347 
348  // Need to hard-code some dimension sizes for now.
349  // TODO: find a better way of configuring this
350  int m_nswbands = 14;
351  int m_nlwbands = 16;
352  int m_nswgpts = 112;
353  int m_nlwgpts = 128;
354 
355  // Rad frequency in number of steps
357 
358  // Whether or not to do subcolumn sampling of cloud state for MCICA
359  bool m_do_subcol_sampling = true;
360 
361  // 1d size (1 or nlay)
363 
364  // 1d size (ncol)
380 
381  // 2d size (ncol, nlay)
400 
401  // 2d size (ncol, nlay+1)
425 
426  // 3d size (ncol, nlay+1, nswbands)
431 
432  // 3d size (ncol, nlay+1, nlwbands)
435 
436  // 2d size (ncol, nswbands)
439 
440  // 2d size (ncol, nlwbands)
442 
443  // 3d size (ncol, nlay, n[sw,lw]bands)
448 
449  // 3d size (ncol, nlay, n[sw,lw]bnds)
452 
453  // 3d size (ncol, nlay, n[sw,lw]gpts)
456 };
457 #endif
Kokkos::View< RealT *, RadDefaultDevice > real1d_k
Definition: ERF_Kokkos_Radiation.H:17
Kokkos::View< RealT ***, layout_t, RadDefaultDevice > real3d_k
Definition: ERF_Kokkos_Radiation.H:19
Kokkos::View< RealT **, layout_t, RadDefaultDevice > real2d_k
Definition: ERF_Kokkos_Radiation.H:18
Definition: ERF_RadiationInterface.H:14
Definition: ERF_Radiation.H:45
real3d_k sw_bnd_flux_dn
Definition: ERF_Radiation.H:428
amrex::MultiFab * m_lon
Definition: ERF_Radiation.H:258
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:300
real3d_k sw_bnd_flux_dir
Definition: ERF_Radiation.H:429
real2d_k sw_clnsky_flux_dn
Definition: ERF_Radiation.H:417
std::string rrtmgp_coeffs_file_sw
Definition: ERF_Radiation.H:277
std::string rrtmgp_coeffs_sw
Definition: ERF_Radiation.H:273
real2d_k d_tint
Definition: ERF_Radiation.H:402
real2d_k lw_clnclrsky_flux_dn
Definition: ERF_Radiation.H:420
int m_rad_freq_in_steps
Definition: ERF_Radiation.H:356
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:379
real2d_k eff_radius_qi
Definition: ERF_Radiation.H:392
void dealloc_buffers()
Definition: ERF_Radiation.cpp:315
real1d_k m_gas_mol_weights
Definition: ERF_Radiation.H:301
amrex::Real m_lon_cons
Definition: ERF_Radiation.H:262
void initialize_impl()
Definition: ERF_Radiation.cpp:833
int m_nswbands
Definition: ERF_Radiation.H:350
bool m_do_aerosol_rad
Definition: ERF_Radiation.H:317
std::string rrtmgp_coeffs_lw
Definition: ERF_Radiation.H:274
bool m_moist
Definition: ERF_Radiation.H:241
std::string rrtmgp_cloud_optics_file_lw
Definition: ERF_Radiation.H:280
amrex::MultiFab * m_cons_in
Definition: ERF_Radiation.H:248
amrex::MultiFab datalog_mf
Definition: ERF_Radiation.H:269
Radiation(const int &lev, SolverChoice &sc)
Definition: ERF_Radiation.cpp:19
real2d_k sw_heating
Definition: ERF_Radiation.H:396
bool m_lsm
Definition: ERF_Radiation.H:245
bool m_do_subcol_sampling
Definition: ERF_Radiation.H:359
real3d_k lw_bnd_flux_dn
Definition: ERF_Radiation.H:434
bool m_rad_write_fluxes
Definition: ERF_Radiation.H:236
void set_lsm_inputs(LandSurfaceModelType *model)
Definition: ERF_Radiation.H:98
real3d_k sw_bnd_flux_up
Definition: ERF_Radiation.H:427
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::MultiFab *qheating_rates, amrex::MultiFab *z_phys, amrex::MultiFab *lat, amrex::MultiFab *lon) override
Definition: ERF_Radiation.H:61
real2d_k qv_lay
Definition: ERF_Radiation.H:387
std::string rrtmgp_cloud_optics_file_sw
Definition: ERF_Radiation.H:279
void run_impl()
Definition: ERF_Radiation.cpp:844
real2d_k sw_clnclrsky_flux_dn_dir
Definition: ERF_Radiation.H:412
real1d_k sfc_flux_dif_vis
Definition: ERF_Radiation.H:373
int m_step
Definition: ERF_Radiation.H:218
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:375
void populateDatalogMF()
Definition: ERF_Radiation.cpp:638
real2d_k qi_lay
Definition: ERF_Radiation.H:389
real2d_k sw_clrsky_flux_up
Definition: ERF_Radiation.H:413
amrex::Real m_orbital_mvelp
Definition: ERF_Radiation.H:337
real2d_k t_lev
Definition: ERF_Radiation.H:404
amrex::Vector< amrex::Real > m_o3vmr
Definition: ERF_Radiation.H:291
void mf_to_kokkos_buffers()
Definition: ERF_Radiation.cpp:423
real1d_k o3_lay
Definition: ERF_Radiation.H:362
std::string rrtmgp_cloud_optics_sw
Definition: ERF_Radiation.H:275
void finalize_impl()
Definition: ERF_Radiation.cpp:1092
real1d_k sfc_alb_dif_vis
Definition: ERF_Radiation.H:369
int m_nlwgpts
Definition: ERF_Radiation.H:353
real1d_k mu0
Definition: ERF_Radiation.H:366
real2d_k cldfrac_tot
Definition: ERF_Radiation.H:390
real3d_k cld_tau_lw_bnd
Definition: ERF_Radiation.H:451
real1d_k sfc_flux_dir_nir
Definition: ERF_Radiation.H:372
amrex::MultiFab * m_z_phys
Definition: ERF_Radiation.H:254
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
std::string rrtmgp_cloud_optics_lw
Definition: ERF_Radiation.H:276
real3d_k sw_bnd_flux_dif
Definition: ERF_Radiation.H:430
real2d_k sfc_alb_dif
Definition: ERF_Radiation.H:438
amrex::Real m_co2vmr
Definition: ERF_Radiation.H:290
bool m_extra_clnsky_diag
Definition: ERF_Radiation.H:320
real2d_k emis_sfc
Definition: ERF_Radiation.H:441
real1d_k sfc_alb_dif_nir
Definition: ERF_Radiation.H:370
amrex::MultiFab * m_lsm_fluxes
Definition: ERF_Radiation.H:265
virtual void Init(const amrex::Geometry &geom, const amrex::BoxArray &ba, amrex::MultiFab *cons_in) override
Definition: ERF_Radiation.H:57
real2d_k lw_clnsky_flux_dn
Definition: ERF_Radiation.H:424
amrex::BoxArray m_ba
Definition: ERF_Radiation.H:230
real2d_k p_lay
Definition: ERF_Radiation.H:383
real2d_k r_lay
Definition: ERF_Radiation.H:382
int m_ncol
Definition: ERF_Radiation.H:310
real2d_k sw_flux_dn_dir
Definition: ERF_Radiation.H:407
const std::vector< amrex::Real > m_mol_weight_gas
Definition: ERF_Radiation.H:286
amrex::Real m_lat_cons
Definition: ERF_Radiation.H:261
amrex::Real m_fixed_total_solar_irradiance
Definition: ERF_Radiation.H:342
std::vector< std::string > gas_names_offset
Definition: ERF_Radiation.H:302
real2d_k sw_clrsky_flux_dn
Definition: ERF_Radiation.H:414
amrex::Real m_dt
Definition: ERF_Radiation.H:224
int m_orbital_mon
Definition: ERF_Radiation.H:327
real3d_k aero_ssa_sw
Definition: ERF_Radiation.H:445
real2d_k sw_clnsky_flux_dn_dir
Definition: ERF_Radiation.H:418
amrex::Real m_o2vmr
Definition: ERF_Radiation.H:295
amrex::MultiFab * m_qheating_rates
Definition: ERF_Radiation.H:251
real2d_k qc_lay
Definition: ERF_Radiation.H:388
virtual void WriteDataLog(const amrex::Real &time) override
Definition: ERF_Radiation.cpp:698
real2d_k lw_clrsky_flux_up
Definition: ERF_Radiation.H:421
amrex::Real m_n2ovmr
Definition: ERF_Radiation.H:292
std::string rrtmgp_file_path
Definition: ERF_Radiation.H:272
real2d_k sw_flux_dn
Definition: ERF_Radiation.H:406
amrex::MultiFab * m_lsm_zenith
Definition: ERF_Radiation.H:266
int m_nlwbands
Definition: ERF_Radiation.H:351
real2d_k tmp2d
Definition: ERF_Radiation.H:393
real2d_k lw_clrsky_flux_dn
Definition: ERF_Radiation.H:422
GasConcsK< amrex::Real, layout_t, RadDefaultDevice > m_gas_concs
Definition: ERF_Radiation.H:304
amrex::Real m_orbital_eccen
Definition: ERF_Radiation.H:335
amrex::Real m_covmr
Definition: ERF_Radiation.H:293
int m_orbital_sec
Definition: ERF_Radiation.H:329
real2d_k z_del
Definition: ERF_Radiation.H:385
real2d_k p_del
Definition: ERF_Radiation.H:386
real3d_k aero_tau_lw
Definition: ERF_Radiation.H:447
real2d_k lw_flux_dn
Definition: ERF_Radiation.H:409
void kokkos_buffers_to_mf()
Definition: ERF_Radiation.cpp:542
real2d_k sw_clrsky_flux_dn_dir
Definition: ERF_Radiation.H:415
void write_rrtmgp_fluxes()
Definition: ERF_Radiation.cpp:604
std::string rrtmgp_coeffs_file_lw
Definition: ERF_Radiation.H:278
int m_lev
Definition: ERF_Radiation.H:215
real1d_k sfc_flux_dif_nir
Definition: ERF_Radiation.H:374
int m_ngas
Definition: ERF_Radiation.H:283
real2d_k lw_heating
Definition: ERF_Radiation.H:397
real1d_k sfc_alb_dir_nir
Definition: ERF_Radiation.H:368
real2d_k lw_clnsky_flux_up
Definition: ERF_Radiation.H:423
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::MultiFab *qheating_rates, amrex::MultiFab *z_phys, amrex::MultiFab *lat, amrex::MultiFab *lon)
Definition: ERF_Radiation.cpp:115
amrex::Real m_ch4vmr
Definition: ERF_Radiation.H:294
bool m_update_rad
Definition: ERF_Radiation.H:233
real2d_k sfc_alb_dir
Definition: ERF_Radiation.H:437
bool m_ice
Definition: ERF_Radiation.H:242
amrex::Real m_fixed_solar_zenith_angle
Definition: ERF_Radiation.H:346
~Radiation()=default
real1d_k lon
Definition: ERF_Radiation.H:376
int m_orbital_day
Definition: ERF_Radiation.H:328
real1d_k sfc_emis
Definition: ERF_Radiation.H:377
amrex::Real m_orbital_obliq
Definition: ERF_Radiation.H:336
amrex::Real m_n2vmr
Definition: ERF_Radiation.H:296
real3d_k cld_tau_sw_bnd
Definition: ERF_Radiation.H:450
bool m_fixed_orbital_year
Definition: ERF_Radiation.H:334
real2d_k sw_clnsky_flux_up
Definition: ERF_Radiation.H:416
int m_nlay
Definition: ERF_Radiation.H:311
real3d_k lw_bnd_flux_up
Definition: ERF_Radiation.H:433
int m_nswgpts
Definition: ERF_Radiation.H:352
bool m_first_step
Definition: ERF_Radiation.H:238
real1d_k sfc_alb_dir_vis
Definition: ERF_Radiation.H:367
const std::vector< std::string > m_gas_names
Definition: ERF_Radiation.H:284
real2d_k sw_clrsky_heating
Definition: ERF_Radiation.H:398
bool m_extra_clnclrsky_diag
Definition: ERF_Radiation.H:321
real2d_k t_lay
Definition: ERF_Radiation.H:384
void rad_run_impl()
Definition: ERF_Radiation.H:194
real1d_k sfc_flux_dir_vis
Definition: ERF_Radiation.H:371
real2d_k iwp
Definition: ERF_Radiation.H:395
amrex::Geometry m_geom
Definition: ERF_Radiation.H:227
amrex::MultiFab * m_lat
Definition: ERF_Radiation.H:257
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:391
real1d_k cosine_zenith
Definition: ERF_Radiation.H:365
amrex::Vector< int > m_col_offsets
Definition: ERF_Radiation.H:314
real1d_k t_sfc
Definition: ERF_Radiation.H:378
void alloc_buffers()
Definition: ERF_Radiation.cpp:193
int m_orbital_year
Definition: ERF_Radiation.H:326
real2d_k sw_clnclrsky_flux_dn
Definition: ERF_Radiation.H:411
amrex::Real m_time
Definition: ERF_Radiation.H:221
Definition: ERF_DataStruct.H:123