| 1 | ! radiation_interface.F90 - Public interface to radiation scheme |
|---|
| 2 | ! |
|---|
| 3 | ! (C) Copyright 2014- ECMWF. |
|---|
| 4 | ! |
|---|
| 5 | ! This software is licensed under the terms of the Apache Licence Version 2.0 |
|---|
| 6 | ! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0. |
|---|
| 7 | ! |
|---|
| 8 | ! In applying this licence, ECMWF does not waive the privileges and immunities |
|---|
| 9 | ! granted to it by virtue of its status as an intergovernmental organisation |
|---|
| 10 | ! nor does it submit to any jurisdiction. |
|---|
| 11 | ! |
|---|
| 12 | ! Author: Robin Hogan |
|---|
| 13 | ! Email: r.j.hogan@ecmwf.int |
|---|
| 14 | ! |
|---|
| 15 | ! Modifications |
|---|
| 16 | ! 2017-04-11 R. Hogan Changes to enable generalized surface description |
|---|
| 17 | ! 2017-09-08 R. Hogan Reverted some changes |
|---|
| 18 | ! |
|---|
| 19 | ! To use the radiation scheme, create a configuration_type object, |
|---|
| 20 | ! call "setup_radiation" on it once to load the look-up-tables and |
|---|
| 21 | ! data describing how gas and hydrometeor absorption/scattering are to |
|---|
| 22 | ! be represented, and call "radiation" multiple times on different |
|---|
| 23 | ! input profiles. |
|---|
| 24 | |
|---|
| 25 | module radiation_interface |
|---|
| 26 | |
|---|
| 27 | implicit none |
|---|
| 28 | |
|---|
| 29 | public :: setup_radiation, set_gas_units, radiation |
|---|
| 30 | private :: radiation_reverse |
|---|
| 31 | |
|---|
| 32 | contains |
|---|
| 33 | |
|---|
| 34 | !--------------------------------------------------------------------- |
|---|
| 35 | ! Load the look-up-tables and data describing how gas and |
|---|
| 36 | ! hydrometeor absorption/scattering are to be represented |
|---|
| 37 | subroutine setup_radiation(config) |
|---|
| 38 | |
|---|
| 39 | use parkind1, only : jprb |
|---|
| 40 | use yomhook, only : lhook, dr_hook, jphook |
|---|
| 41 | use radiation_io, only : nulerr, radiation_abort |
|---|
| 42 | use radiation_config, only : config_type, ISolverMcICA, & |
|---|
| 43 | & IGasModelMonochromatic, IGasModelIFSRRTMG, IGasModelECCKD |
|---|
| 44 | use radiation_spectral_definition, only & |
|---|
| 45 | & : SolarReferenceTemperature, TerrestrialReferenceTemperature |
|---|
| 46 | ! Currently there are two gas absorption models: RRTMG (default) |
|---|
| 47 | ! and monochromatic |
|---|
| 48 | use radiation_monochromatic, only : & |
|---|
| 49 | & setup_gas_optics_mono => setup_gas_optics, & |
|---|
| 50 | & setup_cloud_optics_mono => setup_cloud_optics, & |
|---|
| 51 | & setup_aerosol_optics_mono => setup_aerosol_optics |
|---|
| 52 | use radiation_ifs_rrtm, only : setup_gas_optics_rrtmg => setup_gas_optics |
|---|
| 53 | use radiation_ecckd_interface,only : setup_gas_optics_ecckd => setup_gas_optics |
|---|
| 54 | use radiation_cloud_optics, only : setup_cloud_optics |
|---|
| 55 | use radiation_general_cloud_optics, only : setup_general_cloud_optics |
|---|
| 56 | use radiation_aerosol_optics, only : setup_aerosol_optics |
|---|
| 57 | |
|---|
| 58 | |
|---|
| 59 | type(config_type), intent(inout) :: config |
|---|
| 60 | |
|---|
| 61 | real(jphook) :: hook_handle |
|---|
| 62 | |
|---|
| 63 | if (lhook) call dr_hook('radiation_interface:setup_radiation',0,hook_handle) |
|---|
| 64 | |
|---|
| 65 | ! Consolidate configuration data, including setting data file |
|---|
| 66 | ! names |
|---|
| 67 | call config%consolidate() |
|---|
| 68 | |
|---|
| 69 | ! Load the look-up tables from files in the specified directory |
|---|
| 70 | if (config%i_gas_model == IGasModelMonochromatic) then |
|---|
| 71 | call setup_gas_optics_mono(config, trim(config%directory_name)) |
|---|
| 72 | else if (config%i_gas_model == IGasModelIFSRRTMG) then |
|---|
| 73 | call setup_gas_optics_rrtmg(config, trim(config%directory_name)) |
|---|
| 74 | else if (config%i_gas_model == IGasModelECCKD) then |
|---|
| 75 | call setup_gas_optics_ecckd(config) |
|---|
| 76 | end if |
|---|
| 77 | |
|---|
| 78 | if (config%do_lw_aerosol_scattering & |
|---|
| 79 | & .and. .not. config%do_lw_cloud_scattering) then |
|---|
| 80 | write(nulerr, '(a)') '*** Error: longwave aerosol scattering requires longwave cloud scattering' |
|---|
| 81 | call radiation_abort('Radiation configuration error') |
|---|
| 82 | end if |
|---|
| 83 | |
|---|
| 84 | |
|---|
| 85 | ! Whether or not the "radiation" subroutine needs ssa_lw and g_lw |
|---|
| 86 | ! arrays depends on whether longwave scattering by aerosols is to |
|---|
| 87 | ! be included. If not, one of the array dimensions will be set to |
|---|
| 88 | ! zero. |
|---|
| 89 | if (config%do_lw_aerosol_scattering) then |
|---|
| 90 | config%n_g_lw_if_scattering = config%n_g_lw |
|---|
| 91 | else |
|---|
| 92 | config%n_g_lw_if_scattering = 0 |
|---|
| 93 | end if |
|---|
| 94 | |
|---|
| 95 | ! Whether or not the "radiation" subroutine needs ssa_lw_cloud and |
|---|
| 96 | ! g_lw_cloud arrays depends on whether longwave scattering by |
|---|
| 97 | ! hydrometeors is to be included. If not, one of the array |
|---|
| 98 | ! dimensions will be set to zero. |
|---|
| 99 | if (config%do_lw_cloud_scattering) then |
|---|
| 100 | config%n_bands_lw_if_scattering = config%n_bands_lw |
|---|
| 101 | else |
|---|
| 102 | config%n_bands_lw_if_scattering = 0 |
|---|
| 103 | end if |
|---|
| 104 | |
|---|
| 105 | ! If we have longwave scattering and McICA then even if there is |
|---|
| 106 | ! no aerosol, it is convenient if single scattering albedo and |
|---|
| 107 | ! g factor arrays are allocated before the call to |
|---|
| 108 | ! solver_lw as they will be needed. |
|---|
| 109 | if (config%do_lw_cloud_scattering & |
|---|
| 110 | & .and. config%i_solver_lw == ISolverMcICA) then |
|---|
| 111 | config%n_g_lw_if_scattering = config%n_g_lw |
|---|
| 112 | end if |
|---|
| 113 | |
|---|
| 114 | ! Consolidate the albedo/emissivity intervals with the shortwave |
|---|
| 115 | ! and longwave spectral bands |
|---|
| 116 | if (config%do_sw) then |
|---|
| 117 | call config%consolidate_sw_albedo_intervals |
|---|
| 118 | end if |
|---|
| 119 | if (config%do_lw) then |
|---|
| 120 | call config%consolidate_lw_emiss_intervals |
|---|
| 121 | end if |
|---|
| 122 | |
|---|
| 123 | if (config%do_clouds) then |
|---|
| 124 | if (config%i_gas_model == IGasModelMonochromatic) then |
|---|
| 125 | ! call setup_cloud_optics_mono(config) |
|---|
| 126 | elseif (config%use_general_cloud_optics) then |
|---|
| 127 | call setup_general_cloud_optics(config) |
|---|
| 128 | else |
|---|
| 129 | call setup_cloud_optics(config) |
|---|
| 130 | end if |
|---|
| 131 | end if |
|---|
| 132 | |
|---|
| 133 | if (config%use_aerosols) then |
|---|
| 134 | if (config%i_gas_model == IGasModelMonochromatic) then |
|---|
| 135 | ! call setup_aerosol_optics_mono(config) |
|---|
| 136 | else |
|---|
| 137 | call setup_aerosol_optics(config) |
|---|
| 138 | end if |
|---|
| 139 | end if |
|---|
| 140 | |
|---|
| 141 | ! Load cloud water PDF look-up table for McICA |
|---|
| 142 | if ( config%i_solver_sw == ISolverMcICA & |
|---|
| 143 | & .or. config%i_solver_lw == ISolverMcICA) then |
|---|
| 144 | call config%pdf_sampler%setup(config%cloud_pdf_file_name, & |
|---|
| 145 | & iverbose=config%iverbosesetup) |
|---|
| 146 | end if |
|---|
| 147 | |
|---|
| 148 | if (lhook) call dr_hook('radiation_interface:setup_radiation',1,hook_handle) |
|---|
| 149 | |
|---|
| 150 | end subroutine setup_radiation |
|---|
| 151 | |
|---|
| 152 | |
|---|
| 153 | !--------------------------------------------------------------------- |
|---|
| 154 | ! Scale the gas mixing ratios so that they have the units (and |
|---|
| 155 | ! possibly scale factors) required by the specific gas absorption |
|---|
| 156 | ! model. This subroutine simply passes the gas object on to the |
|---|
| 157 | ! module of the currently active gas model. |
|---|
| 158 | subroutine set_gas_units(config, gas) |
|---|
| 159 | |
|---|
| 160 | use radiation_config |
|---|
| 161 | use radiation_gas, only : gas_type |
|---|
| 162 | use radiation_monochromatic, only : set_gas_units_mono => set_gas_units |
|---|
| 163 | use radiation_ifs_rrtm, only : set_gas_units_ifs => set_gas_units |
|---|
| 164 | use radiation_ecckd_interface, only : set_gas_units_ecckd => set_gas_units |
|---|
| 165 | |
|---|
| 166 | type(config_type), intent(in) :: config |
|---|
| 167 | type(gas_type), intent(inout) :: gas |
|---|
| 168 | |
|---|
| 169 | if (config%i_gas_model == IGasModelMonochromatic) then |
|---|
| 170 | call set_gas_units_mono(gas) |
|---|
| 171 | elseif (config%i_gas_model == IGasModelECCKD) then |
|---|
| 172 | call set_gas_units_ecckd(gas) |
|---|
| 173 | else |
|---|
| 174 | call set_gas_units_ifs(gas) |
|---|
| 175 | end if |
|---|
| 176 | |
|---|
| 177 | end subroutine set_gas_units |
|---|
| 178 | |
|---|
| 179 | |
|---|
| 180 | !--------------------------------------------------------------------- |
|---|
| 181 | ! Run the radiation scheme according to the configuration in the |
|---|
| 182 | ! config object. There are ncol profiles of which only istartcol to |
|---|
| 183 | ! iendcol are to be processed, and there are nlev model levels. The |
|---|
| 184 | ! output fluxes are written to the flux object, and all other |
|---|
| 185 | ! objects contain the input variables. The variables may be defined |
|---|
| 186 | ! either in order of increasing or decreasing pressure, but if in |
|---|
| 187 | ! order of decreasing pressure then radiation_reverse will be called |
|---|
| 188 | ! to reverse the order for the computation and then reverse the |
|---|
| 189 | ! order of the output fluxes to match the inputs. |
|---|
| 190 | subroutine radiation(ncol, nlev, istartcol, iendcol, config, & |
|---|
| 191 | & single_level, thermodynamics, gas, cloud, aerosol, flux) |
|---|
| 192 | |
|---|
| 193 | use parkind1, only : jprb |
|---|
| 194 | use yomhook, only : lhook, dr_hook, jphook |
|---|
| 195 | |
|---|
| 196 | use radiation_io, only : nulout |
|---|
| 197 | use radiation_config, only : config_type, & |
|---|
| 198 | & IGasModelMonochromatic, IGasModelIFSRRTMG, & |
|---|
| 199 | & ISolverMcICA, ISolverSpartacus, ISolverHomogeneous, & |
|---|
| 200 | & ISolverTripleclouds |
|---|
| 201 | use radiation_single_level, only : single_level_type |
|---|
| 202 | use radiation_thermodynamics, only : thermodynamics_type |
|---|
| 203 | use radiation_gas, only : gas_type |
|---|
| 204 | use radiation_cloud, only : cloud_type |
|---|
| 205 | use radiation_aerosol, only : aerosol_type |
|---|
| 206 | use radiation_flux, only : flux_type |
|---|
| 207 | use radiation_spartacus_sw, only : solver_spartacus_sw |
|---|
| 208 | use radiation_spartacus_lw, only : solver_spartacus_lw |
|---|
| 209 | use radiation_tripleclouds_sw,only : solver_tripleclouds_sw |
|---|
| 210 | use radiation_tripleclouds_lw,only : solver_tripleclouds_lw |
|---|
| 211 | use radiation_mcica_sw, only : solver_mcica_sw |
|---|
| 212 | use radiation_mcica_lw, only : solver_mcica_lw |
|---|
| 213 | use radiation_cloudless_sw, only : solver_cloudless_sw |
|---|
| 214 | use radiation_cloudless_lw, only : solver_cloudless_lw |
|---|
| 215 | use radiation_homogeneous_sw, only : solver_homogeneous_sw |
|---|
| 216 | use radiation_homogeneous_lw, only : solver_homogeneous_lw |
|---|
| 217 | use radiation_save, only : save_radiative_properties |
|---|
| 218 | |
|---|
| 219 | ! Treatment of gas and hydrometeor optics |
|---|
| 220 | use radiation_monochromatic, only : & |
|---|
| 221 | & gas_optics_mono => gas_optics, & |
|---|
| 222 | & cloud_optics_mono => cloud_optics, & |
|---|
| 223 | & add_aerosol_optics_mono => add_aerosol_optics |
|---|
| 224 | use radiation_ifs_rrtm, only : gas_optics_rrtmg => gas_optics |
|---|
| 225 | use radiation_ecckd_interface,only : gas_optics_ecckd => gas_optics |
|---|
| 226 | use radiation_cloud_optics, only : cloud_optics |
|---|
| 227 | use radiation_general_cloud_optics, only : general_cloud_optics |
|---|
| 228 | use radiation_aerosol_optics, only : add_aerosol_optics |
|---|
| 229 | USE mod_phys_lmdz_para |
|---|
| 230 | |
|---|
| 231 | |
|---|
| 232 | ! Inputs |
|---|
| 233 | integer, intent(in) :: ncol ! number of columns |
|---|
| 234 | integer, intent(in) :: nlev ! number of model levels |
|---|
| 235 | integer, intent(in) :: istartcol, iendcol ! range of columns to process |
|---|
| 236 | type(config_type), intent(in) :: config |
|---|
| 237 | type(single_level_type), intent(in) :: single_level |
|---|
| 238 | type(thermodynamics_type),intent(in) :: thermodynamics |
|---|
| 239 | type(gas_type), intent(in) :: gas |
|---|
| 240 | type(cloud_type), intent(inout):: cloud |
|---|
| 241 | type(aerosol_type), intent(in) :: aerosol |
|---|
| 242 | ! Output |
|---|
| 243 | type(flux_type), intent(inout):: flux |
|---|
| 244 | |
|---|
| 245 | |
|---|
| 246 | ! Local variables |
|---|
| 247 | |
|---|
| 248 | ! Layer optical depth, single scattering albedo and asymmetry factor of |
|---|
| 249 | ! gases and aerosols at each longwave g-point, where the latter |
|---|
| 250 | ! two variables are only defined if aerosol longwave scattering is |
|---|
| 251 | ! enabled (otherwise both are treated as zero). |
|---|
| 252 | real(jprb), dimension(config%n_g_lw,nlev,istartcol:iendcol) :: od_lw |
|---|
| 253 | real(jprb), dimension(config%n_g_lw_if_scattering,nlev,istartcol:iendcol) :: & |
|---|
| 254 | & ssa_lw, g_lw |
|---|
| 255 | |
|---|
| 256 | ! Layer in-cloud optical depth, single scattering albedo and |
|---|
| 257 | ! asymmetry factor of hydrometeors in each longwave band, where |
|---|
| 258 | ! the latter two variables are only defined if hydrometeor |
|---|
| 259 | ! longwave scattering is enabled (otherwise both are treated as |
|---|
| 260 | ! zero). |
|---|
| 261 | real(jprb), dimension(config%n_bands_lw,nlev,istartcol:iendcol) :: od_lw_cloud |
|---|
| 262 | real(jprb), dimension(config%n_bands_lw_if_scattering,nlev,istartcol:iendcol) :: & |
|---|
| 263 | & ssa_lw_cloud, g_lw_cloud |
|---|
| 264 | |
|---|
| 265 | ! Layer optical depth, single scattering albedo and asymmetry factor of |
|---|
| 266 | ! gases and aerosols at each shortwave g-point |
|---|
| 267 | real(jprb), dimension(config%n_g_sw,nlev,istartcol:iendcol) :: od_sw, ssa_sw, g_sw |
|---|
| 268 | |
|---|
| 269 | ! Layer in-cloud optical depth, single scattering albedo and |
|---|
| 270 | ! asymmetry factor of hydrometeors in each shortwave band |
|---|
| 271 | real(jprb), dimension(config%n_bands_sw,nlev,istartcol:iendcol) :: & |
|---|
| 272 | & od_sw_cloud, ssa_sw_cloud, g_sw_cloud |
|---|
| 273 | |
|---|
| 274 | ! The Planck function (emitted flux from a black body) at half |
|---|
| 275 | ! levels |
|---|
| 276 | real(jprb), dimension(config%n_g_lw,nlev+1,istartcol:iendcol) :: planck_hl |
|---|
| 277 | |
|---|
| 278 | ! The longwave emission from and albedo of the surface in each |
|---|
| 279 | ! longwave g-point; note that these are weighted averages of the |
|---|
| 280 | ! values from individual tiles |
|---|
| 281 | real(jprb), dimension(config%n_g_lw, istartcol:iendcol) :: lw_emission |
|---|
| 282 | real(jprb), dimension(config%n_g_lw, istartcol:iendcol) :: lw_albedo |
|---|
| 283 | |
|---|
| 284 | ! Direct and diffuse shortwave surface albedo in each shortwave |
|---|
| 285 | ! g-point; note that these are weighted averages of the values |
|---|
| 286 | ! from individual tiles |
|---|
| 287 | real(jprb), dimension(config%n_g_sw, istartcol:iendcol) :: sw_albedo_direct |
|---|
| 288 | real(jprb), dimension(config%n_g_sw, istartcol:iendcol) :: sw_albedo_diffuse |
|---|
| 289 | |
|---|
| 290 | ! The incoming shortwave flux into a plane perpendicular to the |
|---|
| 291 | ! incoming radiation at top-of-atmosphere in each of the shortwave |
|---|
| 292 | ! g-points |
|---|
| 293 | real(jprb), dimension(config%n_g_sw,istartcol:iendcol) :: incoming_sw |
|---|
| 294 | |
|---|
| 295 | character(len=100) :: rad_prop_file_name |
|---|
| 296 | character(*), parameter :: rad_prop_base_file_name = "radiative_properties" |
|---|
| 297 | |
|---|
| 298 | real(jphook) :: hook_handle |
|---|
| 299 | integer :: jcol, jlev |
|---|
| 300 | |
|---|
| 301 | if (lhook) call dr_hook('radiation_interface:radiation',0,hook_handle) |
|---|
| 302 | |
|---|
| 303 | if (config%i_solver_sw == ISolverSPARTACUS) then |
|---|
| 304 | print*,'Dans radiation, mpi_rank, omp_rank, size, chape inv_cloud = ',& |
|---|
| 305 | mpi_rank, omp_rank, & |
|---|
| 306 | shape(cloud%inv_cloud_effective_size), & |
|---|
| 307 | size(cloud%inv_cloud_effective_size) |
|---|
| 308 | ! do jcol=istartcol, iendcol |
|---|
| 309 | ! do jlev=1,nlev |
|---|
| 310 | ! print*,'Entree radiation_interf, mpi_rank, omp_rank, jcol, jlev & |
|---|
| 311 | ! & cloud%inv_cloud_effective_size =',mpi_rank, omp_rank, jcol, jlev, & |
|---|
| 312 | ! & cloud%inv_cloud_effective_size(jcol,jlev) |
|---|
| 313 | ! enddo |
|---|
| 314 | ! enddo |
|---|
| 315 | endif |
|---|
| 316 | ! cloud%inv_cloud_effective_size=0.05_jprb |
|---|
| 317 | |
|---|
| 318 | if (thermodynamics%pressure_hl(istartcol,2) & |
|---|
| 319 | & < thermodynamics%pressure_hl(istartcol,1)) then |
|---|
| 320 | ! Input arrays are arranged in order of decreasing pressure / |
|---|
| 321 | ! increasing height: the following subroutine reverses them, |
|---|
| 322 | ! call the radiation scheme and then reverses the returned |
|---|
| 323 | ! fluxes |
|---|
| 324 | call radiation_reverse(ncol, nlev, istartcol, iendcol, config, & |
|---|
| 325 | & single_level, thermodynamics, gas, cloud, aerosol, flux) |
|---|
| 326 | else |
|---|
| 327 | |
|---|
| 328 | ! Input arrays arranged in order of increasing pressure / |
|---|
| 329 | ! decreasing height: progress normally |
|---|
| 330 | |
|---|
| 331 | ! Extract surface albedos at each gridpoint |
|---|
| 332 | call single_level%get_albedos(istartcol, iendcol, config, & |
|---|
| 333 | & sw_albedo_direct, sw_albedo_diffuse, & |
|---|
| 334 | & lw_albedo) |
|---|
| 335 | |
|---|
| 336 | ! Compute gas absorption optical depth in shortwave and |
|---|
| 337 | ! longwave, shortwave single scattering albedo (i.e. fraction of |
|---|
| 338 | ! extinction due to Rayleigh scattering), Planck functions and |
|---|
| 339 | ! incoming shortwave flux at each g-point, for the specified |
|---|
| 340 | ! range of atmospheric columns |
|---|
| 341 | if (config%i_gas_model == IGasModelMonochromatic) then |
|---|
| 342 | call gas_optics_mono(ncol,nlev,istartcol,iendcol, config, & |
|---|
| 343 | & single_level, thermodynamics, gas, lw_albedo, & |
|---|
| 344 | & od_lw, od_sw, ssa_sw, & |
|---|
| 345 | & planck_hl, lw_emission, incoming_sw) |
|---|
| 346 | else if (config%i_gas_model == IGasModelIFSRRTMG) then |
|---|
| 347 | call gas_optics_rrtmg(ncol,nlev,istartcol,iendcol, config, & |
|---|
| 348 | & single_level, thermodynamics, gas, & |
|---|
| 349 | & od_lw, od_sw, ssa_sw, lw_albedo=lw_albedo, & |
|---|
| 350 | & planck_hl=planck_hl, lw_emission=lw_emission, & |
|---|
| 351 | & incoming_sw=incoming_sw) |
|---|
| 352 | else |
|---|
| 353 | call gas_optics_ecckd(ncol,nlev,istartcol,iendcol, config, & |
|---|
| 354 | & single_level, thermodynamics, gas, & |
|---|
| 355 | & od_lw, od_sw, ssa_sw, lw_albedo=lw_albedo, & |
|---|
| 356 | & planck_hl=planck_hl, lw_emission=lw_emission, & |
|---|
| 357 | & incoming_sw=incoming_sw) |
|---|
| 358 | end if |
|---|
| 359 | |
|---|
| 360 | if (config%do_clouds) then |
|---|
| 361 | ! Crop the cloud fraction to remove clouds that have too small |
|---|
| 362 | ! a fraction or water content; after this, we can safely |
|---|
| 363 | ! assume that a cloud is present if cloud%fraction > 0.0. |
|---|
| 364 | call cloud%crop_cloud_fraction(istartcol, iendcol, & |
|---|
| 365 | & config%cloud_fraction_threshold, & |
|---|
| 366 | & config%cloud_mixing_ratio_threshold) |
|---|
| 367 | |
|---|
| 368 | ! Compute hydrometeor absorption/scattering properties in each |
|---|
| 369 | ! shortwave and longwave band |
|---|
| 370 | if (config%i_gas_model == IGasModelMonochromatic) then |
|---|
| 371 | call cloud_optics_mono(nlev, istartcol, iendcol, & |
|---|
| 372 | & config, thermodynamics, cloud, & |
|---|
| 373 | & od_lw_cloud, ssa_lw_cloud, g_lw_cloud, & |
|---|
| 374 | & od_sw_cloud, ssa_sw_cloud, g_sw_cloud) |
|---|
| 375 | elseif (config%use_general_cloud_optics) then |
|---|
| 376 | call general_cloud_optics(nlev, istartcol, iendcol, & |
|---|
| 377 | & config, thermodynamics, cloud, & |
|---|
| 378 | & od_lw_cloud, ssa_lw_cloud, g_lw_cloud, & |
|---|
| 379 | & od_sw_cloud, ssa_sw_cloud, g_sw_cloud) |
|---|
| 380 | else |
|---|
| 381 | call cloud_optics(nlev, istartcol, iendcol, & |
|---|
| 382 | & config, thermodynamics, cloud, & |
|---|
| 383 | & od_lw_cloud, ssa_lw_cloud, g_lw_cloud, & |
|---|
| 384 | & od_sw_cloud, ssa_sw_cloud, g_sw_cloud) |
|---|
| 385 | end if |
|---|
| 386 | end if ! do_clouds |
|---|
| 387 | |
|---|
| 388 | if (config%use_aerosols) then |
|---|
| 389 | if (config%i_gas_model == IGasModelMonochromatic) then |
|---|
| 390 | ! call add_aerosol_optics_mono(nlev,istartcol,iendcol, & |
|---|
| 391 | ! & config, thermodynamics, gas, aerosol, & |
|---|
| 392 | ! & od_lw, ssa_lw, g_lw, od_sw, ssa_sw, g_sw) |
|---|
| 393 | else |
|---|
| 394 | call add_aerosol_optics(nlev,istartcol,iendcol, & |
|---|
| 395 | & config, thermodynamics, gas, aerosol, & |
|---|
| 396 | & od_lw, ssa_lw, g_lw, od_sw, ssa_sw, g_sw) |
|---|
| 397 | end if |
|---|
| 398 | else |
|---|
| 399 | g_sw(:,:,istartcol:iendcol) = 0.0_jprb |
|---|
| 400 | if (config%do_lw_aerosol_scattering) then |
|---|
| 401 | ssa_lw(:,:,istartcol:iendcol) = 0.0_jprb |
|---|
| 402 | g_lw(:,:,istartcol:iendcol) = 0.0_jprb |
|---|
| 403 | end if |
|---|
| 404 | end if |
|---|
| 405 | |
|---|
| 406 | ! For diagnostic purposes, save these intermediate variables to |
|---|
| 407 | ! a NetCDF file |
|---|
| 408 | if (config%do_save_radiative_properties) then |
|---|
| 409 | if (istartcol == 1 .and. iendcol == ncol) then |
|---|
| 410 | rad_prop_file_name = rad_prop_base_file_name // ".nc" |
|---|
| 411 | else |
|---|
| 412 | write(rad_prop_file_name,'(a,a,i4.4,a,i4.4,a)') & |
|---|
| 413 | & rad_prop_base_file_name, '_', istartcol, '-',iendcol,'.nc' |
|---|
| 414 | end if |
|---|
| 415 | call save_radiative_properties(trim(rad_prop_file_name), & |
|---|
| 416 | & nlev, istartcol, iendcol, & |
|---|
| 417 | & config, single_level, thermodynamics, cloud, & |
|---|
| 418 | & planck_hl, lw_emission, lw_albedo, & |
|---|
| 419 | & sw_albedo_direct, sw_albedo_diffuse, incoming_sw, & |
|---|
| 420 | & od_lw, ssa_lw, g_lw, od_sw, ssa_sw, g_sw, & |
|---|
| 421 | & od_lw_cloud, ssa_lw_cloud, g_lw_cloud, & |
|---|
| 422 | & od_sw_cloud, ssa_sw_cloud, g_sw_cloud) |
|---|
| 423 | end if |
|---|
| 424 | |
|---|
| 425 | if (config%do_lw) then |
|---|
| 426 | if (config%iverbose >= 2) then |
|---|
| 427 | write(nulout,'(a)') 'Computing longwave fluxes' |
|---|
| 428 | end if |
|---|
| 429 | |
|---|
| 430 | if (config%i_solver_lw == ISolverMcICA) then |
|---|
| 431 | ! Compute fluxes using the McICA longwave solver |
|---|
| 432 | call solver_mcica_lw(nlev,istartcol,iendcol, & |
|---|
| 433 | & config, single_level, cloud, & |
|---|
| 434 | & od_lw, ssa_lw, g_lw, od_lw_cloud, ssa_lw_cloud, & |
|---|
| 435 | & g_lw_cloud, planck_hl, lw_emission, lw_albedo, flux) |
|---|
| 436 | else if (config%i_solver_lw == ISolverSPARTACUS) then |
|---|
| 437 | ! Compute fluxes using the SPARTACUS longwave solver |
|---|
| 438 | call solver_spartacus_lw(nlev,istartcol,iendcol, & |
|---|
| 439 | & config, thermodynamics, cloud, & |
|---|
| 440 | & od_lw, ssa_lw, g_lw, od_lw_cloud, ssa_lw_cloud, g_lw_cloud, & |
|---|
| 441 | & planck_hl, lw_emission, lw_albedo, flux) |
|---|
| 442 | else if (config%i_solver_lw == ISolverTripleclouds) then |
|---|
| 443 | ! Compute fluxes using the Tripleclouds longwave solver |
|---|
| 444 | call solver_tripleclouds_lw(nlev,istartcol,iendcol, & |
|---|
| 445 | & config, cloud, & |
|---|
| 446 | & od_lw, ssa_lw, g_lw, od_lw_cloud, ssa_lw_cloud, g_lw_cloud, & |
|---|
| 447 | & planck_hl, lw_emission, lw_albedo, flux) |
|---|
| 448 | elseif (config%i_solver_lw == ISolverHomogeneous) then |
|---|
| 449 | ! Compute fluxes using the homogeneous solver |
|---|
| 450 | call solver_homogeneous_lw(nlev,istartcol,iendcol, & |
|---|
| 451 | & config, cloud, & |
|---|
| 452 | & od_lw, ssa_lw, g_lw, od_lw_cloud, ssa_lw_cloud, & |
|---|
| 453 | & g_lw_cloud, planck_hl, lw_emission, lw_albedo, flux) |
|---|
| 454 | else |
|---|
| 455 | ! Compute fluxes using the cloudless solver |
|---|
| 456 | call solver_cloudless_lw(nlev,istartcol,iendcol, & |
|---|
| 457 | & config, od_lw, ssa_lw, g_lw, & |
|---|
| 458 | & planck_hl, lw_emission, lw_albedo, flux) |
|---|
| 459 | end if |
|---|
| 460 | end if |
|---|
| 461 | |
|---|
| 462 | if (config%do_sw) then |
|---|
| 463 | if (config%iverbose >= 2) then |
|---|
| 464 | write(nulout,'(a)') 'Computing shortwave fluxes' |
|---|
| 465 | end if |
|---|
| 466 | |
|---|
| 467 | if (config%i_solver_sw == ISolverMcICA) then |
|---|
| 468 | ! Compute fluxes using the McICA shortwave solver |
|---|
| 469 | call solver_mcica_sw(nlev,istartcol,iendcol, & |
|---|
| 470 | & config, single_level, cloud, & |
|---|
| 471 | & od_sw, ssa_sw, g_sw, od_sw_cloud, ssa_sw_cloud, & |
|---|
| 472 | & g_sw_cloud, sw_albedo_direct, sw_albedo_diffuse, & |
|---|
| 473 | & incoming_sw, flux) |
|---|
| 474 | else if (config%i_solver_sw == ISolverSPARTACUS) then |
|---|
| 475 | ! Compute fluxes using the SPARTACUS shortwave solver |
|---|
| 476 | ! cloud%inv_cloud_effective_size=0.05_jprb |
|---|
| 477 | ! do jcol=istartcol, iendcol |
|---|
| 478 | ! do jlev=1,nlev |
|---|
| 479 | ! print*,'jcol, jlev, dans radiation_interf i & |
|---|
| 480 | ! & cloud%inv_cloud_effective_size =',jcol, jlev, & |
|---|
| 481 | ! cloud%inv_cloud_effective_size(jcol,jlev) |
|---|
| 482 | ! enddo |
|---|
| 483 | ! enddo |
|---|
| 484 | call solver_spartacus_sw(nlev,istartcol,iendcol, & |
|---|
| 485 | & config, single_level, thermodynamics, cloud, & |
|---|
| 486 | & od_sw, ssa_sw, g_sw, od_sw_cloud, ssa_sw_cloud, & |
|---|
| 487 | & g_sw_cloud, sw_albedo_direct, sw_albedo_diffuse, & |
|---|
| 488 | & incoming_sw, flux) |
|---|
| 489 | else if (config%i_solver_sw == ISolverTripleclouds) then |
|---|
| 490 | ! Compute fluxes using the Tripleclouds shortwave solver |
|---|
| 491 | call solver_tripleclouds_sw(nlev,istartcol,iendcol, & |
|---|
| 492 | & config, single_level, cloud, & |
|---|
| 493 | & od_sw, ssa_sw, g_sw, od_sw_cloud, ssa_sw_cloud, & |
|---|
| 494 | & g_sw_cloud, sw_albedo_direct, sw_albedo_diffuse, & |
|---|
| 495 | & incoming_sw, flux) |
|---|
| 496 | elseif (config%i_solver_sw == ISolverHomogeneous) then |
|---|
| 497 | ! Compute fluxes using the homogeneous solver |
|---|
| 498 | call solver_homogeneous_sw(nlev,istartcol,iendcol, & |
|---|
| 499 | & config, single_level, cloud, & |
|---|
| 500 | & od_sw, ssa_sw, g_sw, od_sw_cloud, ssa_sw_cloud, & |
|---|
| 501 | & g_sw_cloud, sw_albedo_direct, sw_albedo_diffuse, & |
|---|
| 502 | & incoming_sw, flux) |
|---|
| 503 | else |
|---|
| 504 | ! Compute fluxes using the cloudless solver |
|---|
| 505 | call solver_cloudless_sw(nlev,istartcol,iendcol, & |
|---|
| 506 | & config, single_level, od_sw, ssa_sw, g_sw, & |
|---|
| 507 | & sw_albedo_direct, sw_albedo_diffuse, & |
|---|
| 508 | & incoming_sw, flux) |
|---|
| 509 | end if |
|---|
| 510 | end if |
|---|
| 511 | |
|---|
| 512 | ! Store surface downwelling, and TOA, fluxes in bands from |
|---|
| 513 | ! fluxes in g points |
|---|
| 514 | call flux%calc_surface_spectral(config, istartcol, iendcol) |
|---|
| 515 | call flux%calc_toa_spectral (config, istartcol, iendcol) |
|---|
| 516 | |
|---|
| 517 | end if |
|---|
| 518 | |
|---|
| 519 | if (lhook) call dr_hook('radiation_interface:radiation',1,hook_handle) |
|---|
| 520 | |
|---|
| 521 | end subroutine radiation |
|---|
| 522 | |
|---|
| 523 | |
|---|
| 524 | !--------------------------------------------------------------------- |
|---|
| 525 | ! If the input arrays are arranged in order of decreasing pressure / |
|---|
| 526 | ! increasing height then this subroutine reverses them, calls the |
|---|
| 527 | ! radiation scheme and then reverses the returned fluxes. Since this |
|---|
| 528 | ! subroutine calls, and is called by "radiation", it must be in this |
|---|
| 529 | ! module to avoid circular dependencies. |
|---|
| 530 | subroutine radiation_reverse(ncol, nlev, istartcol, iendcol, config, & |
|---|
| 531 | & single_level, thermodynamics, gas, cloud, aerosol, flux) |
|---|
| 532 | |
|---|
| 533 | use parkind1, only : jprb |
|---|
| 534 | |
|---|
| 535 | use radiation_io, only : nulout |
|---|
| 536 | use radiation_config, only : config_type |
|---|
| 537 | use radiation_single_level, only : single_level_type |
|---|
| 538 | use radiation_thermodynamics, only : thermodynamics_type |
|---|
| 539 | use radiation_gas, only : gas_type |
|---|
| 540 | use radiation_cloud, only : cloud_type |
|---|
| 541 | use radiation_aerosol, only : aerosol_type |
|---|
| 542 | use radiation_flux, only : flux_type |
|---|
| 543 | |
|---|
| 544 | ! Inputs |
|---|
| 545 | integer, intent(in) :: ncol ! number of columns |
|---|
| 546 | integer, intent(in) :: nlev ! number of model levels |
|---|
| 547 | integer, intent(in) :: istartcol, iendcol ! range of columns to process |
|---|
| 548 | type(config_type), intent(in) :: config |
|---|
| 549 | type(single_level_type), intent(in) :: single_level |
|---|
| 550 | type(thermodynamics_type),intent(in) :: thermodynamics |
|---|
| 551 | type(gas_type), intent(in) :: gas |
|---|
| 552 | type(cloud_type), intent(in) :: cloud |
|---|
| 553 | type(aerosol_type), intent(in) :: aerosol |
|---|
| 554 | ! Output |
|---|
| 555 | type(flux_type), intent(inout):: flux |
|---|
| 556 | |
|---|
| 557 | ! Reversed data structures |
|---|
| 558 | type(thermodynamics_type) :: thermodynamics_rev |
|---|
| 559 | type(gas_type) :: gas_rev |
|---|
| 560 | type(cloud_type) :: cloud_rev |
|---|
| 561 | type(aerosol_type) :: aerosol_rev |
|---|
| 562 | type(flux_type) :: flux_rev |
|---|
| 563 | |
|---|
| 564 | ! Start and end levels for aerosol data |
|---|
| 565 | integer :: istartlev, iendlev |
|---|
| 566 | |
|---|
| 567 | if (config%iverbose >= 2) then |
|---|
| 568 | write(nulout,'(a)') 'Reversing arrays to be in order of increasing pressure' |
|---|
| 569 | end if |
|---|
| 570 | |
|---|
| 571 | ! Allocate reversed arrays |
|---|
| 572 | call thermodynamics_rev%allocate(ncol, nlev) |
|---|
| 573 | call cloud_rev%allocate(ncol, nlev) |
|---|
| 574 | call flux_rev%allocate(config, istartcol, iendcol, nlev) |
|---|
| 575 | if (allocated(aerosol%mixing_ratio)) then |
|---|
| 576 | istartlev = nlev + 1 - aerosol%iendlev |
|---|
| 577 | iendlev = nlev + 1 - aerosol%istartlev |
|---|
| 578 | call aerosol_rev%allocate(ncol, istartlev, iendlev, & |
|---|
| 579 | & config%n_aerosol_types) |
|---|
| 580 | end if |
|---|
| 581 | |
|---|
| 582 | ! Fill reversed thermodynamic arrays |
|---|
| 583 | thermodynamics_rev%pressure_hl(istartcol:iendcol,:) & |
|---|
| 584 | & = thermodynamics%pressure_hl(istartcol:iendcol, nlev+1:1:-1) |
|---|
| 585 | thermodynamics_rev%temperature_hl(istartcol:iendcol,:) & |
|---|
| 586 | & = thermodynamics%temperature_hl(istartcol:iendcol, nlev+1:1:-1) |
|---|
| 587 | |
|---|
| 588 | ! Fill reversed gas arrays |
|---|
| 589 | call gas%reverse(istartcol, iendcol, gas_rev) |
|---|
| 590 | |
|---|
| 591 | if (config%do_clouds) then |
|---|
| 592 | ! Fill reversed cloud arrays |
|---|
| 593 | cloud_rev%q_liq(istartcol:iendcol,:) & |
|---|
| 594 | & = cloud%q_liq(istartcol:iendcol,nlev:1:-1) |
|---|
| 595 | cloud_rev%re_liq(istartcol:iendcol,:) & |
|---|
| 596 | & = cloud%re_liq(istartcol:iendcol,nlev:1:-1) |
|---|
| 597 | cloud_rev%q_ice(istartcol:iendcol,:) & |
|---|
| 598 | & = cloud%q_ice(istartcol:iendcol,nlev:1:-1) |
|---|
| 599 | cloud_rev%re_ice(istartcol:iendcol,:) & |
|---|
| 600 | & = cloud%re_ice(istartcol:iendcol,nlev:1:-1) |
|---|
| 601 | cloud_rev%fraction(istartcol:iendcol,:) & |
|---|
| 602 | & = cloud%fraction(istartcol:iendcol,nlev:1:-1) |
|---|
| 603 | cloud_rev%overlap_param(istartcol:iendcol,:) & |
|---|
| 604 | & = cloud%overlap_param(istartcol:iendcol,nlev-1:1:-1) |
|---|
| 605 | if (allocated(cloud%fractional_std)) then |
|---|
| 606 | cloud_rev%fractional_std(istartcol:iendcol,:) & |
|---|
| 607 | & = cloud%fractional_std(istartcol:iendcol,nlev:1:-1) |
|---|
| 608 | else |
|---|
| 609 | cloud_rev%fractional_std(istartcol:iendcol,:) = 0.0_jprb |
|---|
| 610 | end if |
|---|
| 611 | if (allocated(cloud%inv_cloud_effective_size)) then |
|---|
| 612 | cloud_rev%inv_cloud_effective_size(istartcol:iendcol,:) & |
|---|
| 613 | & = cloud%inv_cloud_effective_size(istartcol:iendcol,nlev:1:-1) |
|---|
| 614 | else |
|---|
| 615 | cloud_rev%inv_cloud_effective_size(istartcol:iendcol,:) = 0.0_jprb |
|---|
| 616 | end if |
|---|
| 617 | end if |
|---|
| 618 | |
|---|
| 619 | if (allocated(aerosol%mixing_ratio)) then |
|---|
| 620 | aerosol_rev%mixing_ratio(:,istartlev:iendlev,:) & |
|---|
| 621 | & = aerosol%mixing_ratio(:,aerosol%iendlev:aerosol%istartlev:-1,:) |
|---|
| 622 | end if |
|---|
| 623 | |
|---|
| 624 | ! Run radiation scheme on reversed profiles |
|---|
| 625 | call radiation(ncol, nlev,istartcol,iendcol, & |
|---|
| 626 | & config, single_level, thermodynamics_rev, gas_rev, & |
|---|
| 627 | & cloud_rev, aerosol_rev, flux_rev) |
|---|
| 628 | |
|---|
| 629 | ! Reorder fluxes |
|---|
| 630 | if (allocated(flux%lw_up)) then |
|---|
| 631 | flux%lw_up(istartcol:iendcol,:) & |
|---|
| 632 | & = flux_rev%lw_up(istartcol:iendcol,nlev+1:1:-1) |
|---|
| 633 | flux%lw_dn(istartcol:iendcol,:) & |
|---|
| 634 | & = flux_rev%lw_dn(istartcol:iendcol,nlev+1:1:-1) |
|---|
| 635 | if (allocated(flux%lw_up_clear)) then |
|---|
| 636 | flux%lw_up_clear(istartcol:iendcol,:) & |
|---|
| 637 | & = flux_rev%lw_up_clear(istartcol:iendcol,nlev+1:1:-1) |
|---|
| 638 | flux%lw_dn_clear(istartcol:iendcol,:) & |
|---|
| 639 | & = flux_rev%lw_dn_clear(istartcol:iendcol,nlev+1:1:-1) |
|---|
| 640 | end if |
|---|
| 641 | end if |
|---|
| 642 | if (allocated(flux%sw_up)) then |
|---|
| 643 | flux%sw_up(istartcol:iendcol,:) & |
|---|
| 644 | & = flux_rev%sw_up(istartcol:iendcol,nlev+1:1:-1) |
|---|
| 645 | flux%sw_dn(istartcol:iendcol,:) & |
|---|
| 646 | & = flux_rev%sw_dn(istartcol:iendcol,nlev+1:1:-1) |
|---|
| 647 | if (allocated(flux%sw_dn_direct)) then |
|---|
| 648 | flux%sw_dn_direct(istartcol:iendcol,:) & |
|---|
| 649 | & = flux_rev%sw_dn_direct(istartcol:iendcol,nlev+1:1:-1) |
|---|
| 650 | end if |
|---|
| 651 | if (allocated(flux%sw_up_clear)) then |
|---|
| 652 | flux%sw_up_clear(istartcol:iendcol,:) & |
|---|
| 653 | & = flux_rev%sw_up_clear(istartcol:iendcol,nlev+1:1:-1) |
|---|
| 654 | flux%sw_dn_clear(istartcol:iendcol,:) & |
|---|
| 655 | & = flux_rev%sw_dn_clear(istartcol:iendcol,nlev+1:1:-1) |
|---|
| 656 | if (allocated(flux%sw_dn_direct_clear)) then |
|---|
| 657 | flux%sw_dn_direct_clear(istartcol:iendcol,:) & |
|---|
| 658 | & = flux_rev%sw_dn_direct_clear(istartcol:iendcol,nlev+1:1:-1) |
|---|
| 659 | end if |
|---|
| 660 | end if |
|---|
| 661 | end if |
|---|
| 662 | |
|---|
| 663 | ! Deallocate reversed arrays |
|---|
| 664 | call thermodynamics_rev%deallocate |
|---|
| 665 | call gas_rev%deallocate |
|---|
| 666 | call cloud_rev%deallocate |
|---|
| 667 | call flux_rev%deallocate |
|---|
| 668 | if (allocated(aerosol%mixing_ratio)) then |
|---|
| 669 | call aerosol_rev%deallocate |
|---|
| 670 | end if |
|---|
| 671 | |
|---|
| 672 | end subroutine radiation_reverse |
|---|
| 673 | |
|---|
| 674 | end module radiation_interface |
|---|
