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, amrex::Vector< amrex::MultiFab * > &lsm_input_ptrs, amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs, amrex::MultiFab *qheating_rates, amrex::MultiFab *rad_fluxes, amrex::MultiFab *z_phys, amrex::MultiFab *lat_ptr, amrex::MultiFab *lon_ptr, const bool updated_lsm) override
 
void set_grids (int &level, int &step, double &time, const double &dt, const amrex::BoxArray &ba, amrex::Geometry &geom, amrex::MultiFab *cons_in, amrex::iMultiFab *lmask, amrex::MultiFab *t_surf, amrex::Vector< amrex::MultiFab * > &lsm_input_ptrs, amrex::MultiFab *qheating_rates, amrex::MultiFab *rad_fluxes, amrex::MultiFab *z_phys, amrex::MultiFab *lat, amrex::MultiFab *lon, const bool updated_lsm)
 
void alloc_buffers ()
 
void dealloc_buffers ()
 
void mf_to_kokkos_buffers (amrex::iMultiFab *lmask, amrex::MultiFab *t_surf, amrex::Vector< amrex::MultiFab * > &lsm_input_ptrs)
 
void kokkos_buffers_to_mf (amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
 
void write_rrtmgp_fluxes ()
 
void initialize_impl ()
 
void run_impl ()
 
void finalize_impl (amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
 
void rad_run_impl (amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
 
virtual amrex::Vector< std::string > get_lsm_input_varnames () override
 
virtual amrex::Vector< std::string > get_lsm_output_varnames () override
 
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::Realm_mol_weight_gas
 
amrex::Real m_co2vmr = amrex::Real(388.717e-6)
 
amrex::Vector< amrex::Realm_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, KokkosDefaultDevicem_gas_concs
 
int m_ncol
 
int m_nlay
 
amrex::Vector< int > m_col_offsets
 
bool m_do_aerosol_rad = false
 
bool m_extra_clnsky_diag = false
 
bool m_extra_clnclrsky_diag = false
 
int m_orbital_year = -9999
 
int m_orbital_mon = -9999
 
int m_orbital_day = -9999
 
int m_orbital_sec = -9999
 
bool m_fixed_orbital_year = false
 
amrex::Real m_orbital_eccen = -amrex::Real(9999.)
 
amrex::Real m_orbital_obliq = -amrex::Real(9999.)
 
amrex::Real m_orbital_mvelp = -amrex::Real(9999.)
 
amrex::Real m_fixed_total_solar_irradiance = -amrex::Real(9999.)
 
amrex::Real m_fixed_solar_zenith_angle = -amrex::Real(9999.)
 
int m_nswbands
 
int m_nlwbands
 
int m_nswgpts
 
int m_nlwgpts
 
int m_rad_freq_in_steps = 1
 
int m_ncol_chunk_requested = 1024
 
int m_ncol_chunk = 1024
 
int m_rad_nvar = 12
 
bool m_do_subcol_sampling = true
 
real1d_k o3_lay
 
real1d_k mu0
 
real1d_k sfc_alb_dir_vis
 
real1d_k sfc_alb_dir_nir
 
real1d_k sfc_alb_dif_vis
 
real1d_k sfc_alb_dif_nir
 
real1d_k sfc_flux_dir_vis
 
real1d_k sfc_flux_dir_nir
 
real1d_k sfc_flux_dif_vis
 
real1d_k sfc_flux_dif_nir
 
real1d_k lat
 
real1d_k lon
 
real1d_k sfc_emis
 
real1d_k t_sfc
 
real1d_k lw_src
 
real2d_k r_lay
 
real2d_k p_lay
 
real2d_k t_lay
 
real2d_k z_del
 
real2d_k qv_lay
 
real2d_k qc_lay
 
real2d_k qi_lay
 
real2d_k cldfrac_tot
 
real2d_k eff_radius_qc
 
real2d_k eff_radius_qi
 
real2d_k lwp
 
real2d_k iwp
 
real2d_k sw_heating
 
real2d_k lw_heating
 
real2d_k sw_clrsky_heating
 
real2d_k lw_clrsky_heating
 
real2d_k d_tint
 
real2d_k p_lev
 
real2d_k t_lev
 
real2d_k sw_flux_up
 
real2d_k sw_flux_dn
 
real2d_k sw_flux_dn_dir
 
real2d_k lw_flux_up
 
real2d_k lw_flux_dn
 
real2d_k sw_clnclrsky_flux_up
 
real2d_k sw_clnclrsky_flux_dn
 
real2d_k sw_clnclrsky_flux_dn_dir
 
real2d_k sw_clrsky_flux_up
 
real2d_k sw_clrsky_flux_dn
 
real2d_k sw_clrsky_flux_dn_dir
 
real2d_k sw_clnsky_flux_up
 
real2d_k sw_clnsky_flux_dn
 
real2d_k sw_clnsky_flux_dn_dir
 
real2d_k lw_clnclrsky_flux_up
 
real2d_k lw_clnclrsky_flux_dn
 
real2d_k lw_clrsky_flux_up
 
real2d_k lw_clrsky_flux_dn
 
real2d_k lw_clnsky_flux_up
 
real2d_k lw_clnsky_flux_dn
 
real3d_k sw_bnd_flux_up
 
real3d_k sw_bnd_flux_dn
 
real3d_k sw_bnd_flux_dir
 
real3d_k sw_bnd_flux_dif
 
real3d_k lw_bnd_flux_up
 
real3d_k lw_bnd_flux_dn
 
real2d_k sfc_alb_dir
 
real2d_k sfc_alb_dif
 
real3d_k aero_tau_sw
 
real3d_k aero_ssa_sw
 
real3d_k aero_g_sw
 
real3d_k aero_tau_lw
 

Additional Inherited Members

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

Constructor & Destructor Documentation

◆ Radiation()

Radiation::Radiation ( const int &  lev,
SolverChoice sc 
)
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 5000). "
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:350
std::string rrtmgp_coeffs_sw
Definition: ERF_Radiation.H:346
int m_rad_freq_in_steps
Definition: ERF_Radiation.H:432
amrex::Real m_lon_cons
Definition: ERF_Radiation.H:339
int m_nswbands
Definition: ERF_Radiation.H:426
bool m_do_aerosol_rad
Definition: ERF_Radiation.H:393
std::string rrtmgp_coeffs_lw
Definition: ERF_Radiation.H:347
bool m_moist
Definition: ERF_Radiation.H:299
std::string rrtmgp_cloud_optics_file_lw
Definition: ERF_Radiation.H:353
int m_ncol_chunk_requested
Definition: ERF_Radiation.H:440
bool m_lsm
Definition: ERF_Radiation.H:304
bool m_do_subcol_sampling
Definition: ERF_Radiation.H:447
bool m_rad_write_fluxes
Definition: ERF_Radiation.H:296
int m_qi_comp
Definition: ERF_Radiation.H:301
std::string rrtmgp_cloud_optics_file_sw
Definition: ERF_Radiation.H:352
amrex::Real m_orbital_mvelp
Definition: ERF_Radiation.H:414
amrex::Vector< amrex::Real > m_o3vmr
Definition: ERF_Radiation.H:364
std::string rrtmgp_cloud_optics_sw
Definition: ERF_Radiation.H:348
int m_nlwgpts
Definition: ERF_Radiation.H:429
std::string rrtmgp_cloud_optics_lw
Definition: ERF_Radiation.H:349
amrex::Real m_co2vmr
Definition: ERF_Radiation.H:363
bool m_extra_clnsky_diag
Definition: ERF_Radiation.H:396
amrex::Real m_lat_cons
Definition: ERF_Radiation.H:338
amrex::Real m_fixed_total_solar_irradiance
Definition: ERF_Radiation.H:419
amrex::Real m_o2vmr
Definition: ERF_Radiation.H:368
amrex::Real m_n2ovmr
Definition: ERF_Radiation.H:365
amrex::Real m_rad_t_sfc
Definition: ERF_Radiation.H:318
std::string rrtmgp_file_path
Definition: ERF_Radiation.H:345
int m_nlwbands
Definition: ERF_Radiation.H:427
amrex::Real m_orbital_eccen
Definition: ERF_Radiation.H:412
amrex::Real m_covmr
Definition: ERF_Radiation.H:366
int m_ncol_chunk
Definition: ERF_Radiation.H:441
std::string rrtmgp_coeffs_file_lw
Definition: ERF_Radiation.H:351
amrex::Real m_ch4vmr
Definition: ERF_Radiation.H:367
bool m_ice
Definition: ERF_Radiation.H:300
amrex::Real m_fixed_solar_zenith_angle
Definition: ERF_Radiation.H:423
amrex::Real m_orbital_obliq
Definition: ERF_Radiation.H:413
amrex::Real m_n2vmr
Definition: ERF_Radiation.H:369
bool m_fixed_orbital_year
Definition: ERF_Radiation.H:411
int m_nswgpts
Definition: ERF_Radiation.H:428
bool m_extra_clnclrsky_diag
Definition: ERF_Radiation.H:397
amrex::Real m_rdOcp
Definition: ERF_Radiation.H:319
int m_rad_nvar
Definition: ERF_Radiation.H:444
int m_orbital_year
Definition: ERF_Radiation.H:402
int qi
cloud ice
Definition: ERF_DataStruct.H:236
MoistureType moisture_type
Moisture or microphysics model.
Definition: ERF_DataStruct.H:2237
amrex::Real rdOcp
Ratio of dry-air gas constant to c_p.
Definition: ERF_DataStruct.H:2062
LandSurfaceType lsm_type
Land-surface model.
Definition: ERF_DataStruct.H:2240
MoistureComponentIndices moisture_indices
Index map of the moisture data carried by the active scheme: conserved-state components for the speci...
Definition: ERF_DataStruct.H:2263
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:97
real3d_k sw_bnd_flux_dn
Definition: ERF_Radiation.H:513
real2d_k lw_flux_up
Definition: ERF_Radiation.H:493
real3d_k aero_tau_sw
Definition: ERF_Radiation.H:532
int m_o3_size
Definition: ERF_Radiation.H:373
real3d_k sw_bnd_flux_dir
Definition: ERF_Radiation.H:514
real2d_k sw_clnsky_flux_dn
Definition: ERF_Radiation.H:502
real2d_k d_tint
Definition: ERF_Radiation.H:487
real2d_k lw_clnclrsky_flux_dn
Definition: ERF_Radiation.H:505
real2d_k lwp
Definition: ERF_Radiation.H:479
real1d_k lw_src
Definition: ERF_Radiation.H:466
real2d_k eff_radius_qi
Definition: ERF_Radiation.H:478
real1d_k m_gas_mol_weights
Definition: ERF_Radiation.H:374
real2d_k sw_heating
Definition: ERF_Radiation.H:481
real3d_k lw_bnd_flux_dn
Definition: ERF_Radiation.H:519
real3d_k sw_bnd_flux_up
Definition: ERF_Radiation.H:512
real2d_k qv_lay
Definition: ERF_Radiation.H:473
real2d_k sw_clnclrsky_flux_dn_dir
Definition: ERF_Radiation.H:497
real1d_k sfc_flux_dif_vis
Definition: ERF_Radiation.H:460
real2d_k lw_clnclrsky_flux_up
Definition: ERF_Radiation.H:504
real1d_k lat
Definition: ERF_Radiation.H:462
real2d_k qi_lay
Definition: ERF_Radiation.H:475
real2d_k sw_clrsky_flux_up
Definition: ERF_Radiation.H:498
real2d_k t_lev
Definition: ERF_Radiation.H:489
real1d_k o3_lay
Definition: ERF_Radiation.H:450
real1d_k sfc_alb_dif_vis
Definition: ERF_Radiation.H:456
real1d_k mu0
Definition: ERF_Radiation.H:453
real2d_k cldfrac_tot
Definition: ERF_Radiation.H:476
real1d_k sfc_flux_dir_nir
Definition: ERF_Radiation.H:459
real2d_k sw_clnclrsky_flux_up
Definition: ERF_Radiation.H:495
real3d_k aero_g_sw
Definition: ERF_Radiation.H:534
real2d_k sw_flux_up
Definition: ERF_Radiation.H:490
real3d_k sw_bnd_flux_dif
Definition: ERF_Radiation.H:515
real2d_k sfc_alb_dif
Definition: ERF_Radiation.H:523
real1d_k sfc_alb_dif_nir
Definition: ERF_Radiation.H:457
real2d_k lw_clnsky_flux_dn
Definition: ERF_Radiation.H:509
real2d_k p_lay
Definition: ERF_Radiation.H:470
real2d_k r_lay
Definition: ERF_Radiation.H:469
int m_ncol
Definition: ERF_Radiation.H:383
real2d_k sw_flux_dn_dir
Definition: ERF_Radiation.H:492
const std::vector< amrex::Real > m_mol_weight_gas
Definition: ERF_Radiation.H:359
std::vector< std::string > gas_names_offset
Definition: ERF_Radiation.H:375
real2d_k sw_clrsky_flux_dn
Definition: ERF_Radiation.H:499
real3d_k aero_ssa_sw
Definition: ERF_Radiation.H:533
real2d_k sw_clnsky_flux_dn_dir
Definition: ERF_Radiation.H:503
real2d_k qc_lay
Definition: ERF_Radiation.H:474
real2d_k lw_clrsky_flux_up
Definition: ERF_Radiation.H:506
real2d_k sw_flux_dn
Definition: ERF_Radiation.H:491
real2d_k lw_clrsky_flux_dn
Definition: ERF_Radiation.H:507
real2d_k z_del
Definition: ERF_Radiation.H:472
real3d_k aero_tau_lw
Definition: ERF_Radiation.H:535
real2d_k lw_flux_dn
Definition: ERF_Radiation.H:494
real2d_k sw_clrsky_flux_dn_dir
Definition: ERF_Radiation.H:500
real1d_k sfc_flux_dif_nir
Definition: ERF_Radiation.H:461
int m_ngas
Definition: ERF_Radiation.H:356
real2d_k lw_heating
Definition: ERF_Radiation.H:482
real1d_k sfc_alb_dir_nir
Definition: ERF_Radiation.H:455
real2d_k lw_clnsky_flux_up
Definition: ERF_Radiation.H:508
real2d_k sfc_alb_dir
Definition: ERF_Radiation.H:522
real1d_k lon
Definition: ERF_Radiation.H:463
real1d_k sfc_emis
Definition: ERF_Radiation.H:464
real2d_k sw_clnsky_flux_up
Definition: ERF_Radiation.H:501
int m_nlay
Definition: ERF_Radiation.H:384
real3d_k lw_bnd_flux_up
Definition: ERF_Radiation.H:518
real1d_k sfc_alb_dir_vis
Definition: ERF_Radiation.H:454
const std::vector< std::string > m_gas_names
Definition: ERF_Radiation.H:357
real2d_k sw_clrsky_heating
Definition: ERF_Radiation.H:483
real2d_k t_lay
Definition: ERF_Radiation.H:471
real1d_k sfc_flux_dir_vis
Definition: ERF_Radiation.H:458
real2d_k iwp
Definition: ERF_Radiation.H:480
real2d_k p_lev
Definition: ERF_Radiation.H:488
real2d_k lw_clrsky_heating
Definition: ERF_Radiation.H:484
real2d_k eff_radius_qc
Definition: ERF_Radiation.H:477
real1d_k t_sfc
Definition: ERF_Radiation.H:465
real2d_k sw_clnclrsky_flux_dn
Definition: ERF_Radiation.H:496
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 ( amrex::Vector< amrex::MultiFab * > &  lsm_output_ptrs)
1521 {
1522  // Reset gas concentrations (k-dist data persists across steps)
1523  m_gas_concs.reset();
1524 
1525  // Fill the AMReX MFs from Kokkos Views
1526  kokkos_buffers_to_mf(lsm_output_ptrs);
1527 
1528  // Write fluxes if requested
1530 
1531  // Fill output data for datalog before deallocating
1532  if (datalog_int > 0) {
1535  Kokkos::fence();
1537  }
1538 
1539  // Deallocate the buffer arrays
1540  dealloc_buffers();
1541 }
void dealloc_buffers()
Definition: ERF_Radiation.cpp:490
void populateDatalogMF()
Definition: ERF_Radiation.cpp:954
void kokkos_buffers_to_mf(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:377
void write_rrtmgp_fluxes()
Definition: ERF_Radiation.cpp:914
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.

239  {
240  return m_lsm_input_names;
241  }
amrex::Vector< std::string > m_lsm_input_names
Definition: ERF_Radiation.H:307

◆ get_lsm_output_varnames()

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

Reimplemented from IRadiation.

247  {
248  return m_lsm_output_names;
249  }
amrex::Vector< std::string > m_lsm_output_names
Definition: ERF_Radiation.H:312

◆ 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:252
bool m_is_nested_patch
Definition: ERF_Radiation.H:272
amrex::Vector< int > m_col_offsets
Definition: ERF_Radiation.H:387
Here is the call graph for this function:

◆ initialize_impl()

void Radiation::initialize_impl ( )
1196 {
1197  // Initialize gas concentrations for this step
1199 
1200  // Load k-distribution and cloud optics data only once.
1201  // These are static lookup tables that never change.
1202  // Size the memory pool for the requested chunk size (not the effective one, and
1203  // not min with the current m_ncol) so that the pool remains valid even if m_ncol
1204  // grows after regridding/load balancing. The pool is created once and never
1205  // resized, whereas the effective chunk size is recomputed at every Init() and is
1206  // bounded above by the request.
1207  if (!rrtmgp::initialized) {
1208  gas_concs_t gas_concs_pool;
1209  gas_concs_pool.init(gas_names_offset, m_ncol_chunk_requested, m_nlay);
1210  rrtmgp::rrtmgp_initialize(gas_concs_pool,
1213  m_rad_nvar);
1214  gas_concs_pool.reset();
1215  }
1216 }
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.

253  {
254  return m_is_nested_patch;
255  }

Referenced by Init().

Here is the caller graph for this function:

◆ kokkos_buffers_to_mf()

void Radiation::kokkos_buffers_to_mf ( amrex::Vector< amrex::MultiFab * > &  lsm_output_ptrs)
825 {
826  // Heating rate, fluxes, zenith, lsm ptrs
827 
828  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))});
829  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))});
830  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))});
831  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))});
832  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))});
833  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))});
834  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))});
835 
836  TableData<Real,1> sfc_flux_sw_dn; sfc_flux_sw_dn.resize({0}, {static_cast<int>(sw_flux_dn.extent(0))});
837  TableData<Real,1> sfc_flux_lw_dn; sfc_flux_lw_dn.resize({0}, {static_cast<int>(lw_flux_dn.extent(0))});
838  Table1D<Real> sfc_flux_sw_dn_tab = sfc_flux_sw_dn.table();
839  Table1D<Real> sfc_flux_lw_dn_tab = sfc_flux_lw_dn.table();
840  Table1D<Real> sfc_flux_sw_dir_vis_tab(sfc_flux_dir_vis.data(), {0}, {static_cast<int>(sfc_flux_dir_vis.extent(0))});
841  Table1D<Real> sfc_flux_sw_dir_nir_tab(sfc_flux_dir_nir.data(), {0}, {static_cast<int>(sfc_flux_dir_nir.extent(0))});
842  Table1D<Real> sfc_flux_sw_dif_vis_tab(sfc_flux_dif_vis.data(), {0}, {static_cast<int>(sfc_flux_dif_vis.extent(0))});
843  Table1D<Real> sfc_flux_sw_dif_nir_tab(sfc_flux_dif_nir.data(), {0}, {static_cast<int>(sfc_flux_dif_nir.extent(0))});
844  Table1D<Real> mu0_tab(mu0.data(), {0}, {static_cast<int>(mu0.extent(0))});
845  Vector<Table1D<Real>> rrtmgp_out_vars = {mu0_tab , sfc_flux_sw_dn_tab ,
846  sfc_flux_sw_dir_vis_tab, sfc_flux_sw_dir_nir_tab,
847  sfc_flux_sw_dif_vis_tab, sfc_flux_sw_dif_nir_tab,
848  sfc_flux_lw_dn_tab };
849 
850  const Real rdOcp = m_rdOcp;
851  for (MFIter mfi(*m_cons_in); mfi.isValid(); ++mfi) {
852  const auto& vbx = mfi.validbox();
853  const auto& sbx = makeSlab(vbx,2,vbx.smallEnd(2));
854  const int nx = vbx.length(0);
855  const int imin = vbx.smallEnd(0);
856  const int jmin = vbx.smallEnd(1);
857  const int ktop = vbx.bigEnd(2);
858  const int offset = m_col_offsets[mfi.index()];
859  const Array4<Real>& q_arr = m_qheating_rates->array(mfi);
860  const Array4<Real>& f_arr = m_rad_fluxes->array(mfi);
861  ParallelFor(vbx, [=]
862  AMREX_GPU_DEVICE (int i, int j, int k)
863  {
864  // map [i,j,k] 0-based to [icol, ilay] 0-based
865  const int icol = (j-jmin)*nx + (i-imin) + offset;
866  const int ilay = k;
867 
868  // Temperature heating rate for SW and LW
869  q_arr(i,j,k,0) = sw_heating_tab(icol,ilay);
870  q_arr(i,j,k,1) = lw_heating_tab(icol,ilay);
871 
872  // Convert the dT/dz to dTheta/dz
873  Real iexner = rrtmgp::inverse_exner(Real(p_lay_tab(icol,ilay)), rdOcp);
874  q_arr(i,j,k,0) *= iexner;
875  q_arr(i,j,k,1) *= iexner;
876 
877  // Populate the fluxes: level ilay is the lower interface of
878  // layer k, and the top-of-atmosphere level nlay goes into the
879  // z-ghost cell above the top layer (rad_fluxes has one).
880  f_arr(i,j,k,0) = sw_flux_up_tab(icol,ilay);
881  f_arr(i,j,k,1) = sw_flux_dn_tab(icol,ilay);
882  f_arr(i,j,k,2) = lw_flux_up_tab(icol,ilay);
883  f_arr(i,j,k,3) = lw_flux_dn_tab(icol,ilay);
884  if (k == ktop) {
885  f_arr(i,j,k+1,0) = sw_flux_up_tab(icol,ilay+1);
886  f_arr(i,j,k+1,1) = sw_flux_dn_tab(icol,ilay+1);
887  f_arr(i,j,k+1,2) = lw_flux_up_tab(icol,ilay+1);
888  f_arr(i,j,k+1,3) = lw_flux_dn_tab(icol,ilay+1);
889  }
890 
891  if (k==0) {
892  sfc_flux_sw_dn_tab(icol) = sw_flux_dn_tab(icol,ilay);
893  sfc_flux_lw_dn_tab(icol) = lw_flux_dn_tab(icol,ilay);
894  }
895  });
896  for (int ivar(0); ivar<lsm_output_ptrs.size(); ivar++) {
897  if (lsm_output_ptrs[ivar]) {
898  auto rrtmgp_for_fill = rrtmgp_out_vars[ivar];
899  const Array4<Real>& lsm_out_arr = lsm_output_ptrs[ivar]->array(mfi);
900  ParallelFor(sbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
901  {
902  // map [i,j,k] 0-based to [icol, ilay] 0-based
903  const int icol = (j-jmin)*nx + (i-imin) + offset;
904 
905  // export the desired variable at surface
906  lsm_out_arr(i,j,k) = rrtmgp_for_fill(icol);
907  });
908  } // valid ptr
909  } // ivar
910  }// mfi
911 }
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:322
amrex::MultiFab * m_rad_fluxes
Definition: ERF_Radiation.H:328
amrex::MultiFab * m_qheating_rates
Definition: ERF_Radiation.H:325
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,
amrex::Vector< 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< Real>& lsm_in_arr = (lsm_input_ptrs[ivar]) ? lsm_input_ptrs[ivar]->array(mfi) :
759  Array4< 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 && (lsm_in_arr(i,j,k) < lsm_undefined);
776  // Match TwoStream: convert SurfaceLayer theta with
777  // the physical surface pressure diagnosed from the
778  // lowest atmospheric cell.
780  is_land,
781  has_lsm_t_sfc,
782  valid_lsm_t_sfc,
783  has_lsm_t_sfc ? lsm_in_arr(i,j,k) : Real(0.),
784  static_cast<bool>(tsurf_arr),
785  tsurf_arr ? tsurf_arr(i,j,k) : Real(0.),
787  cons_arr(i,j,k_surface,Rho_comp),
788  cons_arr(i,j,k_surface,RhoTheta_comp),
789  moist ? std::max(cons_arr(i,j,k_surface,RhoQ1_comp) /
790  cons_arr(i,j,k_surface,Rho_comp), Real(0.)) : Real(0.),
791  z_arr ? Compute_Zrel_AtCellCenter(i,j,k_surface,z_arr) : Real(0.5)*dz),
792  rdOcp,
793  rrtmgp_default_val,
794  rrtmgp_to_fill(icol),
795  has_lsm_t_sfc ? &lsm_in_arr(i,j,k) : nullptr);
796  } else {
797  // Have LSM and are over land.
798  const bool valid_lsm_data =
799  (lsm_in_arr && (lsm_in_arr(i,j,k) < lsm_undefined));
800  if (is_land && valid_lsm_data) {
801  rrtmgp_to_fill(icol) = lsm_in_arr(i,j,k);
802  } else {
803  // Use the default value.
804  rrtmgp_to_fill(icol) = rrtmgp_default_val;
805  if (lsm_in_arr) { lsm_in_arr(i,j,k) = 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:26
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:39
#define RhoTheta_comp
Definition: ERF_IndexDefines.H:40
#define RhoQ2_comp
Definition: ERF_IndexDefines.H:46
#define RhoQ1_comp
Definition: ERF_IndexDefines.H:45
constexpr amrex::Real zero
Definition: ERF_NumericalConstants.H:29
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:751
amrex::MultiFab * m_lon
Definition: ERF_Radiation.H:335
amrex::MultiFab * m_z_phys
Definition: ERF_Radiation.H:331
amrex::Geometry m_geom
Definition: ERF_Radiation.H:287
amrex::MultiFab * m_lat
Definition: ERF_Radiation.H:334
@ 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:30
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, amrex::Real *lsm_t_sfc_out)
Definition: ERF_RRTMGP_SurfaceTemperature.H:40
Here is the call graph for this function:

◆ populateDatalogMF()

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

◆ rad_run_impl()

void Radiation::rad_run_impl ( amrex::Vector< amrex::MultiFab * > &  lsm_output_ptrs)
inline
224  {
225  if (m_update_rad) {
226  amrex::Print() << "Radiation advancing level " << m_lev << " at (YY-MM-DD SS) " << m_orbital_year << '-'
227  << m_orbital_mon << '-' << m_orbital_day << ' ' << m_orbital_sec << " ...";
228  this->initialize_impl();
229  this->run_impl();
230  this->finalize_impl(lsm_output_ptrs);
231  amrex::Print() << "DONE\n";
232  }
233  }
void initialize_impl()
Definition: ERF_Radiation.cpp:1195
void run_impl()
Definition: ERF_Radiation.cpp:1220
int m_orbital_mon
Definition: ERF_Radiation.H:403
void finalize_impl(amrex::Vector< amrex::MultiFab * > &lsm_output_ptrs)
Definition: ERF_Radiation.cpp:1520
int m_orbital_sec
Definition: ERF_Radiation.H:405
int m_lev
Definition: ERF_Radiation.H:275
bool m_update_rad
Definition: ERF_Radiation.H:293
int m_orbital_day
Definition: ERF_Radiation.H:404

Referenced by Run().

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

◆ Run()

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

Implements IRadiation.

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

◆ run_impl()

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

Referenced by Run().

Here is the caller graph for this function:

◆ write_rrtmgp_fluxes()

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

◆ WriteDataLog()

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

Implements IRadiation.

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

Member Data Documentation

◆ aero_g_sw

real3d_k Radiation::aero_g_sw
private

◆ aero_ssa_sw

real3d_k Radiation::aero_ssa_sw
private

◆ aero_tau_lw

real3d_k Radiation::aero_tau_lw
private

◆ aero_tau_sw

real3d_k Radiation::aero_tau_sw
private

◆ cldfrac_tot

real2d_k Radiation::cldfrac_tot
private

◆ d_tint

real2d_k Radiation::d_tint
private

◆ datalog_mf

amrex::MultiFab Radiation::datalog_mf
private

◆ eff_radius_qc

real2d_k Radiation::eff_radius_qc
private

◆ eff_radius_qi

real2d_k Radiation::eff_radius_qi
private

◆ gas_names_offset

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

◆ iwp

real2d_k Radiation::iwp
private

◆ lat

real1d_k Radiation::lat
private

◆ lon

real1d_k Radiation::lon
private

◆ lw_bnd_flux_dn

real3d_k Radiation::lw_bnd_flux_dn
private

◆ lw_bnd_flux_up

real3d_k Radiation::lw_bnd_flux_up
private

◆ lw_clnclrsky_flux_dn

real2d_k Radiation::lw_clnclrsky_flux_dn
private

◆ lw_clnclrsky_flux_up

real2d_k Radiation::lw_clnclrsky_flux_up
private

◆ lw_clnsky_flux_dn

real2d_k Radiation::lw_clnsky_flux_dn
private

◆ lw_clnsky_flux_up

real2d_k Radiation::lw_clnsky_flux_up
private

◆ lw_clrsky_flux_dn

real2d_k Radiation::lw_clrsky_flux_dn
private

◆ lw_clrsky_flux_up

real2d_k Radiation::lw_clrsky_flux_up
private

◆ lw_clrsky_heating

real2d_k Radiation::lw_clrsky_heating
private

◆ lw_flux_dn

real2d_k Radiation::lw_flux_dn
private

◆ lw_flux_up

real2d_k Radiation::lw_flux_up
private

◆ lw_heating

real2d_k Radiation::lw_heating
private

◆ lw_src

real1d_k Radiation::lw_src
private

◆ lwp

real2d_k Radiation::lwp
private

◆ m_ba

amrex::BoxArray Radiation::m_ba
private

◆ m_ch4vmr

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

◆ m_co2vmr

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

◆ m_col_offsets

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

Referenced by Init().

◆ m_cons_in

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

◆ m_covmr

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

◆ m_do_aerosol_rad

bool Radiation::m_do_aerosol_rad = false
private

◆ m_do_subcol_sampling

bool Radiation::m_do_subcol_sampling = true
private

◆ m_dt

double Radiation::m_dt
private

◆ m_extra_clnclrsky_diag

bool Radiation::m_extra_clnclrsky_diag = false
private

◆ m_extra_clnsky_diag

bool Radiation::m_extra_clnsky_diag = false
private

◆ m_fixed_orbital_year

bool Radiation::m_fixed_orbital_year = false
private

◆ m_fixed_solar_zenith_angle

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

◆ m_fixed_total_solar_irradiance

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

◆ m_gas_concs

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

◆ m_gas_mol_weights

real1d_k Radiation::m_gas_mol_weights
private

◆ m_gas_names

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

◆ m_geom

amrex::Geometry Radiation::m_geom
private

◆ m_ice

bool Radiation::m_ice = false
private

◆ m_is_nested_patch

bool Radiation::m_is_nested_patch = false
private

Referenced by Init(), is_nested_patch(), and Run().

◆ 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().

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