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, double &time, const double &dt, const amrex::BoxArray &ba, amrex::Geometry &geom, amrex::MultiFab *cons_in, amrex::iMultiFab *lmask, amrex::MultiFab *t_surf, const amrex::Vector< const amrex::MultiFab * > &lsm_input_ptrs, const amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs, amrex::MultiFab *qheating_rates, amrex::MultiFab *rad_fluxes, amrex::MultiFab *z_phys, amrex::MultiFab *lat_ptr, amrex::MultiFab *lon_ptr, const bool updated_lsm, amrex::Real &declin, amrex::Real &calday) override
 
bool radiation_updated () const override
 Reports whether the last radiation call updated its outputs. More...
 
void set_grids (int &level, int &step, double &time, const double &dt, const amrex::BoxArray &ba, amrex::Geometry &geom, amrex::MultiFab *cons_in, amrex::iMultiFab *lmask, amrex::MultiFab *t_surf, const amrex::Vector< const amrex::MultiFab * > &lsm_input_ptrs, amrex::MultiFab *qheating_rates, amrex::MultiFab *rad_fluxes, amrex::MultiFab *z_phys, amrex::MultiFab *lat, amrex::MultiFab *lon, const bool updated_lsm)
 
void alloc_buffers ()
 
void dealloc_buffers ()
 
void mf_to_kokkos_buffers (amrex::iMultiFab *lmask, amrex::MultiFab *t_surf, const amrex::Vector< const amrex::MultiFab * > &lsm_input_ptrs)
 
void kokkos_buffers_to_mf (const amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
 
void write_rrtmgp_fluxes ()
 
void initialize_impl ()
 
void run_impl ()
 
void finalize_impl (const amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
 
void rad_run_impl (const 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
 
bool is_nested_patch () const override
 
void populateDatalogMF ()
 
virtual void WriteDataLog (const double &time) override
 
- Public Member Functions inherited from IRadiation
virtual ~IRadiation ()=default
 
void setupDataLog ()
 
void setDataLogFrequency (const int nstep)
 
bool hasDatalog ()
 

Private Attributes

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

Additional Inherited Members

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

Constructor & Destructor Documentation

◆ Radiation()

Radiation::Radiation ( const int &  lev,
SolverChoice &  sc 
)
87 {
88  // Note that Kokkos is now initialized in main.cpp
89  m_rdOcp = sc.rdOcp;
90 
91  // Check if we have a valid moisture model
92  if (sc.moisture_type != MoistureType::None) { m_moist = true; }
93 
94  // Cloud-ice support follows the configured moisture-component mapping.
96  m_ice = (m_qi_comp >= 0);
97 
98  // Check if we have a land surface model enabled
99  if (sc.lsm_type != LandSurfaceType::None) { m_lsm = true; }
100 
101  // Construct parser object for following reads
102  ParmParse pp("erf");
103 
104  // Must specify a surface temp (LSM can overwrite)
105  pp.get("rad_t_sfc", m_rad_t_sfc);
106 
107  // Radiation timestep, as a number of atm steps
108  pp.queryAdd("rad_freq_in_steps", m_rad_freq_in_steps);
109 
110  // Get nvar if specified
111  pp.queryAdd("rad_nvar", m_rad_nvar);
113  "erf.rad_nvar must be greater than 0. "
114  "It controls the amount of memory allocated for temporaries with RRTMGP; "
115  "a value of 0 would allocate no memory.");
116 
117  // Number of columns per RRTMGP chunk (controls peak GPU memory)
118  pp.queryAdd("rad_ncol_chunk", m_ncol_chunk_requested);
120  "erf.rad_ncol_chunk must be a positive integer (default 1024). "
121  "It controls the number of columns processed per RRTMGP kernel launch; "
122  "a value of 0 or negative would produce an infinite loop.");
124 
125  // Flag to write fluxes to plt file
126  pp.queryAdd("rad_write_fluxes", m_rad_write_fluxes);
127 
128  // Do MCICA subcolumn sampling
129  pp.queryAdd("rad_do_subcol_sampling", m_do_subcol_sampling);
130 
131  // Determine orbital year. If orbital_year is negative, use current year
132  // from timestamp for orbital year; if non-negative, use provided orbital year
133  // for duration of simulation. Note that this is keyed off the value itself,
134  // not off whether the input was present, so that a negative value defers to
135  // the timestamp as documented.
136  pp.queryAdd("rad_orbital_year", m_orbital_year);
138 
139  // Get orbital parameters from inputs file
140  pp.queryAdd("rad_orbital_eccentricity", m_orbital_eccen);
141  pp.queryAdd("rad_orbital_obliquity" , m_orbital_obliq);
142  pp.queryAdd("rad_orbital_mvelp" , m_orbital_mvelp);
143 
144  // Get a constant lat/lon for idealized simulations
145  pp.queryAdd("rad_cons_lat", m_lat_cons);
146  pp.queryAdd("rad_cons_lon", m_lon_cons);
147 
148  // Both models read these keys and feed them to the same cos-zenith formula in the
149  // same units, so apply the range checks RadChoice::init_params already applies.
150  if (!std::isfinite(m_lat_cons) || m_lat_cons < Real(-90.0) || m_lat_cons > Real(90.0)) {
151  amrex::Abort("erf.rad_cons_lat = " + std::to_string(m_lat_cons) +
152  " must lie in [-90, 90] degrees.");
153  }
154  if (!std::isfinite(m_lon_cons) || m_lon_cons < Real(-180.0) || m_lon_cons > Real(180.0)) {
155  amrex::Abort("erf.rad_cons_lon = " + std::to_string(m_lon_cons) +
156  " must lie in [-180, 180] degrees.");
157  }
158 
159  // Value for prescribing an invariant solar constant (i.e. total solar irradiance at
160  // TOA). Used for idealized experiments such as RCE. Disabled when value is less than zero
161  pp.queryAdd("fixed_total_solar_irradiance", m_fixed_total_solar_irradiance);
162 
163  // Determine whether or not we are using a fixed solar zenith angle (positive value)
164  pp.queryAdd("fixed_solar_zenith_angle", m_fixed_solar_zenith_angle);
165 
166  // The same checks the two-stream model applies to these shared inputs
167  // (RadChoice::init_params): the fixed zenith input is a cosine, and the
168  // surface temperature must be a temperature.
169  if (m_fixed_solar_zenith_angle > Real(1.0)) {
170  amrex::Abort("erf.fixed_solar_zenith_angle = " + std::to_string(m_fixed_solar_zenith_angle) +
171  " is the cosine of the solar zenith angle and cannot exceed 1; 60 degrees is 0.5.");
172  }
173  if (!std::isfinite(m_rad_t_sfc) || m_rad_t_sfc <= Real(0.0)) {
174  amrex::Abort("erf.rad_t_sfc = " + std::to_string(m_rad_t_sfc) +
175  " must be a positive temperature [K].");
176  }
177 
178  // Get prescribed surface values of greenhouse gases
179  pp.queryAdd("co2vmr", m_co2vmr);
180  pp.queryarr("o3vmr" , m_o3vmr );
181  pp.queryAdd("n2ovmr", m_n2ovmr);
182  pp.queryAdd("covmr" , m_covmr );
183  pp.queryAdd("ch4vmr", m_ch4vmr);
184  pp.queryAdd("o2vmr" , m_o2vmr );
185  pp.queryAdd("n2vmr" , m_n2vmr );
186 
187  // Aerosol forcing hook (not implemented). The aerosol arrays that used to be
188  // passed through rrtmgp_main were never populated with real data, so enabling
189  // this flag only ever multiplied radiation by zero aerosol optics. The hook is
190  // kept so a future SPA/prescribed-aerosol scheme can wire in without touching
191  // the ParmParse surface.
192  pp.queryAdd("rad_do_aerosol", m_do_aerosol_rad);
193  if (m_do_aerosol_rad) {
194  amrex::Abort("erf.rad_do_aerosol = true is not supported: aerosol forcing is "
195  "currently not implemented in the ERF RRTMGP interface. The hook "
196  "is retained for a future aerosol coupling; set rad_do_aerosol = "
197  "false (or remove it) to continue.");
198  }
199 
200  // Whether we do extra clean/clear sky calculations
201  pp.queryAdd("rad_extra_clnclrsky_diag", m_extra_clnclrsky_diag);
202  pp.queryAdd("rad_extra_clnsky_diag" , m_extra_clnsky_diag);
203 
204  // Parse path and file names
205  pp.queryAdd("rrtmgp_file_path" , rrtmgp_file_path);
206  pp.queryAdd("rrtmgp_coeffs_sw" , rrtmgp_coeffs_sw );
207  pp.queryAdd("rrtmgp_coeffs_lw" , rrtmgp_coeffs_lw );
208  pp.queryAdd("rrtmgp_cloud_optics_sw", rrtmgp_cloud_optics_sw);
209  pp.queryAdd("rrtmgp_cloud_optics_lw", rrtmgp_cloud_optics_lw);
210 
211  // Fail early, and with a message that names the missing files and the inputs that
212  // control them, rather than letting netCDF report a bare "No such file or directory"
213  // when the coefficients are opened just below (or, for the cloud optics files, much
214  // later inside rrtmgp_initialize). Every rank does this so that all of them abort.
215  check_rrtmgp_data_files(rrtmgp_file_path,
216  {{"rrtmgp_coeffs_sw" , rrtmgp_coeffs_sw },
217  {"rrtmgp_coeffs_lw" , rrtmgp_coeffs_lw },
218  {"rrtmgp_cloud_optics_sw", rrtmgp_cloud_optics_sw},
219  {"rrtmgp_cloud_optics_lw", rrtmgp_cloud_optics_lw}});
220 
221  // Append file names to path
226 
227  // Get dimensions from lookup data
228  if (ParallelDescriptor::IOProcessor()) {
229  auto ncf_sw = ncutils::NCFile::open(rrtmgp_coeffs_file_sw, NC_CLOBBER | NC_NETCDF4);
230  m_nswbands = ncf_sw.dim("bnd").len();
231  m_nswgpts = ncf_sw.dim("gpt").len();
232  ncf_sw.close();
233 
234  auto ncf_lw = ncutils::NCFile::open(rrtmgp_coeffs_file_lw, NC_CLOBBER | NC_NETCDF4);
235  m_nlwbands = ncf_lw.dim("bnd").len();
236  m_nlwgpts = ncf_lw.dim("gpt").len();
237  ncf_lw.close();
238  }
239  int ioproc = ParallelDescriptor::IOProcessorNumber(); // I/O rank
240  ParallelDescriptor::Bcast(&m_nswbands, 1, ioproc);
241  ParallelDescriptor::Bcast(&m_nlwbands, 1, ioproc);
242  ParallelDescriptor::Bcast(&m_nswgpts, 1, ioproc);
243  ParallelDescriptor::Bcast(&m_nlwgpts, 1, ioproc);
244 
245  // Output for user
246  if (lev == 0) {
247  Print() << "Radiation interface constructed:\n";
248  Print() << "========================================================\n";
249  Print() << "Coeff SW file: " << rrtmgp_coeffs_file_sw << "\n";
250  Print() << "Coeff LW file: " << rrtmgp_coeffs_file_lw << "\n";
251  Print() << "Cloud SW file: " << rrtmgp_cloud_optics_file_sw << "\n";
252  Print() << "Cloud LW file: " << rrtmgp_cloud_optics_file_lw << "\n";
253  Print() << "Number of short/longwave bands: "
254  << m_nswbands << " " << m_nlwbands << "\n";
255  Print() << "Number of short/longwave gauss points: "
256  << m_nswgpts << " " << m_nlwgpts << "\n";
257  Print() << "========================================================\n";
258  }
259 }
ParmParse pp("prob")
AMREX_ALWAYS_ASSERT_WITH_MESSAGE(m_cloud_chamber_config.active, "Cloud Chamber: initializer reached without a parsed configuration")
amrex::Real Real
Definition: ERF_ShocInterface.H:19
std::string rrtmgp_coeffs_file_sw
Definition: ERF_Radiation.H:356
std::string rrtmgp_coeffs_sw
Definition: ERF_Radiation.H:352
int m_rad_freq_in_steps
Definition: ERF_Radiation.H:441
amrex::Real m_lon_cons
Definition: ERF_Radiation.H:345
int m_nswbands
Definition: ERF_Radiation.H:435
bool m_do_aerosol_rad
Definition: ERF_Radiation.H:399
std::string rrtmgp_coeffs_lw
Definition: ERF_Radiation.H:353
bool m_moist
Definition: ERF_Radiation.H:305
std::string rrtmgp_cloud_optics_file_lw
Definition: ERF_Radiation.H:359
int m_ncol_chunk_requested
Definition: ERF_Radiation.H:449
bool m_lsm
Definition: ERF_Radiation.H:310
bool m_do_subcol_sampling
Definition: ERF_Radiation.H:456
bool m_rad_write_fluxes
Definition: ERF_Radiation.H:302
int m_qi_comp
Definition: ERF_Radiation.H:307
std::string rrtmgp_cloud_optics_file_sw
Definition: ERF_Radiation.H:358
amrex::Real m_orbital_mvelp
Definition: ERF_Radiation.H:420
amrex::Vector< amrex::Real > m_o3vmr
Definition: ERF_Radiation.H:370
std::string rrtmgp_cloud_optics_sw
Definition: ERF_Radiation.H:354
int m_nlwgpts
Definition: ERF_Radiation.H:438
std::string rrtmgp_cloud_optics_lw
Definition: ERF_Radiation.H:355
amrex::Real m_co2vmr
Definition: ERF_Radiation.H:369
bool m_extra_clnsky_diag
Definition: ERF_Radiation.H:402
amrex::Real m_lat_cons
Definition: ERF_Radiation.H:344
amrex::Real m_fixed_total_solar_irradiance
Definition: ERF_Radiation.H:425
amrex::Real m_o2vmr
Definition: ERF_Radiation.H:374
amrex::Real m_n2ovmr
Definition: ERF_Radiation.H:371
amrex::Real m_rad_t_sfc
Definition: ERF_Radiation.H:324
std::string rrtmgp_file_path
Definition: ERF_Radiation.H:351
int m_nlwbands
Definition: ERF_Radiation.H:436
amrex::Real m_orbital_eccen
Definition: ERF_Radiation.H:418
amrex::Real m_covmr
Definition: ERF_Radiation.H:372
int m_ncol_chunk
Definition: ERF_Radiation.H:450
std::string rrtmgp_coeffs_file_lw
Definition: ERF_Radiation.H:357
amrex::Real m_ch4vmr
Definition: ERF_Radiation.H:373
bool m_ice
Definition: ERF_Radiation.H:306
amrex::Real m_fixed_solar_zenith_angle
Definition: ERF_Radiation.H:429
amrex::Real m_orbital_obliq
Definition: ERF_Radiation.H:419
amrex::Real m_n2vmr
Definition: ERF_Radiation.H:375
bool m_fixed_orbital_year
Definition: ERF_Radiation.H:417
int m_nswgpts
Definition: ERF_Radiation.H:437
bool m_extra_clnclrsky_diag
Definition: ERF_Radiation.H:403
amrex::Real m_rdOcp
Definition: ERF_Radiation.H:325
int m_rad_nvar
Definition: ERF_Radiation.H:453
int m_orbital_year
Definition: ERF_Radiation.H:408
@ None
Definition: ERF_LSM_Types.H:6
int qi
cloud ice
Definition: ERF_DataStruct.H:251
MoistureType moisture_type
Moisture or microphysics model.
Definition: ERF_DataStruct.H:2524
amrex::Real rdOcp
Ratio of dry-air gas constant to c_p.
Definition: ERF_DataStruct.H:2345
LandSurfaceType lsm_type
Land-surface model.
Definition: ERF_DataStruct.H:2535
MoistureComponentIndices moisture_indices
Index map of the moisture data carried by the active scheme: conserved-state components for the speci...
Definition: ERF_DataStruct.H:2575
Here is the call graph for this function:

◆ ~Radiation()

Radiation::~Radiation ( )
inline
55  {
56  // Release k-distribution data and memory pool
57  if (rrtmgp::initialized) {
59  }
60  // Note that Kokkos is now finalized in main.cpp
61  }
void rrtmgp_finalize()
Definition: ERF_RRTMGP_Interface.cpp:275
bool initialized
Definition: ERF_RRTMGP_Interface.cpp:25
Here is the call graph for this function:

Member Function Documentation

◆ alloc_buffers()

void Radiation::alloc_buffers ( )
336 {
337  // 1d size (m_ngas)
338  const Real* mol_weight_gas_p = m_mol_weight_gas.data();
339  const std::string* gas_names_p = m_gas_names.data();
340  m_gas_mol_weights = real1d_k("m_gas_mol_weights", m_ngas);
341  realHost1d_k m_gas_mol_weights_h("m_gas_mol_weights_h", m_ngas);
342  gas_names_offset.clear(); gas_names_offset.resize(m_ngas);
343  std::string* gas_names_offset_p = gas_names_offset.data();
344  Kokkos::parallel_for(Kokkos::RangePolicy<Kokkos::Serial>(0, m_ngas),
345  [&] (int igas)
346  {
347  m_gas_mol_weights_h(igas) = mol_weight_gas_p[igas];
348  gas_names_offset_p[igas] = gas_names_p[igas];
349  });
350  Kokkos::deep_copy(m_gas_mol_weights, m_gas_mol_weights_h);
351 
352  // 1d size (1 or nlay)
353  m_o3_size = m_o3vmr.size();
355  "O3 VMR array must be length 1 or nlay");
356  Real* o3vmr_p = m_o3vmr.data();
357  o3_lay = real1d_k("o3_lay", m_o3_size);
358  realHost1d_k o3_lay_h("o3_lay_h", m_o3_size);
359  Kokkos::parallel_for(Kokkos::RangePolicy<Kokkos::Serial>(0, m_o3_size),
360  [&] (int io3)
361  {
362  o3_lay_h(io3) = o3vmr_p[io3];
363  });
364  Kokkos::deep_copy(o3_lay, o3_lay_h);
365 
366  // 1d size (ncol)
367  mu0 = real1d_k("mu0" , m_ncol);
368  sfc_alb_dir_vis = real1d_k("sfc_alb_dir_vis" , m_ncol);
369  sfc_alb_dir_nir = real1d_k("sfc_alb_dir_nir" , m_ncol);
370  sfc_alb_dif_vis = real1d_k("sfc_alb_dif_vis" , m_ncol);
371  sfc_alb_dif_nir = real1d_k("sfc_alb_dif_nir" , m_ncol);
372  sfc_flux_dir_vis = real1d_k("sfc_flux_dir_vis", m_ncol);
373  sfc_flux_dir_nir = real1d_k("sfc_flux_dir_nir", m_ncol);
374  sfc_flux_dif_vis = real1d_k("sfc_flux_dif_vis", m_ncol);
375  sfc_flux_dif_nir = real1d_k("sfc_flux_dif_nir", m_ncol);
376  lat = real1d_k("lat" , m_ncol);
377  lon = real1d_k("lon" , m_ncol);
378  sfc_emis = real1d_k("sfc_emis" , m_ncol);
379  t_sfc = real1d_k("t_sfc" , m_ncol);
380  lw_src = real1d_k("lw_src" , m_ncol);
381 
382  // 2d size (ncol, nlay)
383  r_lay = real2d_k("r_lay" , m_ncol, m_nlay);
384  p_lay = real2d_k("p_lay" , m_ncol, m_nlay);
385  t_lay = real2d_k("t_lay" , m_ncol, m_nlay);
386  z_del = real2d_k("z_del" , m_ncol, m_nlay);
387  qv_lay = real2d_k("qv" , m_ncol, m_nlay);
388  qc_lay = real2d_k("qc" , m_ncol, m_nlay);
389  qi_lay = real2d_k("qi" , m_ncol, m_nlay);
390  cldfrac_tot = real2d_k("cldfrac_tot" , m_ncol, m_nlay);
391  eff_radius_qc = real2d_k("eff_radius_qc", m_ncol, m_nlay);
392  eff_radius_qi = real2d_k("eff_radius_qi", m_ncol, m_nlay);
393  lwp = real2d_k("lwp" , m_ncol, m_nlay);
394  iwp = real2d_k("iwp" , m_ncol, m_nlay);
395  sw_heating = real2d_k("sw_heating" , m_ncol, m_nlay);
396  lw_heating = real2d_k("lw_heating" , m_ncol, m_nlay);
397  if (datalog_int > 0) {
398  sw_clrsky_heating = real2d_k("sw_clrsky_heating", m_ncol, m_nlay);
399  lw_clrsky_heating = real2d_k("lw_clrsky_heating", m_ncol, m_nlay);
400  }
401 
402  // 2d size (ncol, nlay+1)
403  d_tint = real2d_k("d_tint" , m_ncol, m_nlay+1);
404  p_lev = real2d_k("p_lev" , m_ncol, m_nlay+1);
405  t_lev = real2d_k("t_lev" , m_ncol, m_nlay+1);
406 
407  sw_flux_up = real2d_k("sw_flux_up" , m_ncol, m_nlay+1);
408  sw_flux_dn = real2d_k("sw_flux_dn" , m_ncol, m_nlay+1);
409  sw_flux_dn_dir = real2d_k("sw_flux_dn_dir" , m_ncol, m_nlay+1);
410 
411  lw_flux_up = real2d_k("lw_flux_up" , m_ncol, m_nlay+1);
412  lw_flux_dn = real2d_k("lw_flux_dn" , m_ncol, m_nlay+1);
413 
414  // Clear-sky flux arrays are always needed
415  if (datalog_int > 0) {
416  sw_clrsky_flux_up = real2d_k("sw_clrsky_flux_up" , m_ncol, m_nlay+1);
417  sw_clrsky_flux_dn = real2d_k("sw_clrsky_flux_dn" , m_ncol, m_nlay+1);
418  sw_clrsky_flux_dn_dir = real2d_k("sw_clrsky_flux_dn_dir", m_ncol, m_nlay+1);
419  lw_clrsky_flux_up = real2d_k("lw_clrsky_flux_up" , m_ncol, m_nlay+1);
420  lw_clrsky_flux_dn = real2d_k("lw_clrsky_flux_dn" , m_ncol, m_nlay+1);
421  } else {
422  // Use m_ncol_chunk_requested to prevent pool shrinkage after regrid
423  sw_clrsky_flux_up = real2d_k("sw_clrsky_flux_up" , m_ncol_chunk_requested, m_nlay+1);
424  sw_clrsky_flux_dn = real2d_k("sw_clrsky_flux_dn" , m_ncol_chunk_requested, m_nlay+1);
425  sw_clrsky_flux_dn_dir = real2d_k("sw_clrsky_flux_dn_dir", m_ncol_chunk_requested, m_nlay+1);
426  lw_clrsky_flux_up = real2d_k("lw_clrsky_flux_up" , m_ncol_chunk_requested, m_nlay+1);
427  lw_clrsky_flux_dn = real2d_k("lw_clrsky_flux_dn" , m_ncol_chunk_requested, m_nlay+1);
428  }
429 
430  // Clean-clear-sky diagnostic fluxes (only when enabled)
432  sw_clnclrsky_flux_up = real2d_k("sw_clnclrsky_flux_up" , m_ncol, m_nlay+1);
433  sw_clnclrsky_flux_dn = real2d_k("sw_clnclrsky_flux_dn" , m_ncol, m_nlay+1);
434  sw_clnclrsky_flux_dn_dir = real2d_k("sw_clnclrsky_flux_dn_dir", m_ncol, m_nlay+1);
435  lw_clnclrsky_flux_up = real2d_k("lw_clnclrsky_flux_up" , m_ncol, m_nlay+1);
436  lw_clnclrsky_flux_dn = real2d_k("lw_clnclrsky_flux_dn" , m_ncol, m_nlay+1);
437  } else {
438  sw_clnclrsky_flux_up = real2d_k("sw_clnclrsky_flux_up" , 1, 1);
439  sw_clnclrsky_flux_dn = real2d_k("sw_clnclrsky_flux_dn" , 1, 1);
440  sw_clnclrsky_flux_dn_dir = real2d_k("sw_clnclrsky_flux_dn_dir", 1, 1);
441  lw_clnclrsky_flux_up = real2d_k("lw_clnclrsky_flux_up" , 1, 1);
442  lw_clnclrsky_flux_dn = real2d_k("lw_clnclrsky_flux_dn" , 1, 1);
443  }
444 
445  // Clean-sky diagnostic fluxes (only when enabled)
446  if (m_extra_clnsky_diag) {
447  sw_clnsky_flux_up = real2d_k("sw_clnsky_flux_up" , m_ncol, m_nlay+1);
448  sw_clnsky_flux_dn = real2d_k("sw_clnsky_flux_dn" , m_ncol, m_nlay+1);
449  sw_clnsky_flux_dn_dir = real2d_k("sw_clnsky_flux_dn_dir" , m_ncol, m_nlay+1);
450  lw_clnsky_flux_up = real2d_k("lw_clnsky_flux_up" , m_ncol, m_nlay+1);
451  lw_clnsky_flux_dn = real2d_k("lw_clnsky_flux_dn" , m_ncol, m_nlay+1);
452  } else {
453  sw_clnsky_flux_up = real2d_k("sw_clnsky_flux_up" , 1, 1);
454  sw_clnsky_flux_dn = real2d_k("sw_clnsky_flux_dn" , 1, 1);
455  sw_clnsky_flux_dn_dir = real2d_k("sw_clnsky_flux_dn_dir" , 1, 1);
456  lw_clnsky_flux_up = real2d_k("lw_clnsky_flux_up" , 1, 1);
457  lw_clnsky_flux_dn = real2d_k("lw_clnsky_flux_dn" , 1, 1);
458  }
459 
460  // 3d size (ncol_chunk, nlay+1, nswbands)
461  // Use m_ncol_chunk_requested to prevent pool shrinkage after regrid
466 
467  // 3d size (ncol_chunk, nlay+1, nlwbands)
468  // Use m_ncol_chunk_requested to prevent pool shrinkage after regrid
471 
472  // 2d size (ncol, nswbands)
473  sfc_alb_dir = real2d_k("sfc_alb_dir", m_ncol, m_nswbands);
474  sfc_alb_dif = real2d_k("sfc_alb_dif", m_ncol, m_nswbands);
475 
476  // Aerosol optical properties — allocated only when aerosol coupling is on.
477  // The flag gates allocation so today (coupling not implemented, abort fires
478  // in the constructor) these stay as empty Views and cost nothing. When a
479  // future aerosol scheme populates them, hook up the plumbing into
480  // rrtmgp_main as well.
481  if (m_do_aerosol_rad) {
482  aero_tau_sw = real3d_k("aero_tau_sw", m_ncol, m_nlay, m_nswbands);
483  aero_ssa_sw = real3d_k("aero_ssa_sw", m_ncol, m_nlay, m_nswbands);
484  aero_g_sw = real3d_k("aero_g_sw", m_ncol, m_nlay, m_nswbands);
485  aero_tau_lw = real3d_k("aero_tau_lw", m_ncol, m_nlay, m_nlwbands);
486  }
487 }
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
int datalog_int
Definition: ERF_RadiationInterface.H:105
real3d_k sw_bnd_flux_dn
Definition: ERF_Radiation.H:522
real2d_k lw_flux_up
Definition: ERF_Radiation.H:502
real3d_k aero_tau_sw
Definition: ERF_Radiation.H:541
int m_o3_size
Definition: ERF_Radiation.H:379
real3d_k sw_bnd_flux_dir
Definition: ERF_Radiation.H:523
real2d_k sw_clnsky_flux_dn
Definition: ERF_Radiation.H:511
real2d_k d_tint
Definition: ERF_Radiation.H:496
real2d_k lw_clnclrsky_flux_dn
Definition: ERF_Radiation.H:514
real2d_k lwp
Definition: ERF_Radiation.H:488
real1d_k lw_src
Definition: ERF_Radiation.H:475
real2d_k eff_radius_qi
Definition: ERF_Radiation.H:487
real1d_k m_gas_mol_weights
Definition: ERF_Radiation.H:380
real2d_k sw_heating
Definition: ERF_Radiation.H:490
real3d_k lw_bnd_flux_dn
Definition: ERF_Radiation.H:528
real3d_k sw_bnd_flux_up
Definition: ERF_Radiation.H:521
real2d_k qv_lay
Definition: ERF_Radiation.H:482
real2d_k sw_clnclrsky_flux_dn_dir
Definition: ERF_Radiation.H:506
real1d_k sfc_flux_dif_vis
Definition: ERF_Radiation.H:469
real2d_k lw_clnclrsky_flux_up
Definition: ERF_Radiation.H:513
real1d_k lat
Definition: ERF_Radiation.H:471
real2d_k qi_lay
Definition: ERF_Radiation.H:484
real2d_k sw_clrsky_flux_up
Definition: ERF_Radiation.H:507
real2d_k t_lev
Definition: ERF_Radiation.H:498
real1d_k o3_lay
Definition: ERF_Radiation.H:459
real1d_k sfc_alb_dif_vis
Definition: ERF_Radiation.H:465
real1d_k mu0
Definition: ERF_Radiation.H:462
real2d_k cldfrac_tot
Definition: ERF_Radiation.H:485
real1d_k sfc_flux_dir_nir
Definition: ERF_Radiation.H:468
real2d_k sw_clnclrsky_flux_up
Definition: ERF_Radiation.H:504
real3d_k aero_g_sw
Definition: ERF_Radiation.H:543
real2d_k sw_flux_up
Definition: ERF_Radiation.H:499
real3d_k sw_bnd_flux_dif
Definition: ERF_Radiation.H:524
real2d_k sfc_alb_dif
Definition: ERF_Radiation.H:532
real1d_k sfc_alb_dif_nir
Definition: ERF_Radiation.H:466
real2d_k lw_clnsky_flux_dn
Definition: ERF_Radiation.H:518
real2d_k p_lay
Definition: ERF_Radiation.H:479
real2d_k r_lay
Definition: ERF_Radiation.H:478
int m_ncol
Definition: ERF_Radiation.H:389
real2d_k sw_flux_dn_dir
Definition: ERF_Radiation.H:501
const std::vector< amrex::Real > m_mol_weight_gas
Definition: ERF_Radiation.H:365
std::vector< std::string > gas_names_offset
Definition: ERF_Radiation.H:381
real2d_k sw_clrsky_flux_dn
Definition: ERF_Radiation.H:508
real3d_k aero_ssa_sw
Definition: ERF_Radiation.H:542
real2d_k sw_clnsky_flux_dn_dir
Definition: ERF_Radiation.H:512
real2d_k qc_lay
Definition: ERF_Radiation.H:483
real2d_k lw_clrsky_flux_up
Definition: ERF_Radiation.H:515
real2d_k sw_flux_dn
Definition: ERF_Radiation.H:500
real2d_k lw_clrsky_flux_dn
Definition: ERF_Radiation.H:516
real2d_k z_del
Definition: ERF_Radiation.H:481
real3d_k aero_tau_lw
Definition: ERF_Radiation.H:544
real2d_k lw_flux_dn
Definition: ERF_Radiation.H:503
real2d_k sw_clrsky_flux_dn_dir
Definition: ERF_Radiation.H:509
real1d_k sfc_flux_dif_nir
Definition: ERF_Radiation.H:470
int m_ngas
Definition: ERF_Radiation.H:362
real2d_k lw_heating
Definition: ERF_Radiation.H:491
real1d_k sfc_alb_dir_nir
Definition: ERF_Radiation.H:464
real2d_k lw_clnsky_flux_up
Definition: ERF_Radiation.H:517
real2d_k sfc_alb_dir
Definition: ERF_Radiation.H:531
real1d_k lon
Definition: ERF_Radiation.H:472
real1d_k sfc_emis
Definition: ERF_Radiation.H:473
real2d_k sw_clnsky_flux_up
Definition: ERF_Radiation.H:510
int m_nlay
Definition: ERF_Radiation.H:390
real3d_k lw_bnd_flux_up
Definition: ERF_Radiation.H:527
real1d_k sfc_alb_dir_vis
Definition: ERF_Radiation.H:463
const std::vector< std::string > m_gas_names
Definition: ERF_Radiation.H:363
real2d_k sw_clrsky_heating
Definition: ERF_Radiation.H:492
real2d_k t_lay
Definition: ERF_Radiation.H:480
real1d_k sfc_flux_dir_vis
Definition: ERF_Radiation.H:467
real2d_k iwp
Definition: ERF_Radiation.H:489
real2d_k p_lev
Definition: ERF_Radiation.H:497
real2d_k lw_clrsky_heating
Definition: ERF_Radiation.H:493
real2d_k eff_radius_qc
Definition: ERF_Radiation.H:486
real1d_k t_sfc
Definition: ERF_Radiation.H:474
real2d_k sw_clnclrsky_flux_dn
Definition: ERF_Radiation.H:505
Here is the call graph for this function:

◆ dealloc_buffers()

void Radiation::dealloc_buffers ( )
491 {
492  // 1d size (m_ngas)
494 
495  // 1d size (1 or nlay)
496  o3_lay = real1d_k();
497 
498  // 1d size (ncol)
499  mu0 = real1d_k();
508  lat = real1d_k();
509  lon = real1d_k();
510  sfc_emis = real1d_k();
511  t_sfc = real1d_k();
512  lw_src = real1d_k();
513 
514  // 2d size (ncol, nlay)
515  r_lay = real2d_k();
516  p_lay = real2d_k();
517  t_lay = real2d_k();
518  z_del = real2d_k();
519  qv_lay = real2d_k();
520  qc_lay = real2d_k();
521  qi_lay = real2d_k();
522  cldfrac_tot = real2d_k();
525  lwp = real2d_k();
526  iwp = real2d_k();
527  sw_heating = real2d_k();
528  lw_heating = real2d_k();
531 
532  // 2d size (ncol, nlay+1)
533  d_tint = real2d_k();
534  p_lev = real2d_k();
535  t_lev = real2d_k();
536  sw_flux_up = real2d_k();
537  sw_flux_dn = real2d_k();
539  lw_flux_up = real2d_k();
540  lw_flux_dn = real2d_k();
556 
557  // 3d size (ncol, nlay+1, nswbands)
562 
563  // 3d size (ncol, nlay+1, nlwbands)
566 
567  // 2d size (ncol, nswbands)
568  sfc_alb_dir = real2d_k();
569  sfc_alb_dif = real2d_k();
570 
571  // Aerosol scaffolding (no-op unless m_do_aerosol_rad enabled allocation above)
572  aero_tau_sw = real3d_k();
573  aero_ssa_sw = real3d_k();
574  aero_g_sw = real3d_k();
575  aero_tau_lw = real3d_k();
576 }

◆ finalize_impl()

void Radiation::finalize_impl ( const amrex::Vector< amrex::MultiFab * > &  lsm_output_ptrs)
1528 {
1529  // Reset gas concentrations (k-dist data persists across steps)
1530  m_gas_concs.reset();
1531 
1532  // Fill the AMReX MFs from Kokkos Views
1533  kokkos_buffers_to_mf(lsm_output_ptrs);
1534 
1535  // Write fluxes if requested
1537 
1538  // Fill output data for datalog before deallocating
1539  if (datalog_int > 0) {
1542  Kokkos::fence();
1544  }
1545 
1546  // Deallocate the buffer arrays
1547  dealloc_buffers();
1548 }
void dealloc_buffers()
Definition: ERF_Radiation.cpp:490
void populateDatalogMF()
Definition: ERF_Radiation.cpp:964
void kokkos_buffers_to_mf(const amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
Definition: ERF_Radiation.cpp:824
GasConcsK< amrex::Real, layout_t, KokkosDefaultDevice > m_gas_concs
Definition: ERF_Radiation.H:383
void write_rrtmgp_fluxes()
Definition: ERF_Radiation.cpp:924
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.

245  {
246  return m_lsm_input_names;
247  }
amrex::Vector< std::string > m_lsm_input_names
Definition: ERF_Radiation.H:313

◆ get_lsm_output_varnames()

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

Reimplemented from IRadiation.

253  {
254  return m_lsm_output_names;
255  }
amrex::Vector< std::string > m_lsm_output_names
Definition: ERF_Radiation.H:318

◆ Init()

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

Implements IRadiation.

68  {
69  // Check for nested patches: fine levels that don't reach the model top
70  // These need special handling since radiation requires complete atmospheric columns
71  int klo = geom.Domain().smallEnd(2);
72  int khi = geom.Domain().bigEnd(2);
73 
74  // Detect nested patches: a nested fine level doesn't span the full domain vertically.
75  // Use ba.minimalBox() which gives the bounding box of all grids on this level.
76  // If that box doesn't reach the domain top or bottom, it's a nested patch.
77  amrex::Box minimal_box = ba.minimalBox();
78  bool is_nested_patch = (minimal_box.smallEnd(2) > klo) || (minimal_box.bigEnd(2) < khi);
79 
80  // Update the nested patch flag - must be set/cleared on every Init call
81  // so that levels that stop being nested after regrid will resume radiation computation
83 
84  // Skip radiation on nested patches - heating rates will be interpolated from coarser level
85  if (is_nested_patch) {
86  return;
87  }
88 
89  // Reset vector of offsets for columnar data
90  m_nlay = geom.Domain().length(2);
91 
92  m_ncol = 0;
93  m_col_offsets.clear();
94  m_col_offsets.resize(int(ba.size()));
95  for (amrex::MFIter mfi(*cons_in); mfi.isValid(); ++mfi) {
96  const amrex::Box& vbx = mfi.validbox();
97  // This assertion catches true vertical MPI decomposition (individual boxes that don't
98  // span the full domain height). The nested patch case is handled above by returning early.
99  AMREX_ALWAYS_ASSERT_WITH_MESSAGE((klo == vbx.smallEnd(2)) &&
100  (khi == vbx.bigEnd(2)),
101  "Vertical decomposition with radiation is not allowed.");
102  int nx = vbx.length(0);
103  int ny = vbx.length(1);
104  m_col_offsets[mfi.index()] = m_ncol;
105  m_ncol += nx * ny;
106  }
107 
108  // Recompute the effective chunk size from the user's request rather than
109  // clamping the stored value in place. Init() is called again for every
110  // MakeNewLevelFromScratch / MakeNewLevelFromCoarse / RemakeLevel on the same
111  // Radiation object, so an in-place clamp could only ever shrink: a rank that
112  // owned no boxes at some Init (m_ncol == 0) would pin the chunk size at zero
113  // and then spin forever in the run_impl chunk loop once a regrid handed it
114  // real boxes. Deriving from the immutable request instead lets the chunk size
115  // grow back. It is zero only when m_ncol is zero, and a rank with no columns
116  // does no radiation work.
118 
119  // NOTE: the buffers are deliberately NOT allocated here. set_grids() allocates them
120  // at the top of every Run() that actually updates radiation, and finalize_impl()
121  // frees them again at the end of it, which is what keeps only one level's worth
122  // of RRTMGP working memory resident at a time (see RRTMGP_Memory_Reduction.md).
123  // Allocating here as well would hold every level's buffers from the moment that
124  // level is created until its first Run, so the whole hierarchy would be live
125  // simultaneously at startup -- the peak the chunking exists to avoid.
126  };
const int nx
Definition: ERF_InitCustomPertVels_CloudChamber.H:14
const int ny
Definition: ERF_InitCustomPertVels_CloudChamber.H:15
const int klo
Definition: ERF_InitCustomPert_ABL.H:75
const int khi
Definition: ERF_InitCustomPert_Bubble.H:21
bool is_nested_patch() const override
Definition: ERF_Radiation.H:258
bool m_is_nested_patch
Definition: ERF_Radiation.H:278
amrex::Vector< int > m_col_offsets
Definition: ERF_Radiation.H:393
Here is the call graph for this function:

◆ initialize_impl()

void Radiation::initialize_impl ( )
1206 {
1207  // Initialize gas concentrations for this step
1209 
1210  // Load k-distribution and cloud optics data only once.
1211  // These are static lookup tables that never change.
1212  // Size the memory pool for the requested chunk size (not the effective one, and
1213  // not min with the current m_ncol) so that the pool remains valid even if m_ncol
1214  // grows after regridding/load balancing. The pool is created once and never
1215  // resized, whereas the effective chunk size is recomputed at every Init() and is
1216  // bounded above by the request.
1217  if (!rrtmgp::initialized) {
1218  gas_concs_t gas_concs_pool;
1219  gas_concs_pool.init(gas_names_offset, m_ncol_chunk_requested, m_nlay);
1220  rrtmgp::rrtmgp_initialize(gas_concs_pool,
1223  m_rad_nvar);
1224  gas_concs_pool.reset();
1225  }
1226 }
GasConcsK< RealT, layout_t, KokkosDefaultDevice > gas_concs_t
Definition: ERF_RRTMGP_Interface.H:31
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, const int &nvar)
Definition: ERF_RRTMGP_Interface.cpp:235

Referenced by rad_run_impl().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ is_nested_patch()

bool Radiation::is_nested_patch ( ) const
inlineoverridevirtual

Reimplemented from IRadiation.

259  {
260  return m_is_nested_patch;
261  }

Referenced by Init().

Here is the caller graph for this function:

◆ kokkos_buffers_to_mf()

void Radiation::kokkos_buffers_to_mf ( const amrex::Vector< amrex::MultiFab * > &  lsm_output_ptrs)
825 {
826  // The heating rates and fluxes read below were written by Kokkos kernels in run_impl().
827  // Kokkos launches on its own default instance while the ParallelFors here launch on
828  // amrex::Gpu::gpuStream(), and nothing orders the consumer against the producer: MFIter
829  // synchronizes the streams it used, but only after its own kernels, which is too late.
830  // This is the same cross-stream race that made radqrclw vary run to run in the datalog
831  // path (see section 7 of RRTMGP_Memory_Reduction.md); that fence guarded one call site,
832  // this one guards the hand-off that feeds qheating_rates and rad_fluxes back into the
833  // solution. No-op on a CPU build.
834  Kokkos::fence();
835 
836  // Heating rate, fluxes, zenith, lsm ptrs
837 
838  Table2D<Real,Order::C> p_lay_tab(p_lay.data(), {0,0}, {static_cast<int>(p_lay.extent(0)),static_cast<int>(p_lay.extent(1))});
839  Table2D<Real,Order::C> sw_heating_tab(sw_heating.data(), {0,0}, {static_cast<int>(sw_heating.extent(0)),static_cast<int>(sw_heating.extent(1))});
840  Table2D<Real,Order::C> lw_heating_tab(lw_heating.data(), {0,0}, {static_cast<int>(lw_heating.extent(0)),static_cast<int>(lw_heating.extent(1))});
841  Table2D<Real,Order::C> sw_flux_up_tab(sw_flux_up.data(), {0,0}, {static_cast<int>(sw_flux_up.extent(0)),static_cast<int>(sw_flux_up.extent(1))});
842  Table2D<Real,Order::C> sw_flux_dn_tab(sw_flux_dn.data(), {0,0}, {static_cast<int>(sw_flux_dn.extent(0)),static_cast<int>(sw_flux_dn.extent(1))});
843  Table2D<Real,Order::C> lw_flux_up_tab(lw_flux_up.data(), {0,0}, {static_cast<int>(lw_flux_up.extent(0)),static_cast<int>(lw_flux_up.extent(1))});
844  Table2D<Real,Order::C> lw_flux_dn_tab(lw_flux_dn.data(), {0,0}, {static_cast<int>(lw_flux_dn.extent(0)),static_cast<int>(lw_flux_dn.extent(1))});
845 
846  TableData<Real,1> sfc_flux_sw_dn; sfc_flux_sw_dn.resize({0}, {static_cast<int>(sw_flux_dn.extent(0))});
847  TableData<Real,1> sfc_flux_lw_dn; sfc_flux_lw_dn.resize({0}, {static_cast<int>(lw_flux_dn.extent(0))});
848  Table1D<Real> sfc_flux_sw_dn_tab = sfc_flux_sw_dn.table();
849  Table1D<Real> sfc_flux_lw_dn_tab = sfc_flux_lw_dn.table();
850  Table1D<Real> sfc_flux_sw_dir_vis_tab(sfc_flux_dir_vis.data(), {0}, {static_cast<int>(sfc_flux_dir_vis.extent(0))});
851  Table1D<Real> sfc_flux_sw_dir_nir_tab(sfc_flux_dir_nir.data(), {0}, {static_cast<int>(sfc_flux_dir_nir.extent(0))});
852  Table1D<Real> sfc_flux_sw_dif_vis_tab(sfc_flux_dif_vis.data(), {0}, {static_cast<int>(sfc_flux_dif_vis.extent(0))});
853  Table1D<Real> sfc_flux_sw_dif_nir_tab(sfc_flux_dif_nir.data(), {0}, {static_cast<int>(sfc_flux_dif_nir.extent(0))});
854  Table1D<Real> mu0_tab(mu0.data(), {0}, {static_cast<int>(mu0.extent(0))});
855  Vector<Table1D<Real>> rrtmgp_out_vars = {mu0_tab , sfc_flux_sw_dn_tab ,
856  sfc_flux_sw_dir_vis_tab, sfc_flux_sw_dir_nir_tab,
857  sfc_flux_sw_dif_vis_tab, sfc_flux_sw_dif_nir_tab,
858  sfc_flux_lw_dn_tab };
859 
860  const Real rdOcp = m_rdOcp;
861  for (MFIter mfi(*m_cons_in); mfi.isValid(); ++mfi) {
862  const auto& vbx = mfi.validbox();
863  const auto& sbx = makeSlab(vbx,2,vbx.smallEnd(2));
864  const int nx = vbx.length(0);
865  const int imin = vbx.smallEnd(0);
866  const int jmin = vbx.smallEnd(1);
867  const int ktop = vbx.bigEnd(2);
868  const int offset = m_col_offsets[mfi.index()];
869  const Array4<Real>& q_arr = m_qheating_rates->array(mfi);
870  const Array4<Real>& f_arr = m_rad_fluxes->array(mfi);
871  ParallelFor(vbx, [=]
872  AMREX_GPU_DEVICE (int i, int j, int k)
873  {
874  // map [i,j,k] 0-based to [icol, ilay] 0-based
875  const int icol = (j-jmin)*nx + (i-imin) + offset;
876  const int ilay = k;
877 
878  // Temperature heating rate for SW and LW
879  q_arr(i,j,k,0) = sw_heating_tab(icol,ilay);
880  q_arr(i,j,k,1) = lw_heating_tab(icol,ilay);
881 
882  // Convert the dT/dz to dTheta/dz
883  Real iexner = rrtmgp::inverse_exner(Real(p_lay_tab(icol,ilay)), rdOcp);
884  q_arr(i,j,k,0) *= iexner;
885  q_arr(i,j,k,1) *= iexner;
886 
887  // Populate the fluxes: level ilay is the lower interface of
888  // layer k, and the top-of-atmosphere level nlay goes into the
889  // z-ghost cell above the top layer (rad_fluxes has one).
890  f_arr(i,j,k,0) = sw_flux_up_tab(icol,ilay);
891  f_arr(i,j,k,1) = sw_flux_dn_tab(icol,ilay);
892  f_arr(i,j,k,2) = lw_flux_up_tab(icol,ilay);
893  f_arr(i,j,k,3) = lw_flux_dn_tab(icol,ilay);
894  if (k == ktop) {
895  f_arr(i,j,k+1,0) = sw_flux_up_tab(icol,ilay+1);
896  f_arr(i,j,k+1,1) = sw_flux_dn_tab(icol,ilay+1);
897  f_arr(i,j,k+1,2) = lw_flux_up_tab(icol,ilay+1);
898  f_arr(i,j,k+1,3) = lw_flux_dn_tab(icol,ilay+1);
899  }
900 
901  if (k==0) {
902  sfc_flux_sw_dn_tab(icol) = sw_flux_dn_tab(icol,ilay);
903  sfc_flux_lw_dn_tab(icol) = lw_flux_dn_tab(icol,ilay);
904  }
905  });
906  for (int ivar(0); ivar<lsm_output_ptrs.size(); ivar++) {
907  if (lsm_output_ptrs[ivar]) {
908  auto rrtmgp_for_fill = rrtmgp_out_vars[ivar];
909  const Array4<Real>& lsm_out_arr = lsm_output_ptrs[ivar]->array(mfi);
910  ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
911  {
912  // map [i,j,k] 0-based to [icol, ilay] 0-based
913  const int icol = (j-jmin)*nx + (i-imin) + offset;
914 
915  // export the desired variable at surface
916  lsm_out_arr(i,j,k) = rrtmgp_for_fill(icol);
917  });
918  } // valid ptr
919  } // ivar
920  }// mfi
921 }
const Real rdOcp
Definition: ERF_InitCustomPert_ABL.H:72
ParallelFor(fab_box, [=] AMREX_GPU_DEVICE(int i, int j, int k) { qrcuten_arr(i, j, k)=Real(0);qscuten_arr(i, j, k)=Real(0);qicuten_arr(i, j, k)=Real(0);})
AMREX_FORCE_INLINE IntVect offset(const int face_dir, const int normal)
Definition: ERF_ReadBndryPlanes.cpp:32
amrex::MultiFab * m_cons_in
Definition: ERF_Radiation.H:328
amrex::MultiFab * m_rad_fluxes
Definition: ERF_Radiation.H:334
amrex::MultiFab * m_qheating_rates
Definition: ERF_Radiation.H:331
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real inverse_exner(const amrex::Real pressure, const amrex::Real rdOcp)
Definition: ERF_RRTMGP_SurfaceTemperature.H:12
Here is the call graph for this function:

◆ mf_to_kokkos_buffers()

void Radiation::mf_to_kokkos_buffers ( amrex::iMultiFab *  lmask,
amrex::MultiFab *  t_surf,
const amrex::Vector< const amrex::MultiFab * > &  lsm_input_ptrs 
)
583 {
584  Table2D<Real,Order::C> r_lay_tab(r_lay.data(), {0,0}, {static_cast<int>(r_lay.extent(0)),static_cast<int>(r_lay.extent(1))});
585  Table2D<Real,Order::C> p_lay_tab(p_lay.data(), {0,0}, {static_cast<int>(p_lay.extent(0)),static_cast<int>(p_lay.extent(1))});
586  Table2D<Real,Order::C> t_lay_tab(t_lay.data(), {0,0}, {static_cast<int>(t_lay.extent(0)),static_cast<int>(t_lay.extent(1))});
587  Table2D<Real,Order::C> z_del_tab(z_del.data(), {0,0}, {static_cast<int>(z_del.extent(0)),static_cast<int>(z_del.extent(1))});
588  Table2D<Real,Order::C> qv_lay_tab(qv_lay.data(), {0,0}, {static_cast<int>(qv_lay.extent(0)),static_cast<int>(qv_lay.extent(1))});
589  Table2D<Real,Order::C> qc_lay_tab(qc_lay.data(), {0,0}, {static_cast<int>(qc_lay.extent(0)),static_cast<int>(qc_lay.extent(1))});
590  Table2D<Real,Order::C> qi_lay_tab(qi_lay.data(), {0,0}, {static_cast<int>(qi_lay.extent(0)),static_cast<int>(qi_lay.extent(1))});
591  Table2D<Real,Order::C> cldfrac_tot_tab(cldfrac_tot.data(), {0,0}, {static_cast<int>(cldfrac_tot.extent(0)),static_cast<int>(cldfrac_tot.extent(1))});
592 
593  Table2D<Real,Order::C> lwp_tab(lwp.data(), {0,0}, {static_cast<int>(lwp.extent(0)),static_cast<int>(lwp.extent(1))});
594  Table2D<Real,Order::C> iwp_tab(iwp.data(), {0,0}, {static_cast<int>(iwp.extent(0)),static_cast<int>(iwp.extent(1))});
595  Table2D<Real,Order::C> eff_radius_qc_tab(eff_radius_qc.data(), {0,0}, {static_cast<int>(eff_radius_qc.extent(0)),static_cast<int>(eff_radius_qc.extent(1))});
596  Table2D<Real,Order::C> eff_radius_qi_tab(eff_radius_qi.data(), {0,0}, {static_cast<int>(eff_radius_qi.extent(0)),static_cast<int>(eff_radius_qi.extent(1))});
597 
598  Table2D<Real,Order::C> p_lev_tab(p_lev.data(), {0,0}, {static_cast<int>(p_lev.extent(0)),static_cast<int>(p_lev.extent(1))});
599  Table2D<Real,Order::C> t_lev_tab(t_lev.data(), {0,0}, {static_cast<int>(t_lev.extent(0)),static_cast<int>(t_lev.extent(1))});
600 
601  Table1D<Real> lat_tab(lat.data(), {0}, {static_cast<int>(lat.extent(0))});
602  Table1D<Real> lon_tab(lon.data(), {0}, {static_cast<int>(lon.extent(0))});
603  Table1D<Real> t_sfc_tab(t_sfc.data(), {0}, {static_cast<int>(t_sfc.extent(0))});
604 
605  bool moist = m_moist;
606  bool ice = m_ice;
607  const int qi_comp = m_qi_comp;
608  const bool has_lsm = m_lsm;
609  const bool has_lat = m_lat;
610  const bool has_lon = m_lon;
611  const bool has_surflayer = (t_surf);
612  int ncol = m_ncol;
613  int nlay = m_nlay;
614  Real dz = m_geom.CellSize(2);
615  Real cons_lat = m_lat_cons;
616  Real cons_lon = m_lon_cons;
617  Real rad_t_sfc = m_rad_t_sfc;
618  Real rdOcp = m_rdOcp;
619 
620  for (MFIter mfi(*m_cons_in); mfi.isValid(); ++mfi) {
621  const auto& vbx = mfi.validbox();
622  const int nx = vbx.length(0);
623  const int imin = vbx.smallEnd(0);
624  const int jmin = vbx.smallEnd(1);
625  const int offset = m_col_offsets[mfi.index()];
626  const Array4<const Real>& cons_arr = m_cons_in->const_array(mfi);
627  const Array4<const Real>& z_arr = (m_z_phys) ? m_z_phys->const_array(mfi) :
628  Array4<const Real>{};
629  const Array4<const Real>& lat_arr = (m_lat) ? m_lat->const_array(mfi) :
630  Array4<const Real>{};
631  const Array4<const Real>& lon_arr = (m_lon) ? m_lon->const_array(mfi) :
632  Array4<const Real>{};
633  ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
634  {
635  // map [i,j,k] 0-based to [icol, ilay] 0-based
636  const int icol = (j-jmin)*nx + (i-imin) + offset;
637  const int ilay = k;
638 
639  // EOS input (at CC)
640  Real r = cons_arr(i,j,k,Rho_comp);
641  Real rt = cons_arr(i,j,k,RhoTheta_comp);
642  Real qv = (moist) ? std::max(cons_arr(i,j,k,RhoQ1_comp)/r,Real(0.)) : Real(0.);
643  Real qc = (moist) ? std::max(cons_arr(i,j,k,RhoQ2_comp)/r,Real(0.)) : Real(0.);
644  Real qi = (ice && qi_comp >= 0) ? std::max(cons_arr(i,j,k,qi_comp)/r,Real(0.)) : Real(0.);
645 
646  // EOS avg to z-face
647  Real r_lo = cons_arr(i,j,k-1,Rho_comp);
648  Real rt_lo = cons_arr(i,j,k-1,RhoTheta_comp);
649  Real qv_lo = (moist) ? cons_arr(i,j,k-1,RhoQ1_comp)/r_lo : Real(0.);
650  Real dz_k = (z_arr) ? Real(0.125) * ( (z_arr(i ,j ,k+1) - z_arr(i ,j ,k))
651  + (z_arr(i+1,j ,k+1) - z_arr(i+1,j ,k))
652  + (z_arr(i ,j+1,k+1) - z_arr(i ,j+1,k))
653  + (z_arr(i+1,j+1,k+1) - z_arr(i+1,j+1,k)) ) : Real(0.5)*dz; // Dist from w-face to CC at k
654  Real dz_km1 = (z_arr) ? Real(0.125) * ( (z_arr(i ,j ,k ) - z_arr(i ,j ,k-1))
655  + (z_arr(i+1,j ,k ) - z_arr(i+1,j ,k-1))
656  + (z_arr(i ,j+1,k ) - z_arr(i ,j+1,k-1))
657  + (z_arr(i+1,j+1,k ) - z_arr(i+1,j+1,k-1)) ) : Real(0.5)*dz; // Dist from w-face to CC at k-1
658  // NOTE: Linear interpolation to the w-face weights each CC value by the
659  // distance from the face to the *opposite* CC (inverse distance)
660  Real r_avg = (dz_km1*r + dz_k*r_lo ) / (dz_k + dz_km1);
661  Real rt_avg = (dz_km1*rt + dz_k*rt_lo) / (dz_k + dz_km1);
662  Real qv_avg = (dz_km1*qv + dz_k*qv_lo) / (dz_k + dz_km1);
663 
664  // Views at CC
665  r_lay_tab(icol,ilay) = r;
666 
667  p_lay_tab(icol,ilay) = getPgivenRTh(rt, qv);
668  t_lay_tab(icol,ilay) = getTgivenRandRTh(r, rt, qv);
669  z_del_tab(icol,ilay) = (z_arr) ? Real(0.25) * ( (z_arr(i ,j ,k+1) - z_arr(i ,j ,k))
670  + (z_arr(i+1,j ,k+1) - z_arr(i+1,j ,k))
671  + (z_arr(i ,j+1,k+1) - z_arr(i ,j+1,k))
672  + (z_arr(i+1,j+1,k+1) - z_arr(i+1,j+1,k)) ) : dz;
673  qv_lay_tab(icol,ilay) = qv;
674  qc_lay_tab(icol,ilay) = qc;
675  qi_lay_tab(icol,ilay) = qi;
676  cldfrac_tot_tab(icol,ilay) = ((qc+qi)>Real(0.)) ? Real(1.) : Real(0.);
677 
678  // NOTE: These are populated in 'mixing_ratio_to_cloud_mass'
679  lwp_tab(icol,ilay) = Real(0.);
680  iwp_tab(icol,ilay) = Real(0.);
681 
682  // NOTE: These would be populated from P3 (we use the constants in p3_main_impl.hpp)
683  // NOTE: These are in units of micron!
684  eff_radius_qc_tab(icol,ilay) = (qc>Real(0.)) ? Real(10.0) : Real(0.);
685  eff_radius_qi_tab(icol,ilay) = (qi>Real(0.)) ? Real(25.0) : Real(0.);
686 
687  // Buffers on z-faces (nlay+1)
688  p_lev_tab(icol,ilay) = getPgivenRTh(rt_avg, qv_avg);
689  t_lev_tab(icol,ilay) = getTgivenRandRTh(r_avg, rt_avg, qv_avg);
690  if (ilay==(nlay-1)) {
691  Real r_hi = cons_arr(i,j,k+1,Rho_comp);
692  Real rt_hi = cons_arr(i,j,k+1,RhoTheta_comp);
693  Real qv_hi = (moist) ? std::max(cons_arr(i,j,k+1,RhoQ1_comp)/r_hi,Real(0.)) : Real(0.);
694  Real dz_kp1 = (z_arr) ? Real(0.125) * ( (z_arr(i ,j ,k+2) - z_arr(i ,j ,k+1))
695  + (z_arr(i+1,j ,k+2) - z_arr(i+1,j ,k+1))
696  + (z_arr(i ,j+1,k+2) - z_arr(i ,j+1,k+1))
697  + (z_arr(i+1,j+1,k+2) - z_arr(i+1,j+1,k+1)) ) : Real(0.5)*dz; // Dist from w-face to CC at k+1
698  r_avg = (dz_kp1*r + dz_k*r_hi ) / (dz_k + dz_kp1);
699  rt_avg = (dz_kp1*rt + dz_k*rt_hi) / (dz_k + dz_kp1);
700  qv_avg = (dz_kp1*qv + dz_k*qv_hi) / (dz_k + dz_kp1);
701  p_lev_tab(icol,ilay+1) = getPgivenRTh(rt_avg, qv_avg);
702  t_lev_tab(icol,ilay+1) = getTgivenRandRTh(r_avg, rt_avg, qv_avg);
703  }
704 
705  // 1D data structures
706  if (k==0) {
707  lat_tab(icol) = (has_lat) ? lat_arr(i,j,0) : cons_lat;
708  lon_tab(icol) = (has_lon) ? lon_arr(i,j,0) : cons_lon;
709  }
710 
711  });
712  } // mfi
713 
714  // Populate vars LSM would provide
715  if (!has_lsm && !has_surflayer) {
716  // Parsed surface temp
717  Kokkos::deep_copy(t_sfc, rad_t_sfc);
718 
719  // EAMXX dummy atmos constants
720  Kokkos::deep_copy(sfc_alb_dir_vis, Real(0.06));
721  Kokkos::deep_copy(sfc_alb_dir_nir, Real(0.06));
722  Kokkos::deep_copy(sfc_alb_dif_vis, Real(0.06));
723  Kokkos::deep_copy(sfc_alb_dif_nir, Real(0.06));
724 
725  // AML NOTE: These are not used in current EAMXX, I've left
726  // the code to plug into these if we need it.
727  //
728  // Current EAMXX constants
729  Kokkos::deep_copy(sfc_emis, Real(0.98));
730  Kokkos::deep_copy(lw_src , zero );
731  } else {
732  Vector<real1d_k> rrtmgp_in_vars = {t_sfc, sfc_emis,
735  Vector<Real> rrtmgp_default_vals = {rad_t_sfc, Real(0.98),
736  Real(0.06), Real(0.06),
737  Real(0.06), Real(0.06)};
738  for (int ivar(0); ivar<lsm_input_ptrs.size(); ivar++) {
739  auto rrtmgp_default_val = rrtmgp_default_vals[ivar];
740  auto rrtmgp_to_fill_k = rrtmgp_in_vars[ivar];
741  amrex::Table1D<amrex::Real> rrtmgp_to_fill(rrtmgp_to_fill_k.data(),
742  0, rrtmgp_to_fill_k.extent(0));
743  for (MFIter mfi(*m_cons_in); mfi.isValid(); ++mfi) {
744  const auto& vbx = mfi.validbox();
745  const auto& sbx = makeSlab(vbx,2,vbx.smallEnd(2));
746  const int nx = vbx.length(0);
747  const int imin = vbx.smallEnd(0);
748  const int jmin = vbx.smallEnd(1);
749  const int offset = m_col_offsets[mfi.index()];
750  const int k_surface = vbx.smallEnd(2);
751  const Array4<const Real>& cons_arr = m_cons_in->const_array(mfi);
752  const Array4<const Real>& z_arr = (m_z_phys) ? m_z_phys->const_array(mfi) :
753  Array4<const Real>{};
754  const Array4<const int>& lmask_arr = (lmask) ? lmask->const_array(mfi) :
755  Array4<const int> {};
756  const Array4<const Real>& tsurf_arr = (t_surf) ? t_surf->const_array(mfi) :
757  Array4<const Real> {};
758  const Array4<const Real>& lsm_in_arr = (lsm_input_ptrs[ivar]) ? lsm_input_ptrs[ivar]->const_array(mfi) :
759  Array4<const Real> {};
760  ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
761  {
762  // map [i,j,k] 0-based to [icol, ilay] 0-based
763  const int icol = (j-jmin)*nx + (i-imin) + offset;
764 
765  // Check if over land
766  bool is_land = (lmask_arr) ? lmask_arr(i,j,k) : 1;
767 
768  // Surface temperature has a distinct contract: LSM and
769  // RRTMGP use absolute temperature, while SurfaceLayer
770  // supplies potential temperature. The other LSM inputs
771  // retain their existing validity/default handling.
772  if (ivar == 0) {
773  const bool has_lsm_t_sfc = static_cast<bool>(lsm_in_arr);
774  const bool valid_lsm_t_sfc =
775  has_lsm_t_sfc && amrex::Math::isfinite(lsm_in_arr(i,j,k)) &&
776  (lsm_in_arr(i,j,k) > Real(0.)) &&
777  (lsm_in_arr(i,j,k) < lsm_undefined);
778  // Match TwoStream: convert SurfaceLayer theta with
779  // the physical surface pressure diagnosed from the
780  // lowest atmospheric cell.
782  is_land,
783  has_lsm_t_sfc,
784  valid_lsm_t_sfc,
785  has_lsm_t_sfc ? lsm_in_arr(i,j,k) : Real(0.),
786  static_cast<bool>(tsurf_arr),
787  tsurf_arr ? tsurf_arr(i,j,k) : Real(0.),
789  cons_arr(i,j,k_surface,Rho_comp),
790  cons_arr(i,j,k_surface,RhoTheta_comp),
791  moist ? std::max(cons_arr(i,j,k_surface,RhoQ1_comp) /
792  cons_arr(i,j,k_surface,Rho_comp), Real(0.)) : Real(0.),
793  z_arr ? Compute_Zrel_AtCellCenter(i,j,k_surface,z_arr) : Real(0.5)*dz),
794  rdOcp,
795  rrtmgp_default_val,
796  rrtmgp_to_fill(icol));
797  } else {
798  // Have LSM and are over land.
799  const bool valid_lsm_data =
800  (lsm_in_arr && (lsm_in_arr(i,j,k) < lsm_undefined));
801  if (is_land && valid_lsm_data) {
802  rrtmgp_to_fill(icol) = lsm_in_arr(i,j,k);
803  } else {
804  // Use the default value.
805  rrtmgp_to_fill(icol) = rrtmgp_default_val;
806  }
807  }
808  });
809  } //mfi
810  } // ivar
811  Kokkos::deep_copy(lw_src, zero );
812  } // have lsm
813 
814  // Enforce consistency between t_sfc and t_lev at bottom surface
815  Kokkos::parallel_for(Kokkos::RangePolicy(0, ncol),
816  KOKKOS_LAMBDA (int icol)
817  {
818  t_lev_tab(icol,0) = t_sfc_tab(icol);
819  });
820 }
constexpr amrex::Real lsm_undefined
Definition: ERF_Constants.H:28
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getTgivenRandRTh(const amrex::Real rho, const amrex::Real rhotheta, const amrex::Real qv=amrex::Real(0))
Definition: ERF_EOS.H:46
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getPgivenRTh(const amrex::Real rhotheta, const amrex::Real qv=amrex::Real(0))
Definition: ERF_EOS.H:81
#define Rho_comp
Definition: ERF_IndexDefines.H:42
#define RhoTheta_comp
Definition: ERF_IndexDefines.H:43
#define RhoQ2_comp
Definition: ERF_IndexDefines.H:49
#define RhoQ1_comp
Definition: ERF_IndexDefines.H:48
constexpr amrex::Real zero
Definition: ERF_NumericalConstants.H:36
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real Compute_Zrel_AtCellCenter(const int &i, const int &j, const int &k, const amrex::Array4< const amrex::Real > &z_nd)
Definition: ERF_TerrainMetrics.H:758
amrex::MultiFab * m_lon
Definition: ERF_Radiation.H:341
amrex::MultiFab * m_z_phys
Definition: ERF_Radiation.H:337
amrex::Geometry m_geom
Definition: ERF_Radiation.H:293
amrex::MultiFab * m_lat
Definition: ERF_Radiation.H:340
@ qv
Definition: ERF_Kessler.H:31
@ qc
Definition: ERF_SatAdj.H:42
@ qi
Definition: ERF_WDM6.H:28
@ dz
Definition: ERF_AdvanceWDM6.cpp:272
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real pressure_at_surface(const amrex::Real rho, const amrex::Real rho_theta, const amrex::Real qv, const amrex::Real delta_z)
Definition: ERF_SurfaceTemperature.H:36
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void resolve_surface_temperature(const bool is_land, const bool has_lsm_t_sfc, const bool valid_lsm_t_sfc, const amrex::Real lsm_t_sfc, const bool has_surface_layer, const amrex::Real surface_layer_theta, const amrex::Real surface_pressure, const amrex::Real rdOcp, const amrex::Real default_t_sfc, amrex::Real &t_sfc)
Definition: ERF_RRTMGP_SurfaceTemperature.H:39
Here is the call graph for this function:

◆ populateDatalogMF()

void Radiation::populateDatalogMF ( )
965 {
966  Table2D<Real,Order::C> sw_flux_up_tab(sw_flux_up.data(), {0,0}, {static_cast<int>(sw_flux_up.extent(0)),static_cast<int>(sw_flux_up.extent(1))});
967  Table2D<Real,Order::C> sw_flux_dn_tab(sw_flux_dn.data(), {0,0}, {static_cast<int>(sw_flux_dn.extent(0)),static_cast<int>(sw_flux_dn.extent(1))});
968  Table2D<Real,Order::C> sw_flux_dn_dir_tab(sw_flux_dn_dir.data(), {0,0}, {static_cast<int>(sw_flux_dn_dir.extent(0)),static_cast<int>(sw_flux_dn_dir.extent(1))});
969  Table2D<Real,Order::C> lw_flux_up_tab(lw_flux_up.data(), {0,0}, {static_cast<int>(lw_flux_up.extent(0)),static_cast<int>(lw_flux_up.extent(1))});
970  Table2D<Real,Order::C> lw_flux_dn_tab(lw_flux_dn.data(), {0,0}, {static_cast<int>(lw_flux_dn.extent(0)),static_cast<int>(lw_flux_dn.extent(1))});
971 
972  Table2D<Real,Order::C> sw_clrsky_flux_up_tab(sw_clrsky_flux_up.data(), {0,0},
973  {static_cast<int>(sw_clrsky_flux_up.extent(0)),static_cast<int>(sw_clrsky_flux_up.extent(1))});
974  Table2D<Real,Order::C> sw_clrsky_flux_dn_tab(sw_clrsky_flux_dn.data(), {0,0},
975  {static_cast<int>(sw_clrsky_flux_dn.extent(0)),static_cast<int>(sw_clrsky_flux_dn.extent(1))});
976  Table2D<Real,Order::C> sw_clrsky_flux_dn_dir_tab(sw_clrsky_flux_dn_dir.data(), {0,0},
977  {static_cast<int>(sw_clrsky_flux_dn_dir.extent(0)),static_cast<int>(sw_clrsky_flux_dn_dir.extent(1))});
978  Table2D<Real,Order::C> lw_clrsky_flux_up_tab(lw_clrsky_flux_up.data(), {0,0},
979  {static_cast<int>(lw_clrsky_flux_up.extent(0)),static_cast<int>(lw_clrsky_flux_up.extent(1))});
980  Table2D<Real,Order::C> lw_clrsky_flux_dn_tab(lw_clrsky_flux_dn.data(), {0,0},
981  {static_cast<int>(lw_clrsky_flux_dn.extent(0)),static_cast<int>(lw_clrsky_flux_dn.extent(1))});
982  Table2D<Real,Order::C> sw_clrsky_heating_tab(sw_clrsky_heating.data(), {0,0},
983  {static_cast<int>(sw_clrsky_heating.extent(0)),static_cast<int>(sw_clrsky_heating.extent(1))});
984  Table2D<Real,Order::C> lw_clrsky_heating_tab(lw_clrsky_heating.data(), {0,0},
985  {static_cast<int>(lw_clrsky_heating.extent(0)),static_cast<int>(lw_clrsky_heating.extent(1))});
986  Table2D<Real,Order::C> sw_clnsky_flux_up_tab(sw_clnsky_flux_up.data(), {0,0},
987  {static_cast<int>(sw_clnsky_flux_up.extent(0)),static_cast<int>(sw_clnsky_flux_up.extent(1))});
988  Table2D<Real,Order::C> sw_clnsky_flux_dn_tab(sw_clnsky_flux_dn.data(), {0,0},
989  {static_cast<int>(sw_clnsky_flux_dn.extent(0)),static_cast<int>(sw_clnsky_flux_dn.extent(1))});
990  Table2D<Real,Order::C> sw_clnsky_flux_dn_dir_tab(sw_clnsky_flux_dn_dir.data(), {0,0},
991  {static_cast<int>(sw_clnsky_flux_dn_dir.extent(0)),static_cast<int>(sw_clnsky_flux_dn_dir.extent(1))});
992  Table2D<Real,Order::C> lw_clnsky_flux_up_tab(lw_clnsky_flux_up.data(), {0,0},
993  {static_cast<int>(lw_clnsky_flux_up.extent(0)),static_cast<int>(lw_clnsky_flux_up.extent(1))});
994  Table2D<Real,Order::C> lw_clnsky_flux_dn_tab(lw_clnsky_flux_dn.data(), {0,0},
995  {static_cast<int>(lw_clnsky_flux_dn.extent(0)),static_cast<int>(lw_clnsky_flux_dn.extent(1))});
996  Table2D<Real,Order::C> sw_clnclrsky_flux_up_tab(sw_clnclrsky_flux_up.data(), {0,0},
997  {static_cast<int>(sw_clnclrsky_flux_up.extent(0)),static_cast<int>(sw_clnclrsky_flux_up.extent(1))});
998  Table2D<Real,Order::C> sw_clnclrsky_flux_dn_tab(sw_clnclrsky_flux_dn.data(), {0,0},
999  {static_cast<int>(sw_clnclrsky_flux_dn.extent(0)),static_cast<int>(sw_clnclrsky_flux_dn.extent(1))});
1000  Table2D<Real,Order::C> sw_clnclrsky_flux_dn_dir_tab(sw_clnclrsky_flux_dn_dir.data(), {0,0},
1001  {static_cast<int>(sw_clnclrsky_flux_dn_dir.extent(0)),static_cast<int>(sw_clnclrsky_flux_dn_dir.extent(1))});
1002  Table2D<Real,Order::C> lw_clnclrsky_flux_up_tab(lw_clnclrsky_flux_up.data(), {0,0},
1003  {static_cast<int>(lw_clnclrsky_flux_up.extent(0)),static_cast<int>(lw_clnclrsky_flux_up.extent(1))});
1004  Table2D<Real,Order::C> lw_clnclrsky_flux_dn_tab(lw_clnclrsky_flux_dn.data(), {0,0},
1005  {static_cast<int>(lw_clnclrsky_flux_dn.extent(0)),static_cast<int>(lw_clnclrsky_flux_dn.extent(1))});
1006 
1007  Table1D<Real> mu0_tab(mu0.data(), {0}, {static_cast<int>(mu0.extent(0))});
1008 
1009  auto extra_clnsky_diag = m_extra_clnsky_diag;
1010  auto extra_clnclrsky_diag = m_extra_clnclrsky_diag;
1011 
1012  for (MFIter mfi(datalog_mf); mfi.isValid(); ++mfi) {
1013  const auto& vbx = mfi.validbox();
1014  const int nx = vbx.length(0);
1015  const int imin = vbx.smallEnd(0);
1016  const int jmin = vbx.smallEnd(1);
1017  const int offset = m_col_offsets[mfi.index()];
1018  const Array4<Real>& dst_arr = datalog_mf.array(mfi);
1019  const Array4<Real>& q_arr = m_qheating_rates->array(mfi);
1020  ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
1021  {
1022  // map [i,j,k] 0-based to [icol, ilay] 0-based
1023  const int icol = (j-jmin)*nx + (i-imin) + offset;
1024  const int ilay = k;
1025 
1026  dst_arr(i,j,k,0) = q_arr(i, j, k, 0);
1027  dst_arr(i,j,k,1) = q_arr(i, j, k, 1);
1028 
1029  // SW and LW fluxes
1030  dst_arr(i,j,k,2) = sw_flux_up_tab(icol,ilay);
1031  dst_arr(i,j,k,3) = sw_flux_dn_tab(icol,ilay);
1032  dst_arr(i,j,k,4) = sw_flux_dn_dir_tab(icol,ilay);
1033  dst_arr(i,j,k,5) = lw_flux_up_tab(icol,ilay);
1034  dst_arr(i,j,k,6) = lw_flux_dn_tab(icol,ilay);
1035 
1036  // Cosine zenith angle
1037  dst_arr(i,j,k,7) = mu0_tab(icol);
1038 
1039  // Clear sky heating rates and fluxes:
1040  dst_arr(i,j,k,8) = sw_clrsky_heating_tab(icol, ilay);
1041  dst_arr(i,j,k,9) = lw_clrsky_heating_tab(icol, ilay);
1042 
1043  dst_arr(i,j,k,10) = sw_clrsky_flux_up_tab(icol,ilay);
1044  dst_arr(i,j,k,11) = sw_clrsky_flux_dn_tab(icol,ilay);
1045  dst_arr(i,j,k,12) = sw_clrsky_flux_dn_dir_tab(icol,ilay);
1046  dst_arr(i,j,k,13) = lw_clrsky_flux_up_tab(icol,ilay);
1047  dst_arr(i,j,k,14) = lw_clrsky_flux_dn_tab(icol,ilay);
1048 
1049  // Clean sky fluxes:
1050  if (extra_clnsky_diag) {
1051  dst_arr(i,j,k,15) = sw_clnsky_flux_up_tab(icol,ilay);
1052  dst_arr(i,j,k,16) = sw_clnsky_flux_dn_tab(icol,ilay);
1053  dst_arr(i,j,k,17) = sw_clnsky_flux_dn_dir_tab(icol,ilay);
1054  dst_arr(i,j,k,18) = lw_clnsky_flux_up_tab(icol,ilay);
1055  dst_arr(i,j,k,19) = lw_clnsky_flux_dn_tab(icol,ilay);
1056  }
1057 
1058  // Clean-clear sky fluxes:
1059  if (extra_clnclrsky_diag) {
1060  dst_arr(i,j,k,20) = sw_clnclrsky_flux_up_tab(icol,ilay);
1061  dst_arr(i,j,k,21) = sw_clnclrsky_flux_dn_tab(icol,ilay);
1062  dst_arr(i,j,k,22) = sw_clnclrsky_flux_dn_dir_tab(icol,ilay);
1063  dst_arr(i,j,k,23) = lw_clnclrsky_flux_up_tab(icol,ilay);
1064  dst_arr(i,j,k,24) = lw_clnclrsky_flux_dn_tab(icol,ilay);
1065  }
1066  });
1067  }
1068 }
amrex::MultiFab datalog_mf
Definition: ERF_Radiation.H:348
Here is the call graph for this function:

◆ rad_run_impl()

void Radiation::rad_run_impl ( const amrex::Vector< amrex::MultiFab * > &  lsm_output_ptrs)
inline
230  {
231  if (m_update_rad) {
232  amrex::Print() << "Radiation advancing level " << m_lev << " at (YY-MM-DD SS) " << m_orbital_year << '-'
233  << m_orbital_mon << '-' << m_orbital_day << ' ' << m_orbital_sec << " ...";
234  this->initialize_impl();
235  this->run_impl();
236  this->finalize_impl(lsm_output_ptrs);
237  amrex::Print() << "DONE\n";
238  }
239  }
void initialize_impl()
Definition: ERF_Radiation.cpp:1205
void run_impl()
Definition: ERF_Radiation.cpp:1230
int m_orbital_mon
Definition: ERF_Radiation.H:409
int m_orbital_sec
Definition: ERF_Radiation.H:411
int m_lev
Definition: ERF_Radiation.H:281
bool m_update_rad
Definition: ERF_Radiation.H:299
int m_orbital_day
Definition: ERF_Radiation.H:410
void finalize_impl(const amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
Definition: ERF_Radiation.cpp:1527

Referenced by Run().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ radiation_updated()

bool Radiation::radiation_updated ( ) const
inlineoverridevirtual

Reports whether the last radiation call updated its outputs.

Returns
True when the radiation kernel wrote output fields.

Reimplemented from IRadiation.

172 { return m_update_rad && !m_is_nested_patch; }

◆ Run()

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

Implements IRadiation.

149  {
150  // Skip radiation computation on nested patches - heating rates will be
151  // filled via interpolation from coarser level in ERF::FillPatch
152  if (m_is_nested_patch) {
153  return;
154  }
155 
157  m_nlay == geom.Domain().length(2),
158  "Radiation: m_nlay (" + std::to_string(m_nlay) +
159  ") doesn't match domain vertical extent (" +
160  std::to_string(geom.Domain().length(2)) + ")");
161 
162  set_grids(level, step, time, dt, ba, geom,
163  cons_in, lmask, t_surf,
164  lsm_input_ptrs, qheating_rates,
165  rad_fluxes, z_phys, lat_ptr, lon_ptr,
166  updated_lsm);
167  rad_run_impl(lsm_output_ptrs);
168  declin = m_declin;
169  calday = m_calday;
170  }
void set_grids(int &level, int &step, double &time, const double &dt, const amrex::BoxArray &ba, amrex::Geometry &geom, amrex::MultiFab *cons_in, amrex::iMultiFab *lmask, amrex::MultiFab *t_surf, const amrex::Vector< const amrex::MultiFab * > &lsm_input_ptrs, amrex::MultiFab *qheating_rates, amrex::MultiFab *rad_fluxes, amrex::MultiFab *z_phys, amrex::MultiFab *lat, amrex::MultiFab *lon, const bool updated_lsm)
Definition: ERF_Radiation.cpp:262
void rad_run_impl(const amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
Definition: ERF_Radiation.H:229
amrex::Real m_declin
Definition: ERF_Radiation.H:431
amrex::Real m_calday
Definition: ERF_Radiation.H:432
Here is the call graph for this function:

◆ run_impl()

void Radiation::run_impl ( )
1231 {
1232  // A rank that owns no boxes on this level has no columns and therefore no
1233  // radiation work to do. Bail out before the chunk loop; there are no MPI
1234  // collectives in this routine, so returning early cannot deadlock.
1235  if (m_ncol == 0) { return; }
1236 
1237  // Local copies
1238  const auto ncol = m_ncol;
1239  const auto nlay = m_nlay;
1240  const auto nswbands = m_nswbands;
1241 
1242  // Compute orbital parameters; these are used both for computing
1243  // the solar zenith angle and also for computing total solar
1244  // irradiance scaling (tsi_scaling).
1245  double obliqr, lambm0, mvelpp;
1246  // Any of eccen/obliq/mvelp that the user set (i.e. is non-negative) is used
1247  // as is by orbital_params; the rest are computed from the orbital year.
1248  int orbital_year = m_orbital_year;
1249  double eccen = m_orbital_eccen;
1250  double obliq = m_orbital_obliq;
1251  double mvelp = m_orbital_mvelp;
1252  orbital_params(orbital_year, eccen, obliq,
1253  mvelp, obliqr, lambm0, mvelpp);
1254 
1255  // Use the orbital parameters to calculate the solar declination and eccentricity factor
1256  double delta, eccf;
1257  // Day of the year plus fraction, calday 1 == Jan 1 0Z (leap-aware)
1259  orbital_decl(calday, eccen, mvelpp, lambm0, obliqr, delta, eccf);
1260  m_calday = calday;
1261  m_declin = delta;
1262 
1263  // Overwrite eccf if using a fixed solar constant.
1264  auto fixed_total_solar_irradiance = m_fixed_total_solar_irradiance;
1265  if (fixed_total_solar_irradiance >= 0){
1266  eccf = fixed_total_solar_irradiance/Real(1360.9);
1267  }
1268 
1269  // Precompute volume mixing ratio (VMR) for all gases
1270  //
1271  // H2O is obtained from qv.
1272  // O3 may be a constant or a 1D vector
1273  // All other comps are set to constants for now
1274  Vector<real2d_k> vmr_full_vec(m_ngas);
1275  for (int igas(0); igas < m_ngas; ++igas) {
1276  auto name = m_gas_names[igas];
1277  vmr_full_vec[igas] = real2d_k("vmr_full_" + name, ncol, nlay);
1278  auto tmp2d = vmr_full_vec[igas];
1279  auto gas_mol_weight = m_mol_weight_gas[igas];
1280  if (name == "H2O") {
1281  auto qv_lay_d = qv_lay;
1282  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
1283  KOKKOS_LAMBDA (int icol, int ilay)
1284  {
1285  tmp2d(icol,ilay) = qv_lay_d(icol,ilay) * mwdair/gas_mol_weight;
1286  });
1287  } else if (name == "CO2") {
1288  Kokkos::deep_copy(tmp2d, m_co2vmr);
1289  } else if (name == "O3") {
1290  if (m_o3_size==1) {
1291  Kokkos::deep_copy(tmp2d, m_o3vmr[0] );
1292  } else {
1293  auto o3_lay_d = o3_lay;
1294  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
1295  KOKKOS_LAMBDA (int icol, int ilay)
1296  {
1297  tmp2d(icol,ilay) = o3_lay_d(ilay);
1298  });
1299  }
1300  } else if (name == "N2O") {
1301  Kokkos::deep_copy(tmp2d, m_n2ovmr);
1302  } else if (name == "CO") {
1303  Kokkos::deep_copy(tmp2d, m_covmr );
1304  } else if (name == "CH4") {
1305  Kokkos::deep_copy(tmp2d, m_ch4vmr);
1306  } else if (name == "O2") {
1307  Kokkos::deep_copy(tmp2d, m_o2vmr );
1308  } else if (name == "N2") {
1309  Kokkos::deep_copy(tmp2d, m_n2vmr );
1310  } else {
1311  Abort("Radiation: Unknown gas component.");
1312  }
1313  }
1314 
1315  // Determine the cosine zenith angle.
1316  // Populate mu0 1D array
1317  // This must be done on HOST and copied to device.
1318  auto h_mu0 = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), mu0);
1319  if (m_fixed_solar_zenith_angle > 0) {
1320  Kokkos::deep_copy(h_mu0, m_fixed_solar_zenith_angle);
1321  } else {
1322  auto h_lat = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), lat);
1323  auto h_lon = Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace(), lon);
1324  double dt = double(m_dt);
1325  auto rad_freq_in_steps = m_rad_freq_in_steps;
1326  Kokkos::parallel_for(Kokkos::RangePolicy<Kokkos::Serial>(0, ncol),
1327  [&,PI_d=PI] (int icol)
1328  {
1329  // Convert lat/lon to radians
1330  double lat_col = h_lat(icol)*PI_d/Real(180.0);
1331  double lon_col = h_lon(icol)*PI_d/Real(180.0);
1332  double lcalday = calday;
1333  double ldelta = delta;
1334  double dt_avg = static_cast<double>(rad_freq_in_steps) * dt;
1335  h_mu0(icol) = Real(orbital_cos_zenith(lcalday, lat_col, lon_col, ldelta, dt_avg));
1336  });
1337  }
1338  Kokkos::deep_copy(mu0, h_mu0);
1339 
1340  // Compute layer cloud mass per unit area (populates lwp/iwp)
1343 
1344  // Convert to g/m2 (needed by RRTMGP)
1345  Table2D<Real,Order::C> lwp_tab(lwp.data(), {0,0}, {static_cast<int>(lwp.extent(0)),static_cast<int>(lwp.extent(1))});
1346  Table2D<Real,Order::C> iwp_tab(iwp.data(), {0,0}, {static_cast<int>(iwp.extent(0)),static_cast<int>(iwp.extent(1))});
1347  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol, nlay}),
1348  KOKKOS_LAMBDA (int icol, int ilay)
1349  {
1350  lwp_tab(icol,ilay) *= Real(1.e3);
1351  iwp_tab(icol,ilay) *= Real(1.e3);
1352  });
1353 
1354  // -----------------------------------------------------------------------
1355  // Process radiation in column chunks to limit peak GPU memory.
1356  // Radiation columns are independent (no horizontal coupling), so
1357  // chunking produces bit-identical results.
1358  // -----------------------------------------------------------------------
1359  const int ncol_chunk = std::min(m_ncol_chunk, ncol);
1360  const int kbot = 0;
1361 
1362  for (int col_s = 0; col_s < ncol; col_s += ncol_chunk) {
1363  const int ncol_c = std::min(ncol_chunk, ncol - col_s);
1364 
1365  // --- Chunk subviews: 1D (ncol) ---
1366  real1d_k mu0_c (mu0.data() + col_s, ncol_c);
1367  real1d_k sfc_alb_dir_vis_c (sfc_alb_dir_vis.data() + col_s, ncol_c);
1368  real1d_k sfc_alb_dir_nir_c (sfc_alb_dir_nir.data() + col_s, ncol_c);
1369  real1d_k sfc_alb_dif_vis_c (sfc_alb_dif_vis.data() + col_s, ncol_c);
1370  real1d_k sfc_alb_dif_nir_c (sfc_alb_dif_nir.data() + col_s, ncol_c);
1371  real1d_k sfc_flux_dir_vis_c (sfc_flux_dir_vis.data() + col_s, ncol_c);
1372  real1d_k sfc_flux_dir_nir_c (sfc_flux_dir_nir.data() + col_s, ncol_c);
1373  real1d_k sfc_flux_dif_vis_c (sfc_flux_dif_vis.data() + col_s, ncol_c);
1374  real1d_k sfc_flux_dif_nir_c (sfc_flux_dif_nir.data() + col_s, ncol_c);
1375  real1d_k t_sfc_c (t_sfc.data() + col_s, ncol_c);
1376  real1d_k sfc_emis_c (sfc_emis.data() + col_s, ncol_c);
1377  real1d_k lw_src_c (lw_src.data() + col_s, ncol_c);
1378 
1379  // --- Chunk subviews: 2D (ncol, nlay) via LayoutRight pointer offset ---
1380  const int stride2_nlay = nlay;
1381  const int stride2_nlayp1 = nlay + 1;
1382  real2d_k p_lay_c (p_lay.data() + col_s*stride2_nlay, ncol_c, nlay);
1383  real2d_k t_lay_c (t_lay.data() + col_s*stride2_nlay, ncol_c, nlay);
1384  real2d_k r_lay_c (r_lay.data() + col_s*stride2_nlay, ncol_c, nlay);
1385  real2d_k z_del_c (z_del.data() + col_s*stride2_nlay, ncol_c, nlay);
1386  real2d_k lwp_c (lwp.data() + col_s*stride2_nlay, ncol_c, nlay);
1387  real2d_k iwp_c (iwp.data() + col_s*stride2_nlay, ncol_c, nlay);
1388  real2d_k eff_radius_qc_c(eff_radius_qc.data() + col_s*stride2_nlay, ncol_c, nlay);
1389  real2d_k eff_radius_qi_c(eff_radius_qi.data() + col_s*stride2_nlay, ncol_c, nlay);
1390  real2d_k cldfrac_tot_c (cldfrac_tot.data() + col_s*stride2_nlay, ncol_c, nlay);
1391  real2d_k sw_heating_c (sw_heating.data() + col_s*stride2_nlay, ncol_c, nlay);
1392  real2d_k lw_heating_c (lw_heating.data() + col_s*stride2_nlay, ncol_c, nlay);
1393 
1394  // --- Chunk subviews: 2D (ncol, nlay+1) ---
1395  real2d_k p_lev_c (p_lev.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1396  real2d_k t_lev_c (t_lev.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1397  real2d_k sw_flux_up_c (sw_flux_up.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1398  real2d_k sw_flux_dn_c (sw_flux_dn.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1399  real2d_k sw_flux_dn_dir_c (sw_flux_dn_dir.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1400  real2d_k lw_flux_up_c (lw_flux_up.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1401  real2d_k lw_flux_dn_c (lw_flux_dn.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1402  // Clear-sky flux subviews (always active)
1403  // NOTE: once on m_ncol_chunk if not writing a datalog
1404  real2d_k sw_clrsky_flux_up_c, sw_clrsky_flux_dn_c, sw_clrsky_flux_dn_dir_c;
1405  real2d_k lw_clrsky_flux_up_c, lw_clrsky_flux_dn_c;
1406  if (datalog_int > 0) {
1407  sw_clrsky_flux_up_c = real2d_k(sw_clrsky_flux_up.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1408  sw_clrsky_flux_dn_c = real2d_k(sw_clrsky_flux_dn.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1409  sw_clrsky_flux_dn_dir_c = real2d_k(sw_clrsky_flux_dn_dir.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1410  lw_clrsky_flux_up_c = real2d_k(lw_clrsky_flux_up.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1411  lw_clrsky_flux_dn_c = real2d_k(lw_clrsky_flux_dn.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1412  } else {
1413  sw_clrsky_flux_up_c = real2d_k(sw_clrsky_flux_up.data() , ncol_c, nlay+1);
1414  sw_clrsky_flux_dn_c = real2d_k(sw_clrsky_flux_dn.data() , ncol_c, nlay+1);
1415  sw_clrsky_flux_dn_dir_c = real2d_k(sw_clrsky_flux_dn_dir.data() , ncol_c, nlay+1);
1416  lw_clrsky_flux_up_c = real2d_k(lw_clrsky_flux_up.data() , ncol_c, nlay+1);
1417  lw_clrsky_flux_dn_c = real2d_k(lw_clrsky_flux_dn.data() , ncol_c, nlay+1);
1418  }
1419 
1420  // Diagnostic flux subviews (placeholder when disabled)
1421  real2d_k sw_clnclrsky_flux_up_c, sw_clnclrsky_flux_dn_c, sw_clnclrsky_flux_dn_dir_c;
1422  real2d_k lw_clnclrsky_flux_up_c, lw_clnclrsky_flux_dn_c;
1423  if (m_extra_clnclrsky_diag) {
1424  sw_clnclrsky_flux_up_c = real2d_k(sw_clnclrsky_flux_up.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1425  sw_clnclrsky_flux_dn_c = real2d_k(sw_clnclrsky_flux_dn.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1426  sw_clnclrsky_flux_dn_dir_c = real2d_k(sw_clnclrsky_flux_dn_dir.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1427  lw_clnclrsky_flux_up_c = real2d_k(lw_clnclrsky_flux_up.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1428  lw_clnclrsky_flux_dn_c = real2d_k(lw_clnclrsky_flux_dn.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1429  } else {
1430  sw_clnclrsky_flux_up_c = real2d_k("sw_clnclrsky_flux_up_c" , 1, 1);
1431  sw_clnclrsky_flux_dn_c = real2d_k("sw_clnclrsky_flux_dn_c" , 1, 1);
1432  sw_clnclrsky_flux_dn_dir_c = real2d_k("sw_clnclrsky_flux_dn_dir_c", 1, 1);
1433  lw_clnclrsky_flux_up_c = real2d_k("lw_clnclrsky_flux_up_c" , 1, 1);
1434  lw_clnclrsky_flux_dn_c = real2d_k("lw_clnclrsky_flux_dn_c" , 1, 1);
1435  }
1436 
1437  real2d_k sw_clnsky_flux_up_c, sw_clnsky_flux_dn_c, sw_clnsky_flux_dn_dir_c;
1438  real2d_k lw_clnsky_flux_up_c, lw_clnsky_flux_dn_c;
1439  if (m_extra_clnsky_diag) {
1440  sw_clnsky_flux_up_c = real2d_k(sw_clnsky_flux_up.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1441  sw_clnsky_flux_dn_c = real2d_k(sw_clnsky_flux_dn.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1442  sw_clnsky_flux_dn_dir_c = real2d_k(sw_clnsky_flux_dn_dir.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1443  lw_clnsky_flux_up_c = real2d_k(lw_clnsky_flux_up.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1444  lw_clnsky_flux_dn_c = real2d_k(lw_clnsky_flux_dn.data() + col_s*stride2_nlayp1, ncol_c, nlay+1);
1445  } else {
1446  sw_clnsky_flux_up_c = real2d_k("sw_clnsky_flux_up_c" , 1, 1);
1447  sw_clnsky_flux_dn_c = real2d_k("sw_clnsky_flux_dn_c" , 1, 1);
1448  sw_clnsky_flux_dn_dir_c = real2d_k("sw_clnsky_flux_dn_dir_c", 1, 1);
1449  lw_clnsky_flux_up_c = real2d_k("lw_clnsky_flux_up_c" , 1, 1);
1450  lw_clnsky_flux_dn_c = real2d_k("lw_clnsky_flux_dn_c" , 1, 1);
1451  }
1452 
1453  // --- Chunk subviews: 2D (ncol, nswbands) ---
1454  real2d_k sfc_alb_dir_c(sfc_alb_dir.data() + col_s*nswbands, ncol_c, nswbands);
1455  real2d_k sfc_alb_dif_c(sfc_alb_dif.data() + col_s*nswbands, ncol_c, nswbands);
1456 
1457  // --- Chunk subviews: 3D (ncol, nlay+1, nbands) ---
1458  // NOTE: Allocate these once on m_ncol_chunk and use what we need in the chunk loop
1459  real3d_k sw_bnd_flux_up_c (sw_bnd_flux_up.data() , ncol_c, nlay+1, nswbands);
1460  real3d_k sw_bnd_flux_dn_c (sw_bnd_flux_dn.data() , ncol_c, nlay+1, nswbands);
1461  real3d_k sw_bnd_flux_dir_c(sw_bnd_flux_dir.data(), ncol_c, nlay+1, nswbands);
1462  real3d_k sw_bnd_flux_dif_c(sw_bnd_flux_dif.data(), ncol_c, nlay+1, nswbands);
1463  real3d_k lw_bnd_flux_up_c (lw_bnd_flux_up.data() , ncol_c, nlay+1, m_nlwbands);
1464  real3d_k lw_bnd_flux_dn_c (lw_bnd_flux_dn.data() , ncol_c, nlay+1, m_nlwbands);
1465 
1466  // --- Create chunk gas concentrations by subsetting from pre-fetched VMR ---
1467  gas_concs_t gas_concs_c;
1468  gas_concs_c.init(gas_names_offset, ncol_c, nlay);
1469  for (int igas = 0; igas < m_ngas; ++igas) {
1470  real2d_k vmr_c("vmr_c", ncol_c, nlay);
1471  auto vmr_full = vmr_full_vec[igas];
1472  auto cs = col_s;
1473  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<2>>({0, 0}, {ncol_c, nlay}),
1474  KOKKOS_LAMBDA (int i, int j) {
1475  vmr_c(i, j) = vmr_full(cs + i, j);
1476  });
1477  gas_concs_c.set_vmr(m_gas_names[igas], vmr_c);
1478  }
1479 
1480  // Expand surface albedos along nswbands for this chunk
1482  sfc_alb_dir_vis_c, sfc_alb_dir_nir_c,
1483  sfc_alb_dif_vis_c, sfc_alb_dif_nir_c,
1484  sfc_alb_dir_c , sfc_alb_dif_c);
1485 
1486  // Run RRTMGP driver for this column chunk
1487  rrtmgp::rrtmgp_main(ncol_c, m_nlay,
1488  p_lay_c, t_lay_c,
1489  p_lev_c, t_lev_c,
1490  gas_concs_c,
1491  sfc_alb_dir_c, sfc_alb_dif_c, mu0_c,
1492  t_sfc_c, sfc_emis_c, lw_src_c,
1493  lwp_c, iwp_c, eff_radius_qc_c, eff_radius_qi_c, cldfrac_tot_c,
1494  sw_flux_up_c, sw_flux_dn_c, sw_flux_dn_dir_c,
1495  lw_flux_up_c, lw_flux_dn_c,
1496  sw_clnclrsky_flux_up_c, sw_clnclrsky_flux_dn_c, sw_clnclrsky_flux_dn_dir_c,
1497  sw_clrsky_flux_up_c, sw_clrsky_flux_dn_c, sw_clrsky_flux_dn_dir_c,
1498  sw_clnsky_flux_up_c, sw_clnsky_flux_dn_c, sw_clnsky_flux_dn_dir_c,
1499  lw_clnclrsky_flux_up_c, lw_clnclrsky_flux_dn_c,
1500  lw_clrsky_flux_up_c, lw_clrsky_flux_dn_c,
1501  lw_clnsky_flux_up_c, lw_clnsky_flux_dn_c,
1502  sw_bnd_flux_up_c, sw_bnd_flux_dn_c, sw_bnd_flux_dir_c,
1503  lw_bnd_flux_up_c, lw_bnd_flux_dn_c,
1505 
1506  // Compute heating rates for this chunk
1507  rrtmgp::compute_heating_rate(sw_flux_up_c, sw_flux_dn_c, r_lay_c, z_del_c, sw_heating_c);
1508  rrtmgp::compute_heating_rate(lw_flux_up_c, lw_flux_dn_c, r_lay_c, z_del_c, lw_heating_c);
1509 
1510  // Compute diffuse band fluxes and broadband surface fluxes for this chunk
1511  Kokkos::parallel_for(Kokkos::MDRangePolicy<Kokkos::Rank<3>>({0, 0, 0}, {ncol_c, nlay+1, nswbands}),
1512  KOKKOS_LAMBDA (int icol, int ilay, int ibnd)
1513  {
1514  sw_bnd_flux_dif_c(icol,ilay,ibnd) = sw_bnd_flux_dn_c(icol,ilay,ibnd) - sw_bnd_flux_dir_c(icol,ilay,ibnd);
1515  });
1516  rrtmgp::compute_broadband_surface_fluxes(ncol_c, kbot, nswbands,
1517  sw_bnd_flux_dir_c , sw_bnd_flux_dif_c ,
1518  sfc_flux_dir_vis_c, sfc_flux_dir_nir_c,
1519  sfc_flux_dif_vis_c, sfc_flux_dif_nir_c);
1520 
1521  gas_concs_c.reset();
1522  } // end column chunk loop
1523 }
constexpr amrex::Real mwdair
Definition: ERF_Constants.H:78
constexpr amrex::Real PI
Definition: ERF_NumericalConstants.H:46
AMREX_GPU_HOST AMREX_FORCE_INLINE double orbital_calday(int year, int mon, int day, int sec)
Definition: ERF_OrbCosZenith.H:485
AMREX_GPU_HOST AMREX_FORCE_INLINE double orbital_cos_zenith(double &jday, double &lat, double &lon, double &declin, double dt_avg=-one, double uniform_angle=-one, double constant_zenith_angle_deg=-one)
Definition: ERF_OrbCosZenith.H:620
AMREX_GPU_HOST AMREX_FORCE_INLINE void orbital_decl(double &calday, double &eccen, double &mvelpp, double &lambm0, double &obliqr, double &delta, double &eccf)
Definition: ERF_OrbCosZenith.H:18
AMREX_GPU_HOST AMREX_FORCE_INLINE void orbital_params(int &iyear_AD, double &eccen, double &obliq, double &mvelp, double &obliqr, double &lambm0, double &mvelpp)
Definition: ERF_OrbCosZenith.H:84
std::string name
Definition: ERF_Plotfile2DCatalog.cpp:105
double m_dt
Definition: ERF_Radiation.H:290
real(c_double), private cs
Definition: ERF_module_mp_morr_two_moment.F90:203
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, real1d_k &sfc_emis, real1d_k &lw_src, real2d_k &lwp, real2d_k &iwp, real2d_k &rel, real2d_k &rei, real2d_k &cldfrac, 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:393
void compute_band_by_band_surface_albedos(const int ncol, const int nswbands, real1d_k &sfc_alb_dir_vis, real1d_k &sfc_alb_dir_nir, real1d_k &sfc_alb_dif_vis, real1d_k &sfc_alb_dif_nir, real2d_k &sfc_alb_dir, real2d_k &sfc_alb_dif)
Definition: ERF_RRTMGP_Interface.cpp:291
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:335
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

Referenced by rad_run_impl().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ set_grids()

void Radiation::set_grids ( int &  level,
int &  step,
double &  time,
const double &  dt,
const amrex::BoxArray &  ba,
amrex::Geometry &  geom,
amrex::MultiFab *  cons_in,
amrex::iMultiFab *  lmask,
amrex::MultiFab *  t_surf,
const amrex::Vector< const amrex::MultiFab * > &  lsm_input_ptrs,
amrex::MultiFab *  qheating_rates,
amrex::MultiFab *  rad_fluxes,
amrex::MultiFab *  z_phys,
amrex::MultiFab *  lat,
amrex::MultiFab *  lon,
const bool  updated_lsm 
)
279 {
280  // Set data members that may change
281  m_lev = level;
282  m_step = step;
283  m_time = time;
284  m_dt = dt;
285  m_geom = geom;
286  m_cons_in = cons_in;
287  m_qheating_rates = qheating_rates;
288  m_rad_fluxes = rad_fluxes;
289  m_z_phys = z_phys;
290  m_lat = lat;
291  m_lon = lon;
292 
293  // Update the day and month
294  time_t timestamp = time_t(time);
295  struct tm *timeinfo = gmtime(&timestamp);
296  if (m_fixed_orbital_year) {
297  m_orbital_mon = timeinfo->tm_mon + 1;
298  m_orbital_day = timeinfo->tm_mday;
299  m_orbital_sec = timeinfo->tm_hour*3600 + timeinfo->tm_min*60 + timeinfo->tm_sec;
300  } else {
301  m_orbital_year = timeinfo->tm_year + 1900;
302  m_orbital_mon = timeinfo->tm_mon + 1;
303  m_orbital_day = timeinfo->tm_mday;
304  m_orbital_sec = timeinfo->tm_hour*3600 + timeinfo->tm_min*60 + timeinfo->tm_sec;
305  }
306 
307  // Only allocate and proceed if we are going to update radiation
308  m_update_rad = false;
309  if (m_rad_freq_in_steps > 0) { m_update_rad = ( (m_step == 0) || (m_step % m_rad_freq_in_steps == 0) ); }
310 
311  if (m_update_rad) {
312  // Call to Init() has set the dimensions: ncol & nlay
313 
314  // Allocate the buffer arrays
315  alloc_buffers();
316 
317  // Fill the KOKKOS Views from AMReX MFs
318  mf_to_kokkos_buffers(lmask, t_surf, lsm_input_ptrs);
319 
320  // (Re)define the datalog MF whenever the grids change; this must always
321  // match the layout of cons_in since populateDatalogMF() iterates over it
322  // while indexing m_col_offsets and m_qheating_rates.
323  if (datalog_int > 0) {
324  bool needs_define = ( (datalog_mf.boxArray() != cons_in->boxArray()) ||
325  (datalog_mf.DistributionMap() != cons_in->DistributionMap()) );
326  if (needs_define) {
327  datalog_mf.define(cons_in->boxArray(), cons_in->DistributionMap(), 25, 0);
328  datalog_mf.setVal(0.0);
329  }
330  }
331  }
332 }
void mf_to_kokkos_buffers(amrex::iMultiFab *lmask, amrex::MultiFab *t_surf, const amrex::Vector< const amrex::MultiFab * > &lsm_input_ptrs)
Definition: ERF_Radiation.cpp:580
int m_step
Definition: ERF_Radiation.H:284
void alloc_buffers()
Definition: ERF_Radiation.cpp:335
double m_time
Definition: ERF_Radiation.H:287

Referenced by Run().

Here is the caller graph for this function:

◆ write_rrtmgp_fluxes()

void Radiation::write_rrtmgp_fluxes ( )
925 {
926  Table2D<Real,Order::C> sw_flux_up_tab(sw_flux_up.data(), {0,0}, {static_cast<int>(sw_flux_up.extent(0)),static_cast<int>(sw_flux_up.extent(1))});
927  Table2D<Real,Order::C> sw_flux_dn_tab(sw_flux_dn.data(), {0,0}, {static_cast<int>(sw_flux_dn.extent(0)),static_cast<int>(sw_flux_dn.extent(1))});
928  Table2D<Real,Order::C> sw_flux_dn_dir_tab(sw_flux_dn_dir.data(), {0,0}, {static_cast<int>(sw_flux_dn_dir.extent(0)),static_cast<int>(sw_flux_dn_dir.extent(1))});
929  Table2D<Real,Order::C> lw_flux_up_tab(lw_flux_up.data(), {0,0}, {static_cast<int>(lw_flux_up.extent(0)),static_cast<int>(lw_flux_up.extent(1))});
930  Table2D<Real,Order::C> lw_flux_dn_tab(lw_flux_dn.data(), {0,0}, {static_cast<int>(lw_flux_dn.extent(0)),static_cast<int>(lw_flux_dn.extent(1))});
931 
932  int n_fluxes = 5;
933  MultiFab mf_flux(m_cons_in->boxArray(), m_cons_in->DistributionMap(), n_fluxes, 0);
934 
935  for (MFIter mfi(mf_flux); mfi.isValid(); ++mfi) {
936  const auto& vbx = mfi.validbox();
937  const int nx = vbx.length(0);
938  const int imin = vbx.smallEnd(0);
939  const int jmin = vbx.smallEnd(1);
940  const int offset = m_col_offsets[mfi.index()];
941  const Array4<Real>& dst_arr = mf_flux.array(mfi);
942  ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
943  {
944  // map [i,j,k] 0-based to [icol, ilay] 0-based
945  const int icol = (j-jmin)*nx + (i-imin) + offset;
946  const int ilay = k;
947 
948  // SW and LW fluxes
949  dst_arr(i,j,k,0) = sw_flux_up_tab(icol,ilay);
950  dst_arr(i,j,k,1) = sw_flux_dn_tab(icol,ilay);
951  dst_arr(i,j,k,2) = sw_flux_dn_dir_tab(icol,ilay);
952  dst_arr(i,j,k,3) = lw_flux_up_tab(icol,ilay);
953  dst_arr(i,j,k,4) = lw_flux_dn_tab(icol,ilay);
954  });
955  }
956 
957 
958  std::string plotfilename = amrex::Concatenate("plt_rad", m_step, 5);
959  Vector<std::string> flux_names = {"sw_flux_up", "sw_flux_dn", "sw_flux_dir",
960  "lw_flux_up", "lw_flux_dn"};
961  WriteSingleLevelPlotfile(plotfilename, mf_flux, flux_names, m_geom, static_cast<Real>(m_time), m_step);
962 }
Here is the call graph for this function:

◆ WriteDataLog()

void Radiation::WriteDataLog ( const double &  time)
overridevirtual

Implements IRadiation.

1071 {
1072  constexpr int datwidth = 14;
1073  constexpr int datprecision = 9;
1074  constexpr int timeprecision = 13;
1075 
1076  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;
1077  // Clear sky
1078  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;
1079  // Clean sky
1080  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;
1081  // Clean clear sky
1082  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;
1083 
1084 
1085  auto domain = m_geom.Domain();
1086  h_avg_radqrsw = sumToLine(datalog_mf, 0, 1, domain, 2);
1087  h_avg_radqrlw = sumToLine(datalog_mf, 1, 1, domain, 2);
1088  h_avg_sw_up = sumToLine(datalog_mf, 2, 1, domain, 2);
1089  h_avg_sw_dn = sumToLine(datalog_mf, 3, 1, domain, 2);
1090  h_avg_sw_dn_dir = sumToLine(datalog_mf, 4, 1, domain, 2);
1091  h_avg_lw_up = sumToLine(datalog_mf, 5, 1, domain, 2);
1092  h_avg_lw_dn = sumToLine(datalog_mf, 6, 1, domain, 2);
1093  h_avg_zenith = sumToLine(datalog_mf, 7, 1, domain, 2);
1094 
1095  h_avg_radqrcsw = sumToLine(datalog_mf, 8, 1, domain, 2);
1096  h_avg_radqrclw = sumToLine(datalog_mf, 9, 1, domain, 2);
1097  h_avg_sw_clr_up = sumToLine(datalog_mf, 10, 1, domain, 2);
1098  h_avg_sw_clr_dn = sumToLine(datalog_mf, 11, 1, domain, 2);
1099  h_avg_sw_clr_dn_dir = sumToLine(datalog_mf, 12, 1, domain, 2);
1100  h_avg_lw_clr_up = sumToLine(datalog_mf, 13, 1, domain, 2);
1101  h_avg_lw_clr_dn = sumToLine(datalog_mf, 14, 1, domain, 2);
1102 
1103  if (m_extra_clnsky_diag) {
1104  h_avg_sw_cln_up = sumToLine(datalog_mf, 15, 1, domain, 2);
1105  h_avg_sw_cln_dn = sumToLine(datalog_mf, 16, 1, domain, 2);
1106  h_avg_sw_cln_dn_dir = sumToLine(datalog_mf, 17, 1, domain, 2);
1107  h_avg_lw_cln_up = sumToLine(datalog_mf, 18, 1, domain, 2);
1108  h_avg_lw_cln_dn = sumToLine(datalog_mf, 19, 1, domain, 2);
1109  }
1110 
1111  if (m_extra_clnclrsky_diag) {
1112  h_avg_sw_clnclr_up = sumToLine(datalog_mf, 20, 1, domain, 2);
1113  h_avg_sw_clnclr_dn = sumToLine(datalog_mf, 21, 1, domain, 2);
1114  h_avg_sw_clnclr_dn_dir = sumToLine(datalog_mf, 22, 1, domain, 2);
1115  h_avg_lw_clnclr_up = sumToLine(datalog_mf, 23, 1, domain, 2);
1116  h_avg_lw_clnclr_dn = sumToLine(datalog_mf, 24, 1, domain, 2);
1117  }
1118 
1119  Real area_z = static_cast<Real>(domain.length(0)*domain.length(1));
1120  int nz = domain.length(2);
1121  for (int k = 0; k < nz; k++) {
1122  h_avg_radqrsw[k] /= area_z;
1123  h_avg_radqrlw[k] /= area_z;
1124  h_avg_sw_up[k] /= area_z;
1125  h_avg_sw_dn[k] /= area_z;
1126  h_avg_sw_dn_dir[k] /= area_z;
1127  h_avg_lw_up[k] /= area_z;
1128  h_avg_lw_dn[k] /= area_z;
1129  h_avg_zenith[k] /= area_z;
1130 
1131  h_avg_radqrcsw[k] /= area_z;
1132  h_avg_radqrclw[k] /= area_z;
1133  h_avg_sw_clr_up[k] /= area_z;
1134  h_avg_sw_clr_dn[k] /= area_z;
1135  h_avg_sw_clr_dn_dir[k] /= area_z;
1136  h_avg_lw_clr_up[k] /= area_z;
1137  h_avg_lw_clr_dn[k] /= area_z;
1138  }
1139 
1140  if (m_extra_clnsky_diag) {
1141  for (int k = 0; k < nz; k++) {
1142  h_avg_sw_cln_up[k] /= area_z;
1143  h_avg_sw_cln_dn[k] /= area_z;
1144  h_avg_sw_cln_dn_dir[k] /= area_z;
1145  h_avg_lw_cln_up[k] /= area_z;
1146  h_avg_lw_cln_dn[k] /= area_z;
1147  }
1148  }
1149 
1150  if (m_extra_clnclrsky_diag) {
1151  for (int k = 0; k < nz; k++) {
1152  h_avg_sw_clnclr_up[k] /= area_z;
1153  h_avg_sw_clnclr_dn[k] /= area_z;
1154  h_avg_sw_clnclr_dn_dir[k] /= area_z;
1155  h_avg_lw_clnclr_up[k] /= area_z;
1156  h_avg_lw_clnclr_dn[k] /= area_z;
1157  }
1158  }
1159 
1160  if (ParallelDescriptor::IOProcessor()) {
1161  std::ostream& log = *datalog;
1162  if (log.good()) {
1163 
1164  for (int k = 0; k < nz; k++)
1165  {
1166  Real z = k * m_geom.CellSize(2);
1167  log << std::setw(datwidth) << std::setprecision(timeprecision) << time << " "
1168  << std::setw(datwidth) << std::setprecision(datprecision) << z << " "
1169  << h_avg_radqrsw[k] << " " << h_avg_radqrlw[k] << " " << h_avg_sw_up[k] << " "
1170  << h_avg_sw_dn[k] << " " << h_avg_sw_dn_dir[k] << " " << h_avg_lw_up[k] << " "
1171  << h_avg_lw_dn[k] << " " << h_avg_zenith[k] << " "
1172  << h_avg_radqrcsw[k] << " " << h_avg_radqrclw[k] << " " << h_avg_sw_clr_up[k] << " "
1173  << h_avg_sw_clr_dn[k] << " " << h_avg_sw_clr_dn_dir[k] << " " << h_avg_lw_clr_up[k] << " "
1174  << h_avg_lw_clr_dn[k] << " ";
1175  if (m_extra_clnsky_diag) {
1176  log << h_avg_sw_cln_up[k] << " " << h_avg_sw_cln_dn[k] << " " << h_avg_sw_cln_dn_dir[k] << " "
1177  << h_avg_lw_cln_up[k] << " " << h_avg_lw_cln_dn[k] << " ";
1178  } else {
1179  log << zero << " " << zero << " " << zero << " " << zero << " " << zero << " ";
1180  }
1181 
1182  if (m_extra_clnclrsky_diag) {
1183  log << h_avg_sw_clnclr_up[k] << " " << h_avg_sw_clnclr_dn[k] << " " << h_avg_sw_clnclr_dn_dir[k] << " "
1184  << h_avg_lw_clnclr_up[k] << " " << h_avg_lw_clnclr_dn[k] << std::endl;
1185  } else {
1186  log << zero << " " << zero << " " << zero << " " << zero << " " << zero << std::endl;
1187  }
1188  }
1189  // Write top face values
1190  Real z = nz * m_geom.CellSize(2);
1191  log << std::setw(datwidth) << std::setprecision(timeprecision) << time << " "
1192  << std::setw(datwidth) << std::setprecision(datprecision) << z << " "
1193  << zero << " " << zero << " " << zero << " " << zero << " " << zero << " " << zero << " "
1194  << zero << " " << zero << " "
1195  << zero << " " << zero << " " << zero << " " << zero << " " << zero << " " << zero << " "
1196  << zero << " "
1197  << zero << " " << zero << " " << zero << " " << zero << " " << zero << " "
1198  << zero << " " << zero << " " << zero << " " << zero << " " << zero
1199  << std::endl;
1200  }
1201  }
1202 }
std::unique_ptr< std::fstream > datalog
Definition: ERF_RadiationInterface.H:102

Member Data Documentation

◆ aero_g_sw

real3d_k Radiation::aero_g_sw
private

◆ aero_ssa_sw

real3d_k Radiation::aero_ssa_sw
private

◆ aero_tau_lw

real3d_k Radiation::aero_tau_lw
private

◆ aero_tau_sw

real3d_k Radiation::aero_tau_sw
private

◆ cldfrac_tot

real2d_k Radiation::cldfrac_tot
private

◆ d_tint

real2d_k Radiation::d_tint
private

◆ datalog_mf

amrex::MultiFab Radiation::datalog_mf
private

◆ eff_radius_qc

real2d_k Radiation::eff_radius_qc
private

◆ eff_radius_qi

real2d_k Radiation::eff_radius_qi
private

◆ gas_names_offset

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

◆ iwp

real2d_k Radiation::iwp
private

◆ lat

real1d_k Radiation::lat
private

◆ lon

real1d_k Radiation::lon
private

◆ lw_bnd_flux_dn

real3d_k Radiation::lw_bnd_flux_dn
private

◆ lw_bnd_flux_up

real3d_k Radiation::lw_bnd_flux_up
private

◆ lw_clnclrsky_flux_dn

real2d_k Radiation::lw_clnclrsky_flux_dn
private

◆ lw_clnclrsky_flux_up

real2d_k Radiation::lw_clnclrsky_flux_up
private

◆ lw_clnsky_flux_dn

real2d_k Radiation::lw_clnsky_flux_dn
private

◆ lw_clnsky_flux_up

real2d_k Radiation::lw_clnsky_flux_up
private

◆ lw_clrsky_flux_dn

real2d_k Radiation::lw_clrsky_flux_dn
private

◆ lw_clrsky_flux_up

real2d_k Radiation::lw_clrsky_flux_up
private

◆ lw_clrsky_heating

real2d_k Radiation::lw_clrsky_heating
private

◆ lw_flux_dn

real2d_k Radiation::lw_flux_dn
private

◆ lw_flux_up

real2d_k Radiation::lw_flux_up
private

◆ lw_heating

real2d_k Radiation::lw_heating
private

◆ lw_src

real1d_k Radiation::lw_src
private

◆ lwp

real2d_k Radiation::lwp
private

◆ m_ba

amrex::BoxArray Radiation::m_ba
private

◆ m_calday

amrex::Real Radiation::m_calday = 0.0
private

Referenced by Run().

◆ m_ch4vmr

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

◆ m_co2vmr

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

◆ m_col_offsets

amrex::Vector<int> Radiation::m_col_offsets
private

Referenced by Init().

◆ m_cons_in

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

◆ m_covmr

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

◆ m_declin

amrex::Real Radiation::m_declin = 0.0
private

Referenced by Run().

◆ m_do_aerosol_rad

bool Radiation::m_do_aerosol_rad = false
private

◆ m_do_subcol_sampling

bool Radiation::m_do_subcol_sampling = true
private

◆ m_dt

double Radiation::m_dt
private

◆ m_extra_clnclrsky_diag

bool Radiation::m_extra_clnclrsky_diag = false
private

◆ m_extra_clnsky_diag

bool Radiation::m_extra_clnsky_diag = false
private

◆ m_fixed_orbital_year

bool Radiation::m_fixed_orbital_year = false
private

◆ m_fixed_solar_zenith_angle

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

◆ m_fixed_total_solar_irradiance

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

◆ m_gas_concs

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

◆ m_gas_mol_weights

real1d_k Radiation::m_gas_mol_weights
private

◆ m_gas_names

const std::vector<std::string> Radiation::m_gas_names
private
Initial value:
= {"H2O", "CO2", "O3", "N2O",
"CO" , "CH4", "O2", "N2" }

◆ m_geom

amrex::Geometry Radiation::m_geom
private

◆ m_ice

bool Radiation::m_ice = false
private

◆ m_is_nested_patch

bool Radiation::m_is_nested_patch = false
private

◆ m_lat

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

◆ m_lat_cons

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

◆ m_lev

int Radiation::m_lev
private

Referenced by rad_run_impl().

◆ m_lon

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

◆ m_lon_cons

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

◆ m_lsm

bool Radiation::m_lsm = false
private

◆ m_lsm_input_names

amrex::Vector<std::string> Radiation::m_lsm_input_names
private
Initial value:
= {"t_sfc" , "sfc_emis" ,
"sfc_alb_dir_vis", "sfc_alb_dir_nir",
"sfc_alb_dif_vis", "sfc_alb_dif_nir"}

Referenced by get_lsm_input_varnames().

◆ m_lsm_output_names

amrex::Vector<std::string> Radiation::m_lsm_output_names
private
Initial value:
= {"cos_zenith_angle" , "sw_flux_dn" ,
"sw_flux_dn_dir_vis", "sw_flux_dn_dir_nir",
"sw_flux_dn_dif_vis", "sw_flux_dn_dif_nir",
"lw_flux_dn"}

Referenced by get_lsm_output_varnames().

◆ m_moist

bool Radiation::m_moist = false
private

◆ m_mol_weight_gas

const std::vector<amrex::Real> Radiation::m_mol_weight_gas
private
Initial value:
= {amrex::Real(18.01528), amrex::Real(44.00950), amrex::Real(47.9982), amrex::Real(44.0128),
amrex::Real(28.01010), amrex::Real(16.04246), amrex::Real(31.9980), amrex::Real(28.0134)}

◆ m_n2ovmr

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

◆ m_n2vmr

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

◆ m_ncol

int Radiation::m_ncol
private

Referenced by Init().

◆ m_ncol_chunk

int Radiation::m_ncol_chunk = 1024
private

Referenced by Init().

◆ m_ncol_chunk_requested

int Radiation::m_ncol_chunk_requested = 1024
private

Referenced by Init().

◆ m_ngas

int Radiation::m_ngas = 8
private

◆ m_nlay

int Radiation::m_nlay
private

Referenced by Init(), and Run().

◆ m_nlwbands

int Radiation::m_nlwbands
private

◆ m_nlwgpts

int Radiation::m_nlwgpts
private

◆ m_nswbands

int Radiation::m_nswbands
private

◆ m_nswgpts

int Radiation::m_nswgpts
private

◆ m_o2vmr

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

◆ m_o3_size

int Radiation::m_o3_size
private

◆ m_o3vmr

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

◆ m_orbital_day

int Radiation::m_orbital_day = -9999
private

Referenced by rad_run_impl().

◆ m_orbital_eccen

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

◆ m_orbital_mon

int Radiation::m_orbital_mon = -9999
private

Referenced by rad_run_impl().

◆ m_orbital_mvelp

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

◆ m_orbital_obliq

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

◆ m_orbital_sec

int Radiation::m_orbital_sec = -9999
private

Referenced by rad_run_impl().

◆ m_orbital_year

int Radiation::m_orbital_year = -9999
private

Referenced by rad_run_impl().

◆ m_qheating_rates

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

◆ m_qi_comp

int Radiation::m_qi_comp = -1
private

◆ m_rad_fluxes

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

◆ m_rad_freq_in_steps

int Radiation::m_rad_freq_in_steps = 1
private

◆ m_rad_nvar

int Radiation::m_rad_nvar = 12
private

◆ m_rad_t_sfc

amrex::Real Radiation::m_rad_t_sfc = -1
private

◆ m_rad_write_fluxes

bool Radiation::m_rad_write_fluxes = false
private

◆ m_rdOcp

amrex::Real Radiation::m_rdOcp = RdoCp
private

◆ m_step

int Radiation::m_step
private

◆ m_time

double Radiation::m_time
private

◆ m_update_rad

bool Radiation::m_update_rad = false
private

Referenced by rad_run_impl(), and radiation_updated().

◆ m_z_phys

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

◆ mu0

real1d_k Radiation::mu0
private

◆ o3_lay

real1d_k Radiation::o3_lay
private

◆ p_lay

real2d_k Radiation::p_lay
private

◆ p_lev

real2d_k Radiation::p_lev
private

◆ qc_lay

real2d_k Radiation::qc_lay
private

◆ qi_lay

real2d_k Radiation::qi_lay
private

◆ qv_lay

real2d_k Radiation::qv_lay
private

◆ r_lay

real2d_k Radiation::r_lay
private

◆ rrtmgp_cloud_optics_file_lw

std::string Radiation::rrtmgp_cloud_optics_file_lw
private

◆ rrtmgp_cloud_optics_file_sw

std::string Radiation::rrtmgp_cloud_optics_file_sw
private

◆ rrtmgp_cloud_optics_lw

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

◆ rrtmgp_cloud_optics_sw

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

◆ rrtmgp_coeffs_file_lw

std::string Radiation::rrtmgp_coeffs_file_lw
private

◆ rrtmgp_coeffs_file_sw

std::string Radiation::rrtmgp_coeffs_file_sw
private

◆ rrtmgp_coeffs_lw

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

◆ rrtmgp_coeffs_sw

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

◆ rrtmgp_file_path

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

◆ sfc_alb_dif

real2d_k Radiation::sfc_alb_dif
private

◆ sfc_alb_dif_nir

real1d_k Radiation::sfc_alb_dif_nir
private

◆ sfc_alb_dif_vis

real1d_k Radiation::sfc_alb_dif_vis
private

◆ sfc_alb_dir

real2d_k Radiation::sfc_alb_dir
private

◆ sfc_alb_dir_nir

real1d_k Radiation::sfc_alb_dir_nir
private

◆ sfc_alb_dir_vis

real1d_k Radiation::sfc_alb_dir_vis
private

◆ sfc_emis

real1d_k Radiation::sfc_emis
private

◆ sfc_flux_dif_nir

real1d_k Radiation::sfc_flux_dif_nir
private

◆ sfc_flux_dif_vis

real1d_k Radiation::sfc_flux_dif_vis
private

◆ sfc_flux_dir_nir

real1d_k Radiation::sfc_flux_dir_nir
private

◆ sfc_flux_dir_vis

real1d_k Radiation::sfc_flux_dir_vis
private

◆ sw_bnd_flux_dif

real3d_k Radiation::sw_bnd_flux_dif
private

◆ sw_bnd_flux_dir

real3d_k Radiation::sw_bnd_flux_dir
private

◆ sw_bnd_flux_dn

real3d_k Radiation::sw_bnd_flux_dn
private

◆ sw_bnd_flux_up

real3d_k Radiation::sw_bnd_flux_up
private

◆ sw_clnclrsky_flux_dn

real2d_k Radiation::sw_clnclrsky_flux_dn
private

◆ sw_clnclrsky_flux_dn_dir

real2d_k Radiation::sw_clnclrsky_flux_dn_dir
private

◆ sw_clnclrsky_flux_up

real2d_k Radiation::sw_clnclrsky_flux_up
private

◆ sw_clnsky_flux_dn

real2d_k Radiation::sw_clnsky_flux_dn
private

◆ sw_clnsky_flux_dn_dir

real2d_k Radiation::sw_clnsky_flux_dn_dir
private

◆ sw_clnsky_flux_up

real2d_k Radiation::sw_clnsky_flux_up
private

◆ sw_clrsky_flux_dn

real2d_k Radiation::sw_clrsky_flux_dn
private

◆ sw_clrsky_flux_dn_dir

real2d_k Radiation::sw_clrsky_flux_dn_dir
private

◆ sw_clrsky_flux_up

real2d_k Radiation::sw_clrsky_flux_up
private

◆ sw_clrsky_heating

real2d_k Radiation::sw_clrsky_heating
private

◆ sw_flux_dn

real2d_k Radiation::sw_flux_dn
private

◆ sw_flux_dn_dir

real2d_k Radiation::sw_flux_dn_dir
private

◆ sw_flux_up

real2d_k Radiation::sw_flux_up
private

◆ sw_heating

real2d_k Radiation::sw_heating
private

◆ t_lay

real2d_k Radiation::t_lay
private

◆ t_lev

real2d_k Radiation::t_lev
private

◆ t_sfc

real1d_k Radiation::t_sfc
private

◆ z_del

real2d_k Radiation::z_del
private

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