From 129b829e2d1cd8cd837526d75423172510d37878 Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Tue, 28 May 2019 10:04:25 -0600 Subject: [PATCH] In progress... --- physics/GFS_rrtmgp_post.F90 | 395 +++++--- physics/GFS_rrtmgp_pre.F90 | 1769 ++++++++++++++++++----------------- physics/rrtmgp_lw.F90 | 131 ++- physics/rrtmgp_sw.F90 | 717 +++----------- 4 files changed, 1402 insertions(+), 1610 deletions(-) diff --git a/physics/GFS_rrtmgp_post.F90 b/physics/GFS_rrtmgp_post.F90 index 1d5002b93..150970c8b 100644 --- a/physics/GFS_rrtmgp_post.F90 +++ b/physics/GFS_rrtmgp_post.F90 @@ -9,14 +9,13 @@ module GFS_rrtmgp_post GFS_radtend_type, & GFS_diag_type use module_radiation_aerosols, only: NSPC1 - use module_radsw_parameters, only: cmpfsw_type use module_radlw_parameters, only: topflw_type, sfcflw_type, proflw_type - use module_radsw_parameters, only: topfsw_type, sfcfsw_type, profsw_type + use module_radsw_parameters, only: topfsw_type, sfcfsw_type, profsw_type, cmpfsw_type ! RRTMGP DDT's use mo_gas_optics_rrtmgp, only: ty_gas_optics_rrtmgp use mo_fluxes_byband, only: ty_fluxes_byband use mo_heating_rates, only: compute_heating_rate - use rrtmgp_lw, only: check_error_msg + use rrtmgp_sw, only: check_error_msg implicit none contains @@ -29,46 +28,63 @@ subroutine GFS_rrtmgp_post_init () end subroutine GFS_rrtmgp_post_init !> \section arg_table_GFS_rrtmgp_post_run Argument Table -!! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | -!! |----------------|-----------------------------------------------------------------------------------------------|-------------------------------------------------------------------------------|----------|------|----------------------|-----------|--------|----------| -!! | Model | GFS_control_type_instance | Fortran DDT containing FV3-GFS model control parameters | DDT | 0 | GFS_control_type | | in | F | -!! | Grid | GFS_grid_type_instance | Fortran DDT containing FV3-GFS grid and interpolation related data | DDT | 0 | GFS_grid_type | | in | F | -!! | Diag | GFS_diag_type_instance | Fortran DDT containing FV3-GFS diagnotics data | DDT | 0 | GFS_diag_type | | inout | F | -!! | Radtend | GFS_radtend_type_instance | Fortran DDT containing FV3-GFS radiation tendencies | DDT | 0 | GFS_radtend_type | | inout | F | -!! | Statein | GFS_statein_type_instance | Fortran DDT containing FV3-GFS prognostic state data in from dycore | DDT | 0 | GFS_statein_type | | in | F | -!! | Coupling | GFS_coupling_type_instance | Fortran DDT containing FV3-GFS fields to/from coupling with other components | DDT | 0 | GFS_coupling_type | | inout | F | -!! | scmpsw | components_of_surface_downward_shortwave_fluxes | derived type for special components of surface downward shortwave fluxes | W m-2 | 1 | cmpfsw_type | | in | F | -!! | im | horizontal_loop_extent | horizontal loop extent | count | 0 | integer | | in | F | -!! | lm | vertical_layer_dimension_for_radiation | number of vertical layers for radiation calculation | count | 0 | integer | | in | F | -!! | ltp | extra_top_layer | extra top layers | none | 0 | integer | | in | F | -!! | kt | vertical_index_difference_between_layer_and_upper_bound | vertical index difference between layer and upper bound | index | 0 | integer | | in | F | -!! | kb | vertical_index_difference_between_layer_and_lower_bound | vertical index difference between layer and lower bound | index | 0 | integer | | in | F | -!! | kd | vertical_index_difference_between_inout_and_local | vertical index difference between in/out and local | index | 0 | integer | | in | F | -!! | raddt | time_step_for_radiation | radiation time step | s | 0 | real | kind_phys | in | F | -!! | aerodp | atmosphere_optical_thickness_due_to_ambient_aerosol_particles | vertical integrated optical depth for various aerosol species | none | 2 | real | kind_phys | in | F | -!! | cldsa | cloud_area_fraction_for_radiation | fraction of clouds for low, middle, high, total and BL | frac | 2 | real | kind_phys | in | F | -!! | mtopa | model_layer_number_at_cloud_top | vertical indices for low, middle and high cloud tops | index | 2 | integer | | in | F | -!! | mbota | model_layer_number_at_cloud_base | vertical indices for low, middle and high cloud bases | index | 2 | integer | | in | F | -!! | cloud_fraction | total_cloud_fraction | layer total cloud fraction | frac | 2 | real | kind_phys | in | F | -!! | cldtaulw | cloud_optical_depth_layers_at_10mu_band | approx 10mu band layer cloud optical depth | none | 2 | real | kind_phys | in | F | -!! | cldtausw | cloud_optical_depth_layers_at_0.55mu_band | approx .55mu band layer cloud optical depth | none | 2 | real | kind_phys | in | F | -!! | tsfa | surface_air_temperature_for_radiation | lowest model layer air temperature for radiation | K | 1 | real | kind_phys | in | F | -!! | p_lev | air_pressure_at_interface_for_radiation_in_hPa | air pressure level | hPa | 2 | real | kind_phys | in | F | -!! | fluxLW_allsky | lw_flux_profiles_byband_allsky | Fortran DDT containing RRTMGP 3D fluxes | DDT | 0 | ty_fluxes_byband | | in | F | -!! | fluxLW_clrsky | lw_flux_profiles_byband_clrsky | Fortran DDT containing RRTMGP 3D fluxes | DDT | 0 | ty_fluxes_byband | | in | F | -!! | kdist_lw | K_distribution_file_for_RRTMGP_LW_scheme | DDT containing spectral information for RRTMGP LW radiation scheme | DDT | 0 | ty_gas_optics_rrtmgp | | in | F | -!! | hlwc | tendency_of_air_temperature_due_to_longwave_heating_on_radiation_time_step | longwave total sky heating rate | K s-1 | 2 | real | kind_phys | out | F | -!! | topflx | lw_fluxes_top_atmosphere | longwave total sky fluxes at the top of the atm | W m-2 | 1 | topflw_type | | inout | F | -!! | sfcflx | lw_fluxes_sfc | longwave total sky fluxes at the Earth surface | W m-2 | 1 | sfcflw_type | | inout | F | -!! | hlw0 | tendency_of_air_temperature_due_to_longwave_heating_assuming_clear_sky_on_radiation_time_step | longwave clear sky heating rate | K s-1 | 2 | real | kind_phys | inout | T | -!! | flxprf | lw_fluxes | lw fluxes total sky / csk and up / down at levels | W m-2 | 2 | proflw_type | | inout | T | -!! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | -!! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | +!! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | +!! |-------------------|------------------------------------------------------------------------------------------------|------------------------------------------------------------------------------|----------|------|----------------------|-----------|--------|----------| +!! | Model | GFS_control_type_instance | Fortran DDT containing FV3-GFS model control parameters | DDT | 0 | GFS_control_type | | in | F | +!! | Grid | GFS_grid_type_instance | Fortran DDT containing FV3-GFS grid and interpolation related data | DDT | 0 | GFS_grid_type | | in | F | +!! | Diag | GFS_diag_type_instance | Fortran DDT containing FV3-GFS diagnotics data | DDT | 0 | GFS_diag_type | | inout | F | +!! | Radtend | GFS_radtend_type_instance | Fortran DDT containing FV3-GFS radiation tendencies | DDT | 0 | GFS_radtend_type | | inout | F | +!! | Statein | GFS_statein_type_instance | Fortran DDT containing FV3-GFS prognostic state data in from dycore | DDT | 0 | GFS_statein_type | | in | F | +!! | Coupling | GFS_coupling_type_instance | Fortran DDT containing FV3-GFS fields to/from coupling with other components | DDT | 0 | GFS_coupling_type | | inout | F | +!! | scmpsw | components_of_surface_downward_shortwave_fluxes | derived type for special components of surface downward shortwave fluxes | W m-2 | 1 | cmpfsw_type | | inout | T | +!! | im | horizontal_loop_extent | horizontal loop extent | count | 0 | integer | | in | F | +!! | lm | vertical_layer_dimension_for_radiation | number of vertical layers for radiation calculation | count | 0 | integer | | in | F | +!! | ltp | extra_top_layer | extra top layers | none | 0 | integer | | in | F | +!! | kt | vertical_index_difference_between_layer_and_upper_bound | vertical index difference between layer and upper bound | index | 0 | integer | | in | F | +!! | kb | vertical_index_difference_between_layer_and_lower_bound | vertical index difference between layer and lower bound | index | 0 | integer | | in | F | +!! | kd | vertical_index_difference_between_inout_and_local | vertical index difference between in/out and local | index | 0 | integer | | in | F | +!! | raddt | time_step_for_radiation | radiation time step | s | 0 | real | kind_phys | in | F | +!! | aerodp | atmosphere_optical_thickness_due_to_ambient_aerosol_particles | vertical integrated optical depth for various aerosol species | none | 2 | real | kind_phys | in | F | +!! | cldsa | cloud_area_fraction_for_radiation | fraction of clouds for low, middle, high, total and BL | frac | 2 | real | kind_phys | in | F | +!! | mtopa | model_layer_number_at_cloud_top | vertical indices for low, middle and high cloud tops | index | 2 | integer | | in | F | +!! | mbota | model_layer_number_at_cloud_base | vertical indices for low, middle and high cloud bases | index | 2 | integer | | in | F | +!! | cloud_fraction | total_cloud_fraction | layer total cloud fraction | frac | 2 | real | kind_phys | in | F | +!! | cldtaulw | cloud_optical_depth_layers_at_10mu_band | approx 10mu band layer cloud optical depth | none | 2 | real | kind_phys | in | F | +!! | cldtausw | cloud_optical_depth_layers_at_0.55mu_band | approx .55mu band layer cloud optical depth | none | 2 | real | kind_phys | in | F | +!! | tsfa | surface_air_temperature_for_radiation | lowest model layer air temperature for radiation | K | 1 | real | kind_phys | in | F | +!! | p_lev | air_pressure_at_interface_for_radiation_in_hPa | air pressure level | hPa | 2 | real | kind_phys | in | F | +!! | nday | daytime_points_dimension | daytime points dimension | count | 0 | integer | | in | F | +!! | idxday | daytime_points | daytime points | index | 1 | integer | | in | F | +!! | fluxLW_allsky | lw_flux_profiles_byband_allsky | Fortran DDT containing RRTMGP 3D fluxes | DDT | 0 | ty_fluxes_byband | | in | F | +!! | fluxLW_clrsky | lw_flux_profiles_byband_clrsky | Fortran DDT containing RRTMGP 3D fluxes | DDT | 0 | ty_fluxes_byband | | in | F | +!! | fluxSW_allsky | sw_flux_profiles_byband_allsky | Fortran DDT containing RRTMGP 3D fluxes | DDT | 0 | ty_fluxes_byband | | in | F | +!! | fluxSW_clrsky | sw_flux_profiles_byband_clrsky | Fortran DDT containing RRTMGP 3D fluxes | DDT | 0 | ty_fluxes_byband | | in | F | +!! | kdist_lw | K_distribution_file_for_RRTMGP_LW_scheme | DDT containing spectral information for RRTMGP LW radiation scheme | DDT | 0 | ty_gas_optics_rrtmgp | | in | F | +!! | kdist_sw | K_distribution_file_for_RRTMGP_SW_scheme | DDT containing spectral information for RRTMGP SW radiation scheme | DDT | 0 | ty_gas_optics_rrtmgp | | in | F | +!! | sfc_alb_nir_dir | surface_shortwave_albedo_near_infrared_direct_in_each_band | surface sw near-infrared direct albedo in each SW band | frac | 2 | real | kind_phys | in | F | +!! | sfc_alb_nir_dif | surface_shortwave_albedo_near_infrared_diffuse_in_each_band | surface sw near-infrared diffuse albedo in each SW band | frac | 2 | real | kind_phys | in | F | +!! | sfc_alb_uvvis_dir | surface_shortwave_albedo_uv_visible_direct_in_each_band | surface sw uv-visible direct albedo in each SW band | frac | 2 | real | kind_phys | in | F | +!! | sfc_alb_uvvis_dif | surface_shortwave_albedo_uv_visible_diffuse_in_each_band | surface sw uv-visible diffuse albedo in each SW band | frac | 2 | real | kind_phys | in | F | +!! | hlwc | tendency_of_air_temperature_due_to_longwave_heating_on_radiation_time_step | longwave total sky heating rate | K s-1 | 2 | real | kind_phys | out | F | +!! | hswc | tendency_of_air_temperature_due_to_shortwave_heating_on_radiation_time_step | shortwave total sky heating rate | K s-1 | 2 | real | kind_phys | out | F | +!! | topflx_lw | lw_fluxes_top_atmosphere | longwave total sky fluxes at the top of the atm | W m-2 | 1 | topflw_type | | inout | F | +!! | sfcflx_lw | lw_fluxes_sfc | longwave total sky fluxes at the Earth surface | W m-2 | 1 | sfcflw_type | | inout | F | +!! | topflx_sw | sw_fluxes_top_atmosphere | shortwave total sky fluxes at the top of the atm | W m-2 | 1 | topfsw_type | | inout | F | +!! | sfcflx_sw | sw_fluxes_sfc | shortwave total sky fluxes at the Earth surface | W m-2 | 1 | sfcfsw_type | | inout | F | +!! | flxprf_lw | lw_fluxes | lw fluxes total sky / csk and up / down at levels | W m-2 | 2 | proflw_type | | inout | T | +!! | flxprf_sw | sw_fluxes | sw fluxes total sky / csk and up / down at levels | W m-2 | 2 | profsw_type | | inout | T | +!! | hlw0 | tendency_of_air_temperature_due_to_longwave_heating_assuming_clear_sky_on_radiation_time_step | longwave clear sky heating rate | K s-1 | 2 | real | kind_phys | inout | T | +!! | hsw0 | tendency_of_air_temperature_due_to_shortwave_heating_assuming_clear_sky_on_radiation_time_step | shortwave clear sky heating rate | K s-1 | 2 | real | kind_phys | inout | T | +!! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | +!! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | !! subroutine GFS_rrtmgp_post_run (Model, Grid, Diag, Radtend, Statein, & Coupling, scmpsw, im, lm, ltp, kt, kb, kd, raddt, aerodp, & - cldsa, mtopa, mbota, cloud_fraction, cldtaulw, cldtausw, p_lev, kdist_lw, & - tsfa, fluxLW_allsky, fluxLW_clrsky, hlwc, topflx, sfcflx, hlw0, flxprf,errmsg, errflg) + cldsa, mtopa, mbota, cloud_fraction, cldtaulw, cldtausw, p_lev, kdist_lw, kdist_sw, & + sfc_alb_nir_dir, sfc_alb_nir_dif, sfc_alb_uvvis_dir, & + sfc_alb_uvvis_dif, & + tsfa, nday, idxday, fluxSW_allsky, fluxSW_clrsky,fluxLW_allsky, fluxLW_clrsky,& + hlwc, hswc, topflx_sw, sfcflx_sw, flxprf_sw, topflx_lw, sfcflx_lw, flxprf_lw, hlw0, hsw0, errmsg, errflg) ! Inputs type(GFS_control_type), intent(in) :: & @@ -83,15 +99,17 @@ subroutine GFS_rrtmgp_post_run (Model, Grid, Diag, Radtend, Statein, & Radtend ! Fortran DDT containing FV3-GFS radiation tendencies type(GFS_diag_type), intent(inout) :: & Diag ! Fortran DDT containing FV3-GFS diagnotics data - type(cmpfsw_type), dimension(size(Grid%xlon,1)), intent(in) :: & - scmpsw ! derived type for special components of surface downward shortwave fluxes + integer, intent(in) :: & im, & ! Horizontal loop extent lm, & ! Number of vertical layers for radiation calculation ltp, & ! Extra-top-layers kt, & ! Vertical index difference between layer and upper bound kb, & ! Vertical index difference between layer and upper bound - kd ! Vertical index difference between in/out and local + kd, & ! Vertical index difference between in/out and local + nDay ! Number of daylit columns + integer, intent(in), dimension(nday) :: & + idxday ! Index array for daytime points real(kind_phys), intent(in) :: & raddt ! Radiation time step real(kind_phys), dimension(size(Grid%xlon,1)), intent(in) :: & @@ -108,26 +126,45 @@ subroutine GFS_rrtmgp_post_run (Model, Grid, Diag, Radtend, Statein, & cldtausw, & ! approx .55mu band layer cloud optical depth cldtaulw ! approx 10mu band layer cloud optical depth type(ty_gas_optics_rrtmgp),intent(in) :: & - kdist_lw ! DDT containing LW spectral information + kdist_lw, & ! DDT containing LW spectral information + kdist_sw ! DDT containing SW spectral information real(kind_phys), dimension(size(Grid%xlon,1), Model%levr+LTP+1), intent(in) :: & - p_lev ! Pressure @ model layer-interfaces (hPa) + p_lev ! Pressure @ model layer-interfaces (hPa) type(ty_fluxes_byband),intent(in) :: & - fluxLW_allsky, & ! All-sky flux (W/m2) - fluxLW_clrsky ! Clear-sky flux (W/m2) - - ! Outputs + fluxLW_allsky, & ! Longwave all-sky flux (W/m2) + fluxLW_clrsky, & ! Longwave clear-sky flux (W/m2) + fluxSW_allsky, & ! Shortwave all-sky flux (W/m2) + fluxSW_clrsky ! Shortwave clear-sky flux (W/m2) + real(kind_phys),dimension(kdist_sw%get_nband(),size(Grid%xlon,1)),intent(in) :: & + sfc_alb_nir_dir, & ! Shortwave surface albedo (nIR-direct) + sfc_alb_nir_dif, & ! Shortwave surface albedo (nIR-diffuse) + sfc_alb_uvvis_dir, & ! Shortwave surface albedo (uvvis-direct) + sfc_alb_uvvis_dif ! Shortwave surface albedo (uvvis-diffuse) + + ! Outputs (mandatory) character(len=*), intent(out) :: & errmsg integer, intent(out) :: & errflg real(kind_phys),dimension(size(Grid%xlon,1), Model%levr+LTP),intent(out) :: & - hlwc ! All-sky heating-rate (K/sec) + hlwc, & ! Longwave all-sky heating-rate (K/sec) + hswc ! Shortwave all-sky heating-rate (K/sec) type(topflw_type), dimension(size(Grid%xlon,1)), intent(inout) :: & - topflx ! radiation fluxes at top, components: + topflx_lw ! radiation fluxes at top, components: ! upfxc - total sky upward flux at top (w/m2) ! upfx0 - clear sky upward flux at top (w/m2) type(sfcflw_type), dimension(size(Grid%xlon,1)), intent(inout) :: & - sfcflx ! radiation fluxes at sfc, components: + sfcflx_lw ! radiation fluxes at sfc, components: + ! upfxc - total sky upward flux at sfc (w/m2) + ! upfx0 - clear sky upward flux at sfc (w/m2) + ! dnfxc - total sky downward flux at sfc (w/m2) + ! dnfx0 - clear sky downward flux at sfc (w/m2) + type(topfsw_type), dimension(size(Grid%xlon,1)), intent(inout) :: & + topflx_sw ! radiation fluxes at top, components: + ! upfxc - total sky upward flux at top (w/m2) + ! upfx0 - clear sky upward flux at top (w/m2) + type(sfcfsw_type), dimension(size(Grid%xlon,1)), intent(inout) :: & + sfcflx_sw ! radiation fluxes at sfc, components: ! upfxc - total sky upward flux at sfc (w/m2) ! upfx0 - clear sky upward flux at sfc (w/m2) ! dnfxc - total sky downward flux at sfc (w/m2) @@ -135,20 +172,34 @@ subroutine GFS_rrtmgp_post_run (Model, Grid, Diag, Radtend, Statein, & ! Outputs (optional) real(kind_phys), dimension(size(Grid%xlon,1), Model%levr+LTP), optional, intent(inout) :: & - hlw0 ! Clear-sky heating rate (K/sec) + hlw0, & ! Longwave clear-sky heating rate (K/sec) + hsw0 ! Shortwave clear-sky heating-rate (K/sec) type(proflw_type), dimension(size(Grid%xlon,1), Model%levr+LTP+1), optional, intent(inout) :: & - flxprf ! 2D radiative fluxes, components: - ! upfxc - total sky upward flux (W/m2) - ! dnfxc - total sky dnward flux (W/m2) - ! upfx0 - clear sky upward flux (W/m2) - ! dnfx0 - clear sky dnward flux (W/m2) - + flxprf_lw ! 2D radiative fluxes, components: + ! upfxc - total sky upward flux (W/m2) + ! dnfxc - total sky dnward flux (W/m2) + ! upfx0 - clear sky upward flux (W/m2) + ! dnfx0 - clear sky dnward flux (W/m2) + type(profsw_type), dimension(size(Grid%xlon,1), Model%levr+LTP+1), intent(inout), optional :: & + flxprf_sw ! 2D radiative fluxes, components: + ! upfxc - total sky upward flux (W/m2) + ! dnfxc - total sky dnward flux (W/m2) + ! upfx0 - clear sky upward flux (W/m2) + ! dnfx0 - clear sky dnward flux (W/m2) + type(cmpfsw_type), dimension(size(Grid%xlon,1)), intent(inout), optional :: & + scmpsw ! 2D surface fluxes, components: + ! uvbfc - total sky downward uv-b flux at (W/m2) + ! uvbf0 - clear sky downward uv-b flux at (W/m2) + ! nirbm - downward nir direct beam flux (W/m2) + ! nirdf - downward nir diffused flux (W/m2) + ! visbm - downward uv+vis direct beam flux (W/m2) + ! visdf - downward uv+vis diffused flux (W/m2) ! Local variables integer :: i, j, k, k1, itop, ibtc, iBand, iSFC, iTOA real(kind_phys) :: tem0d, tem1, tem2 - real(kind_phys), dimension(size(Grid%xlon,1), Model%levr+LTP) :: thetaTendClrSky, thetaTendAllSky - logical :: l_ClrSky_HR, l_fluxes2D, top_at_1 - + real(kind_phys), dimension(nDay, Model%levr+LTP) :: thetaTendClrSky, thetaTendAllSky + logical :: l_clrskylw_hr,l_clrskysw_hr, l_fluxeslw2d, l_fluxessw2d, top_at_1, l_sfcFluxessw1D + ! Initialize CCPP error handling variables errmsg = '' errflg = 0 @@ -156,8 +207,12 @@ subroutine GFS_rrtmgp_post_run (Model, Grid, Diag, Radtend, Statein, & if (.not. (Model%lsswr .or. Model%lslwr)) return ! Are any optional outputs requested? - l_ClrSky_HR = present(hlw0) - l_fluxes2D = present(flxprf) + l_clrskylw_hr = present(hlw0) + l_fluxeslw2d = present(flxprf_lw) + l_clrskysw_hr = present(hsw0) + l_fluxessw2d = present(flxprf_sw) + l_sfcfluxessw1D = present(scmpsw) + ! What is vertical ordering? top_at_1 = (p_lev(1,1) .lt. p_lev(1, Model%levr+LTP)) @@ -174,56 +229,50 @@ subroutine GFS_rrtmgp_post_run (Model, Grid, Diag, Radtend, Statein, & ! ####################################################################################### if (Model%lslwr) then ! Clear-sky heating-rate (optional) - if (l_ClrSky_HR) then - call check_error_msg(compute_heating_rate( & + if (l_clrskylw_hr) then + call check_error_msg('GFS_rrtmgp_post',compute_heating_rate( & fluxLW_clrsky%flux_up, & fluxLW_clrsky%flux_dn, & p_lev, & - thetaTendClrSky)) + hlw0)) endif ! All-sky heating-rate (mandatory) - call check_error_msg(compute_heating_rate( & + call check_error_msg('GFS_rrtmgp_post',compute_heating_rate( & fluxLW_allsky%flux_up, & fluxLW_allsky%flux_dn, & p_lev, & - thetaTendAllSky)) + hlwc)) ! Copy fluxes from RRTGMP types into model radiation types. ! Mandatory outputs - topflx%upfxc = fluxLW_allsky%flux_up(:,iTOA) - topflx%upfx0 = fluxLW_clrsky%flux_up(:,iTOA) - sfcflx%upfxc = fluxLW_allsky%flux_up(:,iSFC) - sfcflx%upfx0 = fluxLW_clrsky%flux_up(:,iSFC) - sfcflx%dnfxc = fluxLW_allsky%flux_dn(:,iSFC) - sfcflx%dnfx0 = fluxLW_clrsky%flux_dn(:,iSFC) - hlwc = thetaTendAllSky + topflx_lw%upfxc = fluxLW_allsky%flux_up(:,iTOA) + topflx_lw%upfx0 = fluxLW_clrsky%flux_up(:,iTOA) + sfcflx_lw%upfxc = fluxLW_allsky%flux_up(:,iSFC) + sfcflx_lw%upfx0 = fluxLW_clrsky%flux_up(:,iSFC) + sfcflx_lw%dnfxc = fluxLW_allsky%flux_dn(:,iSFC) + sfcflx_lw%dnfx0 = fluxLW_clrsky%flux_dn(:,iSFC) ! Optional outputs - if(l_fluxes2D) then - flxprf%upfxc = fluxLW_allsky%flux_up - flxprf%dnfxc = fluxLW_allsky%flux_dn - flxprf%upfx0 = fluxLW_clrsky%flux_up - flxprf%dnfx0 = fluxLW_clrsky%flux_dn - endif - if (l_ClrSky_HR) then - hlw0 = thetaTendClrSky + if(l_fluxeslw2d) then + flxprf_lw%upfxc = fluxLW_allsky%flux_up + flxprf_lw%dnfxc = fluxLW_allsky%flux_dn + flxprf_lw%upfx0 = fluxLW_clrsky%flux_up + flxprf_lw%dnfx0 = fluxLW_clrsky%flux_dn endif endif ! ####################################################################################### - ! Save LW heating-rates (Note. This piece was originally in rrtmg_lw_post.F90:_run()) + ! Save LW outputs (Note. This piece was originally in rrtmg_lw_post.F90:_run()) ! ####################################################################################### if (Model%lslwr) then - !> -# Save calculation results - !> - Save surface air temp for diurnal adjustment at model t-steps - + ! Save surface air temp for diurnal adjustment at model t-steps Radtend%tsflw (:) = tsfa(:) do k = 1, LM k1 = k + kd Radtend%htrlw(1:im,k) = hlwc(1:im,k1) enddo - ! --- repopulate the points above levr + ! Repopulate the points above levr if (lm < Model%levs) then do k = lm,Model%levs Radtend%htrlw (1:im,k) = Radtend%htrlw (1:im,LM) @@ -235,7 +284,7 @@ subroutine GFS_rrtmgp_post_run (Model, Grid, Diag, Radtend, Statein, & k1 = k + kd Radtend%lwhc(1:im,k) = hlw0(1:im,k1) enddo - ! --- repopulate the points above levr + ! Repopulate the points above levr if (lm < Model%levs) then do k = lm,Model%levs Radtend%lwhc(1:im,k) = Radtend%lwhc(1:im,LM) @@ -243,19 +292,148 @@ subroutine GFS_rrtmgp_post_run (Model, Grid, Diag, Radtend, Statein, & endif endif - ! --- radiation fluxes for other physics processes + ! Radiation fluxes for other physics processes Coupling%sfcdlw(:) = Radtend%sfcflw(:)%dnfxc - endif ! end_if_lslwr - + endif + + ! ####################################################################################### + ! Compute SW heating-rates + ! ####################################################################################### + ! Initialize outputs + hswc(:,:) = 0. + topflx_sw = topfsw_type ( 0., 0., 0. ) + sfcflx_sw = sfcfsw_type ( 0., 0., 0., 0. ) + if (l_clrskysw_hr) then + hsw0(:,:) = 0. + endif + if (l_fluxessw2D) then + flxprf_sw = profsw_type ( 0., 0., 0., 0. ) + endif + if (l_sfcfluxessw1D) then + scmpsw = cmpfsw_type (0.,0.,0.,0.,0.,0.) + endif + + if (Model%lsswr .and. nDay .gt. 0) then + ! Clear-sky heating-rate (optional) + if (l_clrskysw_HR) then + call check_error_msg('GFS_rrtmgp_post',compute_heating_rate( & + fluxSW_clrsky%flux_up, & + fluxSW_clrsky%flux_dn, & + p_lev(idxday,1:Model%levr+LTP+1), & + thetaTendClrSky)) + hsw0(idxday,:)=thetaTendClrSky + endif + ! All-sky heating-rate (mandatory) + call check_error_msg('GFS_rrtmgp_post',compute_heating_rate( & + fluxSW_allsky%flux_up, & + fluxSW_allsky%flux_dn, & + p_lev(idxday,1:Model%levr+LTP+1), & + thetaTendAllSky)) + hswc(idxday,:) = thetaTendAllSky + print*,'IN POST: ',fluxSW_allsky%flux_up + print*,'IN POSTT: ',fluxSW_allsky%flux_dn + ! Copy fluxes from RRTGMP types into model radiation types. + ! Mandatory outputs + topflx_sw(idxday)%upfxc = fluxSW_allsky%flux_up(:,iTOA) + topflx_sw(idxday)%upfx0 = fluxSW_clrsky%flux_up(:,iTOA) + sfcflx_sw(idxday)%upfxc = fluxSW_allsky%flux_up(:,iSFC) + sfcflx_sw(idxday)%upfx0 = fluxSW_clrsky%flux_up(:,iSFC) + sfcflx_sw(idxday)%dnfxc = fluxSW_allsky%flux_dn(:,iSFC) + sfcflx_sw(idxday)%dnfx0 = fluxSW_clrsky%flux_dn(:,iSFC) + + ! Optional output + if(l_fluxessw2D) then + flxprf_sw(idxday,:)%upfxc = fluxSW_allsky%flux_up + flxprf_sw(idxday,:)%dnfxc = fluxSW_allsky%flux_dn + flxprf_sw(idxday,:)%upfx0 = fluxSW_clrsky%flux_up + flxprf_sw(idxday,:)%dnfx0 = fluxSW_clrsky%flux_dn + endif + endif + + ! ####################################################################################### + ! Save SW outputs + ! ####################################################################################### + if (Model%lsswr) then + if (nday > 0) then + do k = 1, LM + k1 = k + kd + Radtend%htrsw(1:im,k) = hswc(1:im,k1) + enddo + ! Repopulate the points above levr i.e. LM + if (lm < Model%levs) then + do k = lm,Model%levs + Radtend%htrsw (1:im,k) = Radtend%htrsw (1:im,LM) + enddo + endif + + if (Model%swhtr) then + do k = 1, lm + k1 = k + kd + Radtend%swhc(1:im,k) = hsw0(1:im,k1) + enddo + ! Repopulate the points above levr i.e. LM + if (lm < Model%levs) then + do k = lm,Model%levs + Radtend%swhc(1:im,k) = Radtend%swhc(1:im,LM) + enddo + endif + endif + + ! Surface down and up spectral component fluxes + ! - Save two spectral bands' surface downward and upward fluxes for output. + do i=1,im + Coupling%nirbmdi(i) = scmpsw(i)%nirbm + Coupling%nirdfdi(i) = scmpsw(i)%nirdf + Coupling%visbmdi(i) = scmpsw(i)%visbm + Coupling%visdfdi(i) = scmpsw(i)%visdf + + Coupling%nirbmui(i) = scmpsw(i)%nirbm * sfc_alb_nir_dir(1,i) + Coupling%nirdfui(i) = scmpsw(i)%nirdf * sfc_alb_nir_dif(1,i) + Coupling%visbmui(i) = scmpsw(i)%visbm * sfc_alb_uvvis_dir(1,i) + Coupling%visdfui(i) = scmpsw(i)%visdf * sfc_alb_uvvis_dif(1,i) + enddo + else ! if_nday_block + Radtend%htrsw(:,:) = 0.0 + Radtend%sfcfsw = sfcfsw_type( 0.0, 0.0, 0.0, 0.0 ) + Diag%topfsw = topfsw_type( 0.0, 0.0, 0.0 ) + scmpsw = cmpfsw_type( 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 ) + + do i=1,im + Coupling%nirbmdi(i) = 0.0 + Coupling%nirdfdi(i) = 0.0 + Coupling%visbmdi(i) = 0.0 + Coupling%visdfdi(i) = 0.0 + + Coupling%nirbmui(i) = 0.0 + Coupling%nirdfui(i) = 0.0 + Coupling%visbmui(i) = 0.0 + Coupling%visdfui(i) = 0.0 + enddo + + if (Model%swhtr) then + Radtend%swhc(:,:) = 0 + endif + endif ! end_if_nday + + ! Radiation fluxes for other physics processes + do i=1,im + Coupling%sfcnsw(i) = Radtend%sfcfsw(i)%dnfxc - Radtend%sfcfsw(i)%upfxc + Coupling%sfcdsw(i) = Radtend%sfcfsw(i)%dnfxc + enddo + endif ! end_if_lsswr + + ! ####################################################################################### + ! ####################################################################################### !> - For time averaged output quantities (including total-sky and !! clear-sky SW and LW fluxes at TOA and surface; conventional !! 3-domain cloud amount, cloud top and base pressure, and cloud top !! temperature; aerosols AOD, etc.), store computed results in !! corresponding slots of array fluxr with appropriate time weights. - ! --- ... collect the fluxr data for wrtsfc - + ! ####################################################################################### + ! ####################################################################################### + if (Model%lssav) then if (Model%lsswr) then do i=1,im @@ -268,47 +446,42 @@ subroutine GFS_rrtmgp_post_run (Model, Grid, Diag, Radtend, Statein, & enddo endif - ! --- save lw toa and sfc fluxes + ! Save LW TOA and SFC fluxes if (Model%lslwr) then do i=1,im - ! --- lw total-sky fluxes + ! LW all-sky fluxes Diag%fluxr(i,1 ) = Diag%fluxr(i,1 ) + Model%fhlwr * Diag%topflw(i)%upfxc ! total sky top lw up Diag%fluxr(i,19) = Diag%fluxr(i,19) + Model%fhlwr * Radtend%sfcflw(i)%dnfxc ! total sky sfc lw dn Diag%fluxr(i,20) = Diag%fluxr(i,20) + Model%fhlwr * Radtend%sfcflw(i)%upfxc ! total sky sfc lw up - ! --- lw clear-sky fluxes + ! LW clear-sky fluxes Diag%fluxr(i,28) = Diag%fluxr(i,28) + Model%fhlwr * Diag%topflw(i)%upfx0 ! clear sky top lw up Diag%fluxr(i,30) = Diag%fluxr(i,30) + Model%fhlwr * Radtend%sfcflw(i)%dnfx0 ! clear sky sfc lw dn Diag%fluxr(i,33) = Diag%fluxr(i,33) + Model%fhlwr * Radtend%sfcflw(i)%upfx0 ! clear sky sfc lw up enddo endif - ! --- save sw toa and sfc fluxes with proper diurnal sw wgt. coszen=mean cosz over daylight - ! part of sw calling interval, while coszdg= mean cosz over entire interval + ! Save sw toa and sfc fluxes with proper diurnal sw wgt. coszen=mean cosz over daylight + ! part of sw calling interval, while coszdg= mean cosz over entire interval if (Model%lsswr) then do i = 1, IM if (Radtend%coszen(i) > 0.) then - ! --- sw total-sky fluxes - ! ------------------- + ! SW all-sky fluxes tem0d = Model%fhswr * Radtend%coszdg(i) / Radtend%coszen(i) Diag%fluxr(i,2 ) = Diag%fluxr(i,2) + Diag%topfsw(i)%upfxc * tem0d ! total sky top sw up Diag%fluxr(i,3 ) = Diag%fluxr(i,3) + Radtend%sfcfsw(i)%upfxc * tem0d ! total sky sfc sw up Diag%fluxr(i,4 ) = Diag%fluxr(i,4) + Radtend%sfcfsw(i)%dnfxc * tem0d ! total sky sfc sw dn - ! --- sw uv-b fluxes - ! -------------- + ! SW uv-b fluxes Diag%fluxr(i,21) = Diag%fluxr(i,21) + scmpsw(i)%uvbfc * tem0d ! total sky uv-b sw dn Diag%fluxr(i,22) = Diag%fluxr(i,22) + scmpsw(i)%uvbf0 * tem0d ! clear sky uv-b sw dn - ! --- sw toa incoming fluxes - ! ---------------------- + ! SW TOA incoming fluxes Diag%fluxr(i,23) = Diag%fluxr(i,23) + Diag%topfsw(i)%dnfxc * tem0d ! top sw dn - ! --- sw sfc flux components - ! ---------------------- + ! SW SFC flux components Diag%fluxr(i,24) = Diag%fluxr(i,24) + scmpsw(i)%visbm * tem0d ! uv/vis beam sw dn Diag%fluxr(i,25) = Diag%fluxr(i,25) + scmpsw(i)%visdf * tem0d ! uv/vis diff sw dn Diag%fluxr(i,26) = Diag%fluxr(i,26) + scmpsw(i)%nirbm * tem0d ! nir beam sw dn Diag%fluxr(i,27) = Diag%fluxr(i,27) + scmpsw(i)%nirdf * tem0d ! nir diff sw dn - ! --- sw clear-sky fluxes - ! ------------------- - Diag%fluxr(i,29) = Diag%fluxr(i,29) + Diag%topfsw(i)%upfx0 * tem0d ! clear sky top sw up + ! SW clear-sky fluxes + Diag%fluxr(i,29) = Diag%fluxr(i,29) + Diag%topfsw(i)%upfx0 * tem0d ! clear sky top sw up Diag%fluxr(i,31) = Diag%fluxr(i,31) + Radtend%sfcfsw(i)%upfx0 * tem0d ! clear sky sfc sw up Diag%fluxr(i,32) = Diag%fluxr(i,32) + Radtend%sfcfsw(i)%dnfx0 * tem0d ! clear sky sfc sw dn endif diff --git a/physics/GFS_rrtmgp_pre.F90 b/physics/GFS_rrtmgp_pre.F90 index 07cc4d9de..fcb0e0fb0 100644 --- a/physics/GFS_rrtmgp_pre.F90 +++ b/physics/GFS_rrtmgp_pre.F90 @@ -1,851 +1,944 @@ !> \file GFS_rrtmgp_pre.f90 !! This file contains - module GFS_rrtmgp_pre - use machine, only: kind_phys - - real(kind_phys), parameter :: & - amd = 28.9644_kind_phys, & ! Molecular weight of dry-air (g/mol) - amw = 18.0154_kind_phys, & ! Molecular weight of water vapor (g/mol) - amo3 = 47.9982_kind_phys, & ! Modelular weight of ozone (g/mol) - amdw = amd/amw, & ! Molecular weight of dry air / water vapor - amdo3 = amd/amo3 ! Molecular weight of dry air / ozone - - public GFS_rrtmgp_pre_run,GFS_rrtmgp_pre_init,GFS_rrtmgp_pre_finalize - - contains - +module GFS_rrtmgp_pre + use physparam + use machine, only: & + kind_phys ! Working type + use GFS_typedefs, only: & + GFS_statein_type, & ! Prognostic state data in from dycore + GFS_stateout_type, & ! Prognostic state or tendencies return to dycore + GFS_sfcprop_type, & ! Surface fields + GFS_coupling_type, & ! Fields to/from coupling with other components (e.g. land/ice/ocean/etc.) + GFS_control_type, & ! Model control parameters + GFS_grid_type, & ! Grid and interpolation related data + GFS_tbd_type, & ! To-Be-Determined data that doesn't fit in any one container + GFS_radtend_type, & ! Radiation tendencies needed in physics + GFS_diag_type ! Fields targetted for diagnostic output + use physcons, only: & + eps => con_eps, & ! Rd/Rv + epsm1 => con_epsm1, & ! Rd/Rv-1 + fvirt => con_fvirt, & ! Rv/Rd-1 + rog => con_rog, & ! Rd/g + rocp => con_rocp ! Rd/cp + use radcons, only: & + itsfc, & ! Flag for LW sfc. temp. + ltp, & ! 1-add extra-top layer; 0-no extra layer + lextop, & ! ltp > 0 + qmin,qme5, qme6, epsq ! Minimum vlaues for varius calculations + use funcphys, only: & + fpvs ! Function ot compute sat. vapor pressure over liq. + use module_radiation_astronomy,only: & + coszmn ! Function to compute cos(SZA) + use module_radiation_gases, only: & + NF_VGAS, & ! Number of active gas species + getgases, & ! Routine to setup trace gases + getozn ! Routine to setup ozone + use module_radiation_aerosols, only: & + NF_AESW, & ! Number of optical-fields in SW output (3=tau+g+omega) + NF_AELW, & ! Number of optical-fields in LW output (3=tau+g+omega) + setaer, & ! Routine to compute aerosol radiative properties (tau,g,omega) + NSPC1 ! Number of species for vertically integrated aerosol optical-depth + use module_radiation_clouds, only: & + NF_CLDS, & ! Number of fields in "clouds" array (e.g. (cloud(1)=lwp,clouds(2)=ReffLiq,...) + progcld1, & ! Zhao/Moorthi's prognostic cloud scheme + progcld3, & ! Zhao/Moorthi's prognostic cloud+pdfcld + progcld4, & ! GFDL cloud scheme + progcld5, & ! Thompson / WSM6 cloud micrphysics scheme + progclduni ! Unified cloud-scheme + use surface_perturbation, only: & + cdfnor ! Routine to compute CDF (used to compute percentiles) + use module_radiation_surface, only: & + setemis, & ! Routine to compute surface-emissivity + NF_ALBD, & ! Number of surface albedo categories (4; nir-direct, nir-diffuse, uvvis-direct, uvvis-diffuse) + setalb ! Routine to compute surface albedo + + use rrtmgp_lw, only: nrghice_lw => nrghice, ipsdlw0 + use rrtmgp_sw, only: nrghice_sw => nrghice, ipsdsw0, check_error_msg + use mersenne_twister, only: & + random_setseed, & + random_number, & + random_stat + ! RRTMGP types + use mo_gas_optics_rrtmgp, only: ty_gas_optics_rrtmgp + use mo_cloud_optics, only: ty_cloud_optics + use mo_optical_props, only: ty_optical_props_1scl, ty_optical_props_2str + use mo_cloud_sampling, only: sampled_mask_max_ran, sampled_mask_exp_ran, draw_samples + use mo_gas_concentrations, only: ty_gas_concs + + real(kind_phys), parameter :: & + amd = 28.9644_kind_phys, & ! Molecular weight of dry-air (g/mol) + amw = 18.0154_kind_phys, & ! Molecular weight of water vapor (g/mol) + amo3 = 47.9982_kind_phys, & ! Modelular weight of ozone (g/mol) + amdw = amd/amw, & ! Molecular weight of dry air / water vapor + amdo3 = amd/amo3 ! Molecular weight of dry air / ozone + + public GFS_rrtmgp_pre_run,GFS_rrtmgp_pre_init,GFS_rrtmgp_pre_finalize + +contains + !> \defgroup GFS_rrtmgp_pre GFS RRTMGP Scheme Pre !! @{ !! \section arg_table_GFS_rrtmgp_pre_init Argument Table !! - subroutine GFS_rrtmgp_pre_init () - open(58,file='GFS_rrtmgp_aux_dump.txt',status='unknown') - end subroutine GFS_rrtmgp_pre_init + subroutine GFS_rrtmgp_pre_init () + end subroutine GFS_rrtmgp_pre_init !> \section arg_table_GFS_rrtmgp_pre_run Argument Table -!! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | -!! |-----------------------|---------------------------------------------------------------|-------------------------------------------------------------------------------|----------|------|------------------|-----------|--------|----------| -!! | Model | GFS_control_type_instance | Fortran DDT containing FV3-GFS model control parameters | DDT | 0 | GFS_control_type | | in | F | -!! | Grid | GFS_grid_type_instance | Fortran DDT containing FV3-GFS grid and interpolation related data | DDT | 0 | GFS_grid_type | | in | F | -!! | Sfcprop | GFS_sfcprop_type_instance | Fortran DDT containing FV3-GFS surface fields | DDT | 0 | GFS_sfcprop_type | | in | F | -!! | Statein | GFS_statein_type_instance | Fortran DDT containing FV3-GFS prognostic state data in from dycore | DDT | 0 | GFS_statein_type | | in | F | -!! | Tbd | GFS_tbd_type_instance | Fortran DDT containing FV3-GFS data not yet assigned to a defined container | DDT | 0 | GFS_tbd_type | | in | F | -!! | Cldprop | GFS_cldprop_type_instance | Fortran DDT containing FV3-GFS cloud fields needed by radiation from physics | DDT | 0 | GFS_cldprop_type | | in | F | -!! | Coupling | GFS_coupling_type_instance | Fortran DDT containing FV3-GFS fields needed for coupling | DDT | 0 | GFS_coupling_type| | in | F | -!! | Radtend | GFS_radtend_type_instance | Fortran DDT containing FV3-GFS radiation tendencies | DDT | 0 | GFS_radtend_type | | inout | F | -!! | lm | vertical_layer_dimension_for_radiation | number of vertical layers for radiation calculation | count | 0 | integer | | in | F | -!! | im | horizontal_loop_extent | horizontal loop extent | count | 0 | integer | | in | F | -!! | lmk | adjusted_vertical_layer_dimension_for_radiation | number of vertical layers for radiation | count | 0 | integer | | in | F | -!! | lmp | adjusted_vertical_level_dimension_for_radiation | number of vertical levels for radiation | count | 0 | integer | | in | F | -!! | kd | vertical_index_difference_between_inout_and_local | vertical index difference between in/out and local | index | 0 | integer | | out | F | -!! | kt | vertical_index_difference_between_layer_and_upper_bound | vertical index difference between layer and upper bound | index | 0 | integer | | out | F | -!! | kb | vertical_index_difference_between_layer_and_lower_bound | vertical index difference between layer and lower bound | index | 0 | integer | | out | F | -!! | raddt | time_step_for_radiation | radiation time step | s | 0 | real | kind_phys | out | F | -!! | delp | layer_pressure_thickness_for_radiation | layer pressure thickness on radiation levels | hPa | 2 | real | kind_phys | out | F | -!! | dz | layer_thickness_for_radiation | layer thickness on radiation levels | km | 2 | real | kind_phys | out | F | -!! | plvl | air_pressure_at_interface_for_radiation_in_hPa | air pressure at vertical interface for radiation calculation | hPa | 2 | real | kind_phys | out | F | -!! | plyr | air_pressure_at_layer_for_radiation_in_hPa | air pressure at vertical layer for radiation calculation | hPa | 2 | real | kind_phys | out | F | -!! | tlvl | air_temperature_at_interface_for_radiation | air temperature at vertical interface for radiation calculation | K | 2 | real | kind_phys | out | F | -!! | tlyr | air_temperature_at_layer_for_radiation | air temperature at vertical layer for radiation calculation | K | 2 | real | kind_phys | out | F | -!! | tsfg | surface_ground_temperature_for_radiation | surface ground temperature for radiation | K | 1 | real | kind_phys | out | F | -!! | tsfa | surface_air_temperature_for_radiation | lowest model layer air temperature for radiation | K | 1 | real | kind_phys | out | F | -!! | qlyr | water_vapor_specific_humidity_at_layer_for_radiation | water vapor specific humidity at vertical layer for radiation calculation | kg kg-1 | 2 | real | kind_phys | out | F | -!! | olyr | ozone_concentration_at_layer_for_radiation | ozone concentration | kg kg-1 | 2 | real | kind_phys | out | F | -!! | icseed | seed_random_numbers_lw | seed for random number generation for longwave radiation | none | 1 | integer | | in | F | -!! | gasvmr_co2 | volume_mixing_ratio_co2 | CO2 volume mixing ratio | kg kg-1 | 2 | real | kind_phys | out | F | -!! | gasvmr_n2o | volume_mixing_ratio_n2o | N2O volume mixing ratio | kg kg-1 | 2 | real | kind_phys | out | F | -!! | gasvmr_ch4 | volume_mixing_ratio_ch4 | CH4 volume mixing ratio | kg kg-1 | 2 | real | kind_phys | out | F | -!! | gasvmr_o2 | volume_mixing_ratio_o2 | O2 volume mixing ratio | kg kg-1 | 2 | real | kind_phys | out | F | -!! | gasvmr_co | volume_mixing_ratio_co | CO volume mixing ratio | kg kg-1 | 2 | real | kind_phys | out | F | -!! | gasvmr_cfc11 | volume_mixing_ratio_cfc11 | CFC11 volume mixing ratio | kg kg-1 | 2 | real | kind_phys | out | F | -!! | gasvmr_cfc12 | volume_mixing_ratio_cfc12 | CFC12 volume mixing ratio | kg kg-1 | 2 | real | kind_phys | out | F | -!! | gasvmr_cfc22 | volume_mixing_ratio_cfc22 | CFC22 volume mixing ratio | kg kg-1 | 2 | real | kind_phys | out | F | -!! | gasvmr_ccl4 | volume_mixing_ratio_ccl4 | CCL4 volume mixing ratio | kg kg-1 | 2 | real | kind_phys | out | F | -!! | gasvmr_cfc113 | volume_mixing_ratio_cfc113 | CFC113 volume mixing ratio | kg kg-1 | 2 | real | kind_phys | out | F | -!! | faersw1 | aerosol_optical_depth_for_shortwave_bands_01-16 | aerosol optical depth for shortwave bands 01-16 | none | 3 | real | kind_phys | out | F | -!! | faersw2 | aerosol_single_scattering_albedo_for_shortwave_bands_01-16 | aerosol single scattering albedo for shortwave bands 01-16 | frac | 3 | real | kind_phys | out | F | -!! | faersw3 | aerosol_asymmetry_parameter_for_shortwave_bands_01-16 | aerosol asymmetry parameter for shortwave bands 01-16 | none | 3 | real | kind_phys | out | F | -!! | faerlw1 | aerosol_optical_depth_for_longwave_bands_01-16 | aerosol optical depth for longwave bands 01-16 | none | 3 | real | kind_phys | out | F | -!! | faerlw2 | aerosol_single_scattering_albedo_for_longwave_bands_01-16 | aerosol single scattering albedo for longwave bands 01-16 | frac | 3 | real | kind_phys | out | F | -!! | faerlw3 | aerosol_asymmetry_parameter_for_longwave_bands_01-16 | aerosol asymmetry parameter for longwave bands 01-16 | none | 3 | real | kind_phys | out | F | -!! | aerodp | atmosphere_optical_thickness_due_to_ambient_aerosol_particles | vertical integrated optical depth for various aerosol species | none | 2 | real | kind_phys | out | F | -!! | clouds1 | total_cloud_fraction | layer total cloud fraction | frac | 2 | real | kind_phys | out | F | -!! | clouds2 | cloud_liquid_water_path | layer cloud liquid water path | g m-2 | 2 | real | kind_phys | out | F | -!! | clouds3 | mean_effective_radius_for_liquid_cloud | mean effective radius for liquid cloud | micron | 2 | real | kind_phys | out | F | -!! | clouds4 | cloud_ice_water_path | layer cloud ice water path | g m-2 | 2 | real | kind_phys | out | F | -!! | clouds5 | mean_effective_radius_for_ice_cloud | mean effective radius for ice cloud | micron | 2 | real | kind_phys | out | F | -!! | clouds6 | cloud_rain_water_path | cloud rain water path | g m-2 | 2 | real | kind_phys | out | F | -!! | clouds7 | mean_effective_radius_for_rain_drop | mean effective radius for rain drop | micron | 2 | real | kind_phys | out | F | -!! | clouds8 | cloud_snow_water_path | cloud snow water path | g m-2 | 2 | real | kind_phys | out | F | -!! | clouds9 | mean_effective_radius_for_snow_flake | mean effective radius for snow flake | micron | 2 | real | kind_phys | out | F | -!! | cldsa | cloud_area_fraction_for_radiation | fraction of clouds for low, middle,high, total and BL | frac | 2 | real | kind_phys | out | F | -!! | mtopa | model_layer_number_at_cloud_top | vertical indices for low, middle and high cloud tops | index | 2 | integer | | out | F | -!! | mbota | model_layer_number_at_cloud_base | vertical indices for low, middle and high cloud bases | index | 2 | integer | | out | F | -!! | de_lgth | cloud_decorrelation_length | cloud decorrelation length | km | 1 | real | kind_phys | out | F | -!! | alb1d | surface_albedo_perturbation | surface albedo perturbation | frac | 1 | real | kind_phys | out | F | -!! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | -!! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | -!! | kdist_lw | K_distribution_file_for_RRTMGP_LW_scheme | DDT containing spectral information for RRTMGP LW radiation scheme | DDT | 0 | ty_gas_optics_rrtmgp | | in | F | -!! | kdist_sw | K_distribution_file_for_RRTMGP_SW_scheme | DDT containing spectral information for RRTMGP SW radiation scheme | DDT | 0 | ty_gas_optics_rrtmgp | | in | F | -!! | kdist_cldy_lw | K_distribution_file_for_cloudy_RRTMGP_LW_scheme | DDT containing spectral information for cloudy RRTMGP LW radiation scheme | DDT | 0 | ty_cloud_optics | | in | F | -!! | kdist_cldy_sw | K_distribution_file_for_cloudy_RRTMGP_SW_scheme | DDT containing spectral information for cloudy RRTMGP SW radiation scheme | DDT | 0 | ty_cloud_optics | | in | F | -!! | optical_props_clouds | optical_properties_for_cloudy_atmosphere | Fortran DDT containing RRTMGP optical properties | DDT | 0 | ty_optical_props_1scl | | out | F | -!! | optical_props_aerosol | optical_properties_for_aerosols | Fortran DDT containing RRTMGP optical properties | DDT | 0 | ty_optical_props_1scl | | out | F | -!! | gas_concentrations | Gas_concentrations_for_RRTMGP_suite | DDT containing gas concentrations for RRTMGP radiation scheme | DDT | 0 | ty_gas_concs | | out | F | -!! | sfc_emiss_byband | surface_longwave_emissivity_in_each_band | surface lw emissivity in fraction in each LW band | frac | 2 | real | kind_phys | out | F | +!! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | +!! |-------------------------|---------------------------------------------------------------|-------------------------------------------------------------------------------|----------|------|-----------------------|-----------|--------|----------| +!! | Model | GFS_control_type_instance | Fortran DDT containing FV3-GFS model control parameters | DDT | 0 | GFS_control_type | | in | F | +!! | Grid | GFS_grid_type_instance | Fortran DDT containing FV3-GFS grid and interpolation related data | DDT | 0 | GFS_grid_type | | in | F | +!! | Sfcprop | GFS_sfcprop_type_instance | Fortran DDT containing FV3-GFS surface fields | DDT | 0 | GFS_sfcprop_type | | in | F | +!! | Statein | GFS_statein_type_instance | Fortran DDT containing FV3-GFS prognostic state data in from dycore | DDT | 0 | GFS_statein_type | | in | F | +!! | Tbd | GFS_tbd_type_instance | Fortran DDT containing FV3-GFS data not yet assigned to a defined container | DDT | 0 | GFS_tbd_type | | in | F | +!! | Coupling | GFS_coupling_type_instance | Fortran DDT containing FV3-GFS fields needed for coupling | DDT | 0 | GFS_coupling_type | | in | F | +!! | Radtend | GFS_radtend_type_instance | Fortran DDT containing FV3-GFS radiation tendencies | DDT | 0 | GFS_radtend_type | | inout | F | +!! | lm | vertical_layer_dimension_for_radiation | number of vertical layers for radiation calculation | count | 0 | integer | | in | F | +!! | im | horizontal_loop_extent | horizontal loop extent | count | 0 | integer | | in | F | +!! | lmk | adjusted_vertical_layer_dimension_for_radiation | number of vertical layers for radiation | count | 0 | integer | | in | F | +!! | lmp | adjusted_vertical_level_dimension_for_radiation | number of vertical levels for radiation | count | 0 | integer | | in | F | +!! | kd | vertical_index_difference_between_inout_and_local | vertical index difference between in/out and local | index | 0 | integer | | out | F | +!! | kt | vertical_index_difference_between_layer_and_upper_bound | vertical index difference between layer and upper bound | index | 0 | integer | | out | F | +!! | kb | vertical_index_difference_between_layer_and_lower_bound | vertical index difference between layer and lower bound | index | 0 | integer | | out | F | +!! | raddt | time_step_for_radiation | radiation time step | s | 0 | real | kind_phys | out | F | +!! | delp | layer_pressure_thickness_for_radiation | layer pressure thickness on radiation levels | hPa | 2 | real | kind_phys | out | F | +!! | dz | layer_thickness_for_radiation | layer thickness on radiation levels | km | 2 | real | kind_phys | out | F | +!! | plvl | air_pressure_at_interface_for_radiation_in_hPa | air pressure at vertical interface for radiation calculation | hPa | 2 | real | kind_phys | out | F | +!! | plyr | air_pressure_at_layer_for_radiation_in_hPa | air pressure at vertical layer for radiation calculation | hPa | 2 | real | kind_phys | out | F | +!! | tlvl | air_temperature_at_interface_for_radiation | air temperature at vertical interface for radiation calculation | K | 2 | real | kind_phys | out | F | +!! | tlyr | air_temperature_at_layer_for_radiation | air temperature at vertical layer for radiation calculation | K | 2 | real | kind_phys | out | F | +!! | tsfg | surface_ground_temperature_for_radiation | surface ground temperature for radiation | K | 1 | real | kind_phys | out | F | +!! | tsfa | surface_air_temperature_for_radiation | lowest model layer air temperature for radiation | K | 1 | real | kind_phys | out | F | +!! | qlyr | water_vapor_specific_humidity_at_layer_for_radiation | water vapor specific humidity at vertical layer for radiation calculation | kg kg-1 | 2 | real | kind_phys | out | F | +!! | olyr | ozone_concentration_at_layer_for_radiation | ozone concentration | kg kg-1 | 2 | real | kind_phys | out | F | +!! | icseed | seed_random_numbers_lw | seed for random number generation for longwave radiation | none | 1 | integer | | in | F | +!! | aerodp | atmosphere_optical_thickness_due_to_ambient_aerosol_particles | vertical integrated optical depth for various aerosol species | none | 2 | real | kind_phys | out | F | +!! | cldsa | cloud_area_fraction_for_radiation | fraction of clouds for low, middle,high, total and BL | frac | 2 | real | kind_phys | out | F | +!! | mtopa | model_layer_number_at_cloud_top | vertical indices for low, middle and high cloud tops | index | 2 | integer | | out | F | +!! | mbota | model_layer_number_at_cloud_base | vertical indices for low, middle and high cloud bases | index | 2 | integer | | out | F | +!! | de_lgth | cloud_decorrelation_length | cloud decorrelation length | km | 1 | real | kind_phys | out | F | +!! | alb1d | surface_albedo_perturbation | surface albedo perturbation | frac | 1 | real | kind_phys | out | F | +!! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | +!! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | +!! | kdist_lw | K_distribution_file_for_RRTMGP_LW_scheme | DDT containing spectral information for RRTMGP LW radiation scheme | DDT | 0 | ty_gas_optics_rrtmgp | | in | F | +!! | kdist_sw | K_distribution_file_for_RRTMGP_SW_scheme | DDT containing spectral information for RRTMGP SW radiation scheme | DDT | 0 | ty_gas_optics_rrtmgp | | in | F | +!! | kdist_cldy_lw | K_distribution_file_for_cloudy_RRTMGP_LW_scheme | DDT containing spectral information for cloudy RRTMGP LW radiation scheme | DDT | 0 | ty_cloud_optics | | in | F | +!! | kdist_cldy_sw | K_distribution_file_for_cloudy_RRTMGP_SW_scheme | DDT containing spectral information for cloudy RRTMGP SW radiation scheme | DDT | 0 | ty_cloud_optics | | in | F | +!! | optical_propsLW_clouds | longwave_optical_properties_for_cloudy_atmosphere | Fortran DDT containing RRTMGP optical properties | DDT | 0 | ty_optical_props_1scl | | out | F | +!! | optical_propsLW_aerosol | longwave_optical_properties_for_aerosols | Fortran DDT containing RRTMGP optical properties | DDT | 0 | ty_optical_props_1scl | | out | F | +!! | optical_propsSW_clouds | shortwave_optical_properties_for_cloudy_atmosphere | Fortran DDT containing RRTMGP optical properties | DDT | 0 | ty_optical_props_2str | | out | F | +!! | optical_propsSW_aerosol | shortwave_optical_properties_for_aerosols | Fortran DDT containing RRTMGP optical properties | DDT | 0 | ty_optical_props_2str | | out | F | +!! | gas_concentrations_lw | Gas_concentrations_for_RRTMGP_suite_lw | DDT containing gas concentrations for RRTMGP radiation scheme | DDT | 0 | ty_gas_concs | | out | F | +!! | gas_concentrations_sw | Gas_concentrations_for_RRTMGP_suite_sw | DDT containing gas concentrations for RRTMGP radiation scheme | DDT | 0 | ty_gas_concs | | out | F | +!! | sfc_emiss_byband | surface_longwave_emissivity_in_each_band | surface lw emissivity in fraction in each LW band | frac | 2 | real | kind_phys | out | F | +!! | sfc_alb_nir_dir | surface_shortwave_albedo_near_infrared_direct_in_each_band | surface sw near-infrared direct albedo in each SW band | frac | 2 | real | kind_phys | out | F | +!! | sfc_alb_nir_dif | surface_shortwave_albedo_near_infrared_diffuse_in_each_band | surface sw near-infrared diffuse albedo in each SW band | frac | 2 | real | kind_phys | out | F | +!! | sfc_alb_uvvis_dir | surface_shortwave_albedo_uv_visible_direct_in_each_band | surface sw uv-visible direct albedo in each SW band | frac | 2 | real | kind_phys | out | F | +!! | sfc_alb_uvvis_dif | surface_shortwave_albedo_uv_visible_diffuse_in_each_band | surface sw uv-visible diffuse albedo in each SW band | frac | 2 | real | kind_phys | out | F | +!! | nday | daytime_points_dimension | daytime points dimension | count | 0 | integer | | out | F | +!! | idxday | daytime_points | daytime points | index | 1 | integer | | out | F | !! - ! Attention - the output arguments lm, im, lmk, lmp must not be set - ! in the CCPP version - they are defined in the interstitial_create routine - ! ######################################################################################### - subroutine GFS_rrtmgp_pre_run (Model, Grid, Sfcprop, Statein, Tbd, Cldprop, Coupling, & ! IN - Radtend, & ! INOUT - lm, im, lmk, lmp, kdist_lw, kdist_sw, kdist_cldy_lw, kdist_cldy_sw, & ! IN - kd, kt, kb, raddt, delp, dz, plvl, plyr, tlvl, tlyr, tsfg, tsfa, qlyr, olyr, icseed, & ! OUT - gasvmr_co2, gasvmr_n2o, gasvmr_ch4, gasvmr_o2, gasvmr_co, gasvmr_cfc11, & ! OUT - gasvmr_cfc12, gasvmr_cfc22, gasvmr_ccl4, gasvmr_cfc113, faersw1, faersw2, faersw3, & ! OUT - faerlw1, faerlw2, faerlw3, aerodp, clouds1, clouds2, clouds3, clouds4, clouds5, & ! OUT - clouds6, clouds7, clouds8, clouds9, cldsa, mtopa, mbota, de_lgth, alb1d, & ! OUT - optical_props_clouds, optical_props_aerosol, gas_concentrations, sfc_emiss_byband, errmsg, errflg) - - use physparam - use machine, only: & - kind_phys ! Working type - use GFS_typedefs, only: & - GFS_statein_type, & ! Prognostic state data in from dycore - GFS_stateout_type, & ! Prognostic state or tendencies return to dycore - GFS_sfcprop_type, & ! Surface fields - GFS_coupling_type, & ! Fields to/from coupling with other components (e.g. land/ice/ocean/etc.) - GFS_control_type, & ! Model control parameters - GFS_grid_type, & ! Grid and interpolation related data - GFS_tbd_type, & ! To-Be-Determined data that doesn't fit in any one container - GFS_cldprop_type, & ! Cloud fields needed by radiation from physics - GFS_radtend_type, & ! Radiation tendencies needed in physics - GFS_diag_type ! Fields targetted for diagnostic output - use physcons, only: & - eps => con_eps, & ! Rd/Rv - epsm1 => con_epsm1, & ! Rd/Rv-1 - fvirt => con_fvirt, & ! Rv/Rd-1 - rog => con_rog, & ! Rd/g - rocp => con_rocp ! Rd/cp - use radcons, only: & - itsfc, & ! Flag for LW sfc. temp. - ltp, & ! 1-add extra-top layer; 0-no extra layer - lextop, & ! ltp > 0 - qmin,qme5, qme6, epsq ! Minimum vlaues for varius calculations - use funcphys, only: & - fpvs ! Function ot compute sat. vapor pressure over liq. - use module_radiation_astronomy,only: & - coszmn ! Function to compute cos(SZA) - use module_radiation_gases, only: & - NF_VGAS, & ! Number of active gas species - getgases, & ! Routine to setup trace gases - getozn ! Routine to setup ozone - use module_radiation_aerosols, only: & - NF_AESW, & ! Number of optical-fields in SW output (3=tau+g+omega) - NF_AELW, & ! Number of optical-fields in LW output (3=tau+g+omega) - setaer, & ! Routine to compute aerosol radiative properties (tau,g,omega) - NSPC1 ! Number of species for vertically integrated aerosol optical-depth - use module_radiation_clouds, only: & - NF_CLDS, & ! Number of fields in "clouds" array (e.g. (cloud(1)=lwp,clouds(2)=ReffLiq,...) - progcld1, & ! Zhao/Moorthi's prognostic cloud scheme - progcld3, & ! Zhao/Moorthi's prognostic cloud+pdfcld - progcld4, & ! GFDL cloud scheme - progcld5, & ! Thompson / WSM6 cloud micrphysics scheme - progclduni ! Unified cloud-scheme - use surface_perturbation, only: & - cdfnor ! Routine to compute CDF (used to compute percentiles) - use module_radiation_surface, only: & - setemis ! Routine to compute surface-emissivity - use rrtmgp_lw, only: check_error_msg, nrghice, ipsdlw0 - use mersenne_twister, only: & - random_setseed, & - random_number, & - random_stat - ! RRTMGP types - use mo_gas_optics_rrtmgp, only: ty_gas_optics_rrtmgp - use mo_cloud_optics, only: ty_cloud_optics - use mo_optical_props, only: ty_optical_props_1scl - use mo_cloud_sampling, only: sampled_mask_max_ran, sampled_mask_exp_ran, draw_samples - use mo_gas_concentrations, only: ty_gas_concs - - implicit none - - ! Inputs - type(GFS_control_type), intent(in) :: Model - type(GFS_grid_type), intent(in) :: Grid - type(GFS_sfcprop_type), intent(in) :: Sfcprop - type(GFS_statein_type), intent(in) :: Statein - type(GFS_radtend_type), intent(inout) :: Radtend - type(GFS_tbd_type), intent(in) :: Tbd - type(GFS_cldprop_type), intent(in) :: Cldprop - type(GFS_coupling_type), intent(in) :: Coupling - integer, intent(in) :: im, lm, lmk, lmp - type(ty_gas_optics_rrtmgp),intent(in) :: & - kdist_lw, & ! RRTMGP DDT containing spectral information for LW calculation - kdist_sw ! RRTMGP DDT containing spectral information for SW calculation - type(ty_cloud_optics),intent(in) :: & - kdist_cldy_lw, & - kdist_cldy_sw - type(ty_gas_concs),intent(out) :: & - gas_concentrations - integer,intent(in),dimension(IM) :: & - icseed ! auxiliary special cloud related array when module - ! variable isubclw=2, it provides permutation seed - ! for each column profile that are used for generating - ! random numbers. when isubclw /=2, it will not be used. - - ! Outputs - integer, intent(out) :: kd, kt, kb - real(kind_phys), intent(out) :: raddt - real(kind_phys), dimension(size(Grid%xlon,1)), intent(out) :: & - tsfg, tsfa - real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP), intent(out) :: & - delp, dz, plyr, tlyr, qlyr, olyr, gasvmr_co2, gasvmr_n2o, gasvmr_ch4, & - gasvmr_o2, gasvmr_co, gasvmr_cfc11, gasvmr_cfc12, gasvmr_cfc22, & - gasvmr_ccl4, gasvmr_cfc113, clouds1, clouds2, clouds3, clouds4, clouds5, & - clouds6, clouds7, clouds8, clouds9 - real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+1+LTP), intent(out) :: & - plvl, tlvl - real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,kdist_sw%get_nband()), intent(out) :: & - faersw1, faersw2, faersw3 - real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,kdist_lw%get_nband()), intent(out) :: & - faerlw1, faerlw2, faerlw3 - real(kind_phys), dimension(size(Grid%xlon,1),NSPC1), intent(out) :: & - aerodp - - real(kind_phys), dimension(size(Grid%xlon,1),5), intent(out) :: cldsa - integer, dimension(size(Grid%xlon,1),3), intent(out) :: mbota,mtopa - real(kind_phys), dimension(size(Grid%xlon,1)), intent(out) :: de_lgth,alb1d - character(len=*), intent(out) :: errmsg - integer, intent(out) :: errflg - type(ty_optical_props_1scl),intent(out) :: & - optical_props_clouds, & - optical_props_aerosol - real(kind_phys),dimension(kdist_lw%get_nband(),Model%levr+LTP),intent(out) :: sfc_emiss_byband - - ! Local variables - integer :: me, nfxr, ntrac, ntcw, ntiw, ncld, ntrw, ntsw, ntgl,i, j, k, k1, k2, lsk, & - lv, n, itop, ibtc, LP1, lla, llb, lya, lyb, iCol, iBand - integer,dimension(IM) :: ipseed - logical,dimension(IM,Model%levr+LTP) :: & - liqmask,icemask - real(kind_phys),dimension(IM,Model%levr+LTP) :: & - cld_ref_ice2,cld_ref_liq2, vmr_o3, vmr_h2o - real(kind_phys) :: es, qs, delt, tem0d - real(kind_phys), dimension(size(Grid%xlon,1)) :: cvt1, cvb1, tem1d, tskn - real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP) :: htswc, htlwc, & - gcice, grain, grime, htsw0, htlw0, rhly, tvly,qstl, vvel, clw, ciw, prslk1, & - tem2da, cldcov, deltaq, cnvc, cnvw, effrl, effri, effrr, effrs - real (kind_phys) :: clwmin, clwm, clwt, onemrh, value, tem1, tem2, tem3 - real (kind_phys), parameter :: xrc3 = 100. - real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP+1) :: tem2db - real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,Model%ncnd) :: ccnd - real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,2:Model%ntrac) :: tracer1 - real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,NF_CLDS) :: clouds - real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,NF_VGAS) :: gasvmr - real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,kdist_sw%get_nband(),NF_AESW)::faersw - real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,kdist_lw%get_nband(),NF_AELW)::faerlw - type(ty_optical_props_1scl) :: optical_props_clear, optical_props_cloudsByBand - real(kind_phys), dimension(kdist_lw%get_ngpt(),Model%levr+LTP,IM) :: & - rng3D - real(kind_phys), dimension(kdist_lw%get_ngpt()*(Model%levr+LTP)) :: & - rng1D - logical, dimension(IM,Model%levr+LTP,kdist_lw%get_ngpt()) :: & - cldfracMCICA - type(random_stat) :: rng_stat - - ! Initialize CCPP error handling variables - errmsg = '' - errflg = 0 - - if (.not. (Model%lsswr .or. Model%lslwr)) return - - ! Define some commonly used integers - me = Model%me ! MPI rank designator - NFXR = Model%nfxr ! second dimension for fluxr diagnostic variable (radiation) - NTRAC = Model%ntrac ! Number of tracers - ntcw = Model%ntcw ! Tracer index for cloud condensate (or liquid water) - ntiw = Model%ntiw ! Tracer index for ice - ncld = Model%ncld ! Cloud scheme - ntrw = Model%ntrw ! Tracer index for rain - ntsw = Model%ntsw ! Tracer index for snow - ntgl = Model%ntgl ! Tracer index for groupel - LP1 = LM + 1 ! num of in/out levels - - ! Set local /level/layer indexes corresponding to in/out variables - if ( lextop ) then - if ( ivflip == 1 ) then ! vertical from sfc upward - kd = 0 ! index diff between in/out and local - kt = 1 ! index diff between lyr and upper bound - kb = 0 ! index diff between lyr and lower bound - lla = LMK ! local index at the 2nd level from top - llb = LMP ! local index at toa level - lya = LM ! local index for the 2nd layer from top - lyb = LP1 ! local index for the top layer - else ! vertical from toa downward - kd = 1 ! index diff between in/out and local - kt = 0 ! index diff between lyr and upper bound - kb = 1 ! index diff between lyr and lower bound - lla = 2 ! local index at the 2nd level from top - llb = 1 ! local index at toa level - lya = 2 ! local index for the 2nd layer from top - lyb = 1 ! local index for the top layer - endif ! end if_ivflip_block - else - kd = 0 - if ( ivflip == 1 ) then ! vertical from sfc upward - kt = 1 ! index diff between lyr and upper bound - kb = 0 ! index diff between lyr and lower bound - else ! vertical from toa downward - kt = 0 ! index diff between lyr and upper bound - kb = 1 ! index diff between lyr and lower bound - endif ! end if_ivflip_block - endif ! end if_lextop_block - - ! Radiation time step (output) - raddt = min(Model%fhswr, Model%fhlwr) - - ! Setup surface ground temperature and ground/air skin temperature if required. - if ( itsfc == 0 ) then ! use same sfc skin-air/ground temp - tskn(1:IM) = Sfcprop%tsfc(1:IM) - tsfg(1:IM) = Sfcprop%tsfc(1:IM) - else ! use diff sfc skin-air/ground temp - tskn(1:IM) = Sfcprop%tsfc(1:IM) - tsfg(1:IM) = Sfcprop%tsfc(1:IM) - endif - - ! Prepare atmospheric profiles for radiation input. - lsk = 0 - if (ivflip == 0 .and. lm < Model%levs) lsk = Model%levs - lm - - ! Copy over state fields into fields, compute some needed quantities. - do k = 1, LM - k1 = k + kd - k2 = k + lsk - do i = 1, IM - plvl(i,k1+kb) = Statein%prsi(i,k2+kb) - plyr(i,k1) = Statein%prsl(i,k2) - tlyr(i,k1) = Statein%tgrs(i,k2) - prslk1(i,k1) = Statein%prslk(i,k2) - - ! Compute relative humidity. - es = min( Statein%prsl(i,k2), fpvs( Statein%tgrs(i,k2) ) ) ! fpvs and prsl in pa - qs = max( QMIN, eps * es / (Statein%prsl(i,k2) + epsm1*es) ) - rhly(i,k1) = max( 0.0, min( 1.0, max(QMIN, Statein%qgrs(i,k2,1))/qs ) ) - qstl(i,k1) = qs - enddo - plvl(i,k1+kb) = kdist_lw%get_press_min() - enddo - - ! Recast remaining all tracers (except sphum) forcing them all to be positive - do j = 2, NTRAC - do k = 1, LM - k1 = k + kd - k2 = k + lsk - tracer1(:,k1,j) = max(0.0, Statein%qgrs(:,k2,j)) - enddo - enddo - - ! Input data from toa to sfc - if (ivflip == 0) then - plvl(1:IM-1,1+kd) = Statein%prsi(1:IM-1,1) - plvl(IM,1+kd) = kdist_lw%get_press_min() - if (lsk /= 0) then - plvl(1:IM-1,1+kd) = 0.5 * (plvl(1:IM-1,2+kd) + plvl(1:IM-1,1+kd)) - plvl(IM,1+kd) = kdist_lw%get_press_min() - endif - ! Input data from sfc to top - else - plvl(1:IM-1,LP1+kd) = Statein%prsi(1:IM-1,LP1+lsk) - plvl(IM,LP1+kd) = kdist_lw%get_press_min() - if (lsk /= 0) then - plvl(1:IM-1,LM+kd) = 0.5 * (plvl(1:IM-1,LP1+kd) + plvl(1:IM-1,LM+kd)) - plvl(IM,LM+kd) = kdist_lw%get_press_min() - endif - endif - - if ( lextop ) then ! values for extra top layer - do i = 1, IM - plvl(i,llb) = kdist_lw%get_press_min() - if ( plvl(i,lla) <= kdist_lw%get_press_min() ) plvl(i,lla) = 2.0* kdist_lw%get_press_min() - plyr(i,lyb) = 0.5 * plvl(i,lla) - tlyr(i,lyb) = tlyr(i,lya) - prslk1(i,lyb) = (plyr(i,lyb)*0.001) ** rocp ! plyr in Pa - rhly(i,lyb) = rhly(i,lya) - qstl(i,lyb) = qstl(i,lya) - enddo - - ! note: may need to take care the top layer amount - tracer1(:,lyb,:) = tracer1(:,lya,:) - endif - - ! Get layer ozone mass mixing ratio - if (Model%ntoz > 0) then - do k=1,lmk - do i=1,im - olyr(i,k) = max( QMIN, tracer1(i,k,Model%ntoz) ) - enddo - enddo - ! Use climatological ozone data - else - call getozn (prslk1, Grid%xlat, IM, LMK, olyr) - endif - - ! Compute cosine of zenith angle (only when SW is called) - if (Model%lsswr) then - call coszmn (Grid%xlon, Grid%sinlat, Grid%coslat, Model%solhr, IM, me, & - Radtend%coszen, Radtend%coszdg) - endif - - ! Call getgases(), to set up non-prognostic gas volume mixing ratios (gasvmr). - call getgases (plvl/100., Grid%xlon, Grid%xlat, IM, LMK, gasvmr) - - ! Assign to gasvmr_XXXX - gasvmr_co2 (1:IM,1:LMK) = gasvmr(1:IM,1:LMK,1) - gasvmr_n2o (1:IM,1:LMK) = gasvmr(1:IM,1:LMK,2) - gasvmr_ch4 (1:IM,1:LMK) = gasvmr(1:IM,1:LMK,3) - gasvmr_o2 (1:IM,1:LMK) = gasvmr(1:IM,1:LMK,4) - gasvmr_co (1:IM,1:LMK) = gasvmr(1:IM,1:LMK,5) - gasvmr_cfc11 (1:IM,1:LMK) = gasvmr(1:IM,1:LMK,6) - gasvmr_cfc12 (1:IM,1:LMK) = gasvmr(1:IM,1:LMK,7) - gasvmr_cfc22 (1:IM,1:LMK) = gasvmr(1:IM,1:LMK,8) - gasvmr_ccl4 (1:IM,1:LMK) = gasvmr(1:IM,1:LMK,9) - gasvmr_cfc113 (1:IM,1:LMK) = gasvmr(1:IM,1:LMK,10) - - ! Get temperature at layer interface, and layer moisture. - tem2da(1:IM,2:LMK) = log( plyr(1:IM,2:LMK) ) - tem2db(1:IM,2:LMK) = log( plvl(1:IM,2:LMK) ) - - if (ivflip == 0) then ! input data from toa to sfc - do i = 1, IM - tem1d (i) = QME6 - tem2da(i,1) = log( plyr(i,1) ) - tem2db(i,1) = log( max( kdist_lw%get_press_min(), plvl(i,1)) ) - tem2db(i,LMP) = log( plvl(i,LMP) ) - tsfa (i) = tlyr(i,LMK) ! sfc layer air temp - tlvl(i,1) = tlyr(i,1) - tlvl(i,LMP) = tskn(i) - enddo - do k = 1, LM - k1 = k + kd - do i = 1, IM - qlyr(i,k1) = max( tem1d(i), Statein%qgrs(i,k,1) ) - tem1d(i) = min( QME5, qlyr(i,k1) ) - tvly(i,k1) = Statein%tgrs(i,k) * (1.0 + fvirt*qlyr(i,k1)) ! virtual T (K) - delp(i,k1) = plvl(i,k1+1) - plvl(i,k1) - enddo - enddo - - if ( lextop ) then - do i = 1, IM - qlyr(i,lyb) = qlyr(i,lya) - tvly(i,lyb) = tvly(i,lya) - delp(i,lyb) = plvl(i,lla) - plvl(i,llb) - enddo - endif - - do k = 2, LMK - do i = 1, IM - tlvl(i,k) = tlyr(i,k) + (tlyr(i,k-1) - tlyr(i,k)) * (tem2db(i,k) - tem2da(i,k)) / & - (tem2da(i,k-1) - tem2da(i,k)) - enddo - enddo - - ! Comput lvel height and layer thickness (km) - tem0d = 0.001 * rog - do i = 1, IM - do k = 1, LMK - dz(i,k) = tem0d * (tem2db(i,k+1) - tem2db(i,k)) * tvly(i,k) - enddo - enddo - - else ! input data from sfc to toa - do i = 1, IM - tem1d (i) = QME6 - tem2da(i,1) = log( plyr(i,1) ) - tem2db(i,1) = log( plvl(i,1) ) - tem2db(i,LMP) = log( max( kdist_lw%get_press_min(), plvl(i,LMP)) ) - tsfa (i) = tlyr(i,1) ! sfc layer air temp - tlvl(i,1) = tskn(i) - tlvl(i,LMP) = tlyr(i,LMK) - enddo - - do k = LM, 1, -1 - do i = 1, IM - qlyr(i,k) = max( tem1d(i), Statein%qgrs(i,k,1) ) - tem1d(i) = min( QME5, qlyr(i,k) ) - tvly(i,k) = Statein%tgrs(i,k) * (1.0 + fvirt*qlyr(i,k)) ! virtual T (K) - delp(i,k) = plvl(i,k) - plvl(i,k+1) - enddo - enddo - - if ( lextop ) then - do i = 1, IM - qlyr(i,lyb) = qlyr(i,lya) - tvly(i,lyb) = tvly(i,lya) - delp(i,lyb) = plvl(i,lla) - plvl(i,llb) - enddo - endif - - do k = 1, LMK-1 - do i = 1, IM - tlvl(i,k+1) = tlyr(i,k) + (tlyr(i,k+1) - tlyr(i,k)) * (tem2db(i,k+1) - tem2da(i,k)) / & - (tem2da(i,k+1) - tem2da(i,k)) - enddo - enddo - - ! Compute level height and layer thickness (km) - tem0d = 0.001 * rog - do i = 1, IM - do k = LMK, 1, -1 - dz(i,k) = tem0d * (tem2db(i,k) - tem2db(i,k+1)) * tvly(i,k) - enddo - enddo - endif ! end_if_ivflip - - ! Call module_radiation_aerosols::setaer(),to setup aerosols property profile for radiation. - call setaer (plvl, plyr, prslk1, tvly, rhly, Sfcprop%slmsk, tracer1, Grid%xlon, & - Grid%xlat, IM, LMK, LMP, Model%lsswr, Model%lslwr, faersw, faerlw, aerodp) - - ! Store aerosol optical properties - ! SW. - ! For RRTMGP SW the bands are now ordered from [IR(band) -> nIR -> UV], in RRTMG the - ! band ordering was [nIR -> UV -> IR(band)] - faersw1(1:IM,1:LMK,1) = faersw(1:IM,1:LMK,kdist_sw%get_nband(),1) - faersw2(1:IM,1:LMK,1) = faersw(1:IM,1:LMK,kdist_sw%get_nband(),2) - faersw3(1:IM,1:LMK,1) = faersw(1:IM,1:LMK,kdist_sw%get_nband(),3) - faersw1(1:IM,1:LMK,2:kdist_sw%get_nband()) = faersw(1:IM,1:LMK,1:kdist_sw%get_nband()-1,1) - faersw2(1:IM,1:LMK,2:kdist_sw%get_nband()) = faersw(1:IM,1:LMK,1:kdist_sw%get_nband()-1,2) - faersw3(1:IM,1:LMK,2:kdist_sw%get_nband()) = faersw(1:IM,1:LMK,1:kdist_sw%get_nband()-1,3) - ! LW - faerlw1(1:IM,1:LMK,1:kdist_lw%get_nband()) = faerlw(1:IM,1:LMK,1:kdist_lw%get_nband(),1) - faerlw2(1:IM,1:LMK,1:kdist_lw%get_nband()) = faerlw(1:IM,1:LMK,1:kdist_lw%get_nband(),2) - faerlw3(1:IM,1:LMK,1:kdist_lw%get_nband()) = faerlw(1:IM,1:LMK,1:kdist_lw%get_nband(),3) - - ! Obtain cloud information for radiation calculations - ! (clouds,cldsa,mtopa,mbota) - ! for prognostic cloud: - ! - For Zhao/Moorthi's prognostic cloud scheme, - ! call module_radiation_clouds::progcld1() - ! - For Zhao/Moorthi's prognostic cloud+pdfcld, - ! call module_radiation_clouds::progcld3() - ! call module_radiation_clouds::progclduni() for unified cloud and ncld=2 - ccnd = 0.0_kind_phys - if (Model%ncnd == 1) then ! Zhao_Carr_Sundqvist - ccnd(1:IM,1:LMK,1) = tracer1(1:IM,1:LMK,ntcw) ! -liquid water/ice - elseif (Model%ncnd == 2) then ! MG - ccnd(1:IM,1:LMK,1) = tracer1(1:IM,1:LMK,ntcw) ! -liquid water - ccnd(1:IM,1:LMK,2) = tracer1(1:IM,1:LMK,ntiw) ! -ice water - elseif (Model%ncnd == 4) then ! MG2 - ccnd(1:IM,1:LMK,1) = tracer1(1:IM,1:LMK,ntcw) ! -liquid water - ccnd(1:IM,1:LMK,2) = tracer1(1:IM,1:LMK,ntiw) ! -ice water - ccnd(1:IM,1:LMK,3) = tracer1(1:IM,1:LMK,ntrw) ! -rain water - ccnd(1:IM,1:LMK,4) = tracer1(1:IM,1:LMK,ntsw) ! -snow water - elseif (Model%ncnd == 5) then ! GFDL MP, Thompson, MG3 - ccnd(1:IM,1:LMK,1) = tracer1(1:IM,1:LMK,ntcw) ! -liquid water - ccnd(1:IM,1:LMK,2) = tracer1(1:IM,1:LMK,ntiw) ! -ice water - ccnd(1:IM,1:LMK,3) = tracer1(1:IM,1:LMK,ntrw) ! -rain water - ccnd(1:IM,1:LMK,4) = tracer1(1:IM,1:LMK,ntsw) + & ! -snow + grapuel - tracer1(1:IM,1:LMK,ntgl) - endif - where(ccnd < epsq) ccnd = 0.0 - - if (Model%imp_physics == 11 ) then - if (.not. Model%lgfdlmprad) then - ccnd(:,:,1) = tracer1(:,1:LMK,ntcw) - ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:LMK,ntrw) - ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:LMK,ntiw) - ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:LMK,ntsw) - ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:LMK,ntgl) - endif - do k=1,LMK - do i=1,IM - if (ccnd(i,k,1) < EPSQ ) ccnd(i,k,1) = 0.0 - enddo - enddo - endif - - ! Add suspended convective cloud water to grid-scale cloud water - ! only for cloud fraction & radiation computation it is to enhance - ! cloudiness due to suspended convec cloud water for zhao/moorthi's - ! (imp_phys=99) & ferrier's (imp_phys=5) microphysics schemes - if ((Model%num_p3d == 4) .and. (Model%npdf3d == 3)) then ! same as Model%imp_physics = 99 - deltaq(1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,5) - cnvw (1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,6) - cnvc (1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,7) - elseif ((Model%npdf3d == 0) .and. (Model%ncnvcld3d == 1)) then ! same as MOdel%imp_physics=98 - deltaq(1:im,1+kd:lm+kd) = 0.0 - cnvw (1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,Model%num_p3d+1) - cnvc (1:im,1+kd:lm+kd) = 0.0 - else ! all the rest - deltaq(1:im,1:lmk) = 0.0 - cnvw (1:im,1:lmk) = 0.0 - cnvc (1:im,1:lmk) = 0.0 - endif - - if (lextop) then - cldcov(1:im,lyb) = cldcov(1:im,lya) - deltaq(1:im,lyb) = deltaq(1:im,lya) - cnvw (1:im,lyb) = cnvw (1:im,lya) - cnvc (1:im,lyb) = cnvc (1:im,lya) - if (Model%effr_in) then - effrl(1:im,lyb) = effrl(1:im,lya) - effri(1:im,lyb) = effri(1:im,lya) - effrr(1:im,lyb) = effrr(1:im,lya) - effrs(1:im,lyb) = effrs(1:im,lya) - endif - endif - - if (Model%imp_physics == 99) then - ccnd(1:IM,1:LMK,1) = ccnd(1:IM,1:LMK,1) + cnvw(1:IM,1:LMK) - endif - - if (Model%imp_physics == 10) then - ccnd(1:IM,1:LMK,1) = ccnd(1:IM,1:LMK,1) + cnvw(1:IM,1:LMK) + ccnd(1:IM,1:LMK,2) - endif - - ! DJS2019: START - ! Compute layer cloud fraction. - clwmin = 0.0 - cldcov(:,:) = 0.0 - if (.not. Model%lmfshal) then - do k = 1, LMK - do i = 1, IM - clwt = 1.0e-6 * (plyr(i,k)*0.1) - if (ccnd(i,k,1) > 0.) then - onemrh= max( 1.e-10, 1.0-rhly(i,k) ) - clwm = clwmin / max( 0.01, plyr(i,k)*0.1 ) - tem1 = min(max(sqrt(sqrt(onemrh*qstl(i,k))),0.0001),1.0) - tem1 = 2000.0 / tem1 - value = max( min( tem1*(ccnd(i,k,1)-clwm), 50.0 ), 0.0 ) - tem2 = sqrt( sqrt(rhly(i,k)) ) - cldcov(i,k) = max( tem2*(1.0-exp(-value)), 0.0 ) - endif - enddo - enddo - else - do k = 1, LMK - do i = 1, IM - clwt = 1.0e-6 * (plyr(i,k)*0.1) - if (ccnd(i,k,1) .gt. 0) then - onemrh= max( 1.e-10, 1.0-rhly(i,k) ) - clwm = clwmin / max( 0.01, plyr(i,k)*0.1 ) - tem1 = min(max((onemrh*qstl(i,k))**0.49,0.0001),1.0) !jhan - if (Model%lmfdeep2) then - tem1 = xrc3 / tem1 - else - tem1 = 100.0 / tem1 - endif - value = max( min( tem1*(ccnd(i,k,1)-clwm), 50.0 ), 0.0 ) - tem2 = sqrt( sqrt(rhly(i,k)) ) - cldcov(i,k) = max( tem2*(1.0-exp(-value)), 0.0 ) - endif - enddo - enddo - endif - ! DJS2019: END - - if (Model%uni_cld) then - if (Model%effr_in) then - cldcov(1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,Model%indcld) - effrl(1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,2) - effri(1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,3) - effrr(1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,4) - effrs(1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,5) - else - do k=1,lm - k1 = k + kd - do i=1,im - !cldcov(i,k1) = Tbd%phy_f3d(i,k,Model%indcld) - !if (tracer1(i,k,ntcw) .gt. 0 .or. tracer1(i,k,ntiw) .gt. 0) then - ! cldcov(i,k1) = 0.1 - !else - ! cldcov(i,k1) = 0.0 - !endif - enddo - enddo - endif - elseif (Model%imp_physics == Model%imp_physics_gfdl) then ! GFDL MP - cldcov(1:IM,1+kd:LM+kd) = tracer1(1:IM,1:LM,Model%ntclamt) - else ! neither of the other two cases - ! cldcov = 0.0 - endif - - ! MICROPHYSICS - ! *) zhao/moorthi's prognostic cloud scheme or unified cloud and/or with MG microphysics - if (Model%imp_physics == 99 .or. Model%imp_physics == 10) then - if (Model%uni_cld .and. Model%ncld >= 2) then - call progclduni (plyr/100., plvl/100., tlyr, tvly, ccnd, Model%ncnd, & ! IN - Grid%xlat, Grid%xlon, Sfcprop%slmsk, dz, delp/100.,IM, & ! IN - LMK, LMP, cldcov, effrl, effri, effrr, effrs, Model%effr_in, & ! IN - clouds, cldsa, mtopa, mbota, de_lgth) ! OUT - else - call progcld1 (plyr/100. ,plvl/100., tlyr, tvly, qlyr, qstl, rhly, & ! IN - ccnd(1:IM,1:LMK,1), Grid%xlat,Grid%xlon,Sfcprop%slmsk, dz, & ! IN - delp/100., IM, LMK, LMP, Model%uni_cld, Model%lmfshal, & ! IN - Model%lmfdeep2, cldcov, effrl, effri, effrr, effrs, & ! IN - Model%effr_in, & ! IN - clouds, cldsa, mtopa, mbota, de_lgth) ! OUT - endif - ! *) zhao/moorthi's prognostic cloud+pdfcld - elseif(Model%imp_physics == 98) then - call progcld3 (plyr/100., plvl/100., tlyr, tvly, qlyr, qstl, rhly, & ! IN - ccnd(1:IM,1:LMK,1), cnvw, cnvc, Grid%xlat, Grid%xlon, & ! IN - Sfcprop%slmsk, dz, delp/100., im, lmk, lmp, deltaq, Model%sup, & ! IN - Model%kdt, me, & ! IN - clouds, cldsa, mtopa, mbota, de_lgth) ! OUT - ! *) GFDL cloud scheme - elseif (Model%imp_physics == 11) then - if (.not.Model%lgfdlmprad) then - call progcld4 (plyr/100., plvl/100., tlyr, tvly, qlyr, qstl, rhly, & ! IN - ccnd(1:IM,1:LMK,1), cnvw, cnvc,Grid%xlat, Grid%xlon, & ! IN - Sfcprop%slmsk, cldcov, dz, delp/100., im, lmk, lmp, & ! IN - clouds, cldsa, mtopa, mbota, de_lgth) ! OUT - else - call progclduni (plyr/100., plvl/100., tlyr, tvly, ccnd, Model%ncnd, & ! IN - Grid%xlat, Grid%xlon, Sfcprop%slmsk, dz,delp/100., IM, LMK, & ! IN - LMP, cldcov, effrl, effri, effrr, effrs, Model%effr_in, & ! IN - clouds, cldsa, mtopa, mbota, de_lgth) ! OUT - endif - ! *) Thompson / WSM6 cloud micrphysics scheme - elseif(Model%imp_physics == 8 .or. Model%imp_physics == 6) then - if (Model%kdt == 1) then - Tbd%phy_f3d(:,:,1) = 10. - Tbd%phy_f3d(:,:,2) = 50. - Tbd%phy_f3d(:,:,3) = 250. - endif - - call progcld5 (plyr/100., plvl/100., tlyr, qlyr, qstl, rhly, tracer1, Grid%xlat, & ! IN - Grid%xlon,Sfcprop%slmsk,dz,delp/100., ntrac-1, ntcw-1, ntiw-1, & ! IN - ntrw-1, ntsw-1, ntgl-1, im, lmk, lmp, Model%uni_cld, Model%lmfshal,& ! IN - Model%lmfdeep2, cldcov(:,1:LMK),Tbd%phy_f3d(:,:,1), & ! IN - Tbd%phy_f3d(:,:,2), Tbd%phy_f3d(:,:,3), & ! IN - clouds,cldsa,mtopa,mbota, de_lgth) ! OUT - endif ! end if_imp_physics - - ! CCPP - clouds1(1:IM,1:LMK) = clouds(1:IM,1:LMK,1) - clouds2(1:IM,1:LMK) = clouds(1:IM,1:LMK,2) - clouds3(1:IM,1:LMK) = clouds(1:IM,1:LMK,3) - clouds4(1:IM,1:LMK) = clouds(1:IM,1:LMK,4) - clouds5(1:IM,1:LMK) = clouds(1:IM,1:LMK,5) - clouds6(1:IM,1:LMK) = clouds(1:IM,1:LMK,6) - clouds7(1:IM,1:LMK) = clouds(1:IM,1:LMK,7) - clouds8(1:IM,1:LMK) = clouds(1:IM,1:LMK,8) - clouds9(1:IM,1:LMK) = clouds(1:IM,1:LMK,9) - - ! mg, sfc-perts - ! --- scale random patterns for surface perturbations with - ! perturbation size - ! --- turn vegetation fraction pattern into percentile pattern - alb1d(:) = 0. - if (Model%do_sfcperts) then - if (Model%pertalb(1) > 0.) then - do i=1,im - call cdfnor(Coupling%sfc_wts(i,5),alb1d(i)) - enddo - endif - endif - ! mg, sfc-perts - - ! ####################################################################################### - ! Call module_radiation_surface::setemis(),to setup surface emissivity for LW radiation. - ! ####################################################################################### - if (Model%lslwr) then - call setemis (Grid%xlon, Grid%xlat, Sfcprop%slmsk, Sfcprop%snowd, Sfcprop%sncovr, & - Sfcprop%zorl, tsfg, tsfa, Sfcprop%hprim, IM, Radtend%semis) - do iBand=1,kdist_lw%get_nband() - sfc_emiss_byband(iBand,1:IM) = Radtend%semis(1:IM) - enddo - endif - - ! ####################################################################################### - ! Compute radiative properties needed for RRTMGP - ! ####################################################################################### - - ! Change random number seed value for each radiation invocation (isubclw =1 or 2). - if(isubclw == 1) then ! advance prescribed permutation seed - do iCol = 1, IM - ipseed(iCol) = ipsdlw0 + iCol - enddo - elseif (isubclw == 2) then ! use input array of permutaion seeds - do iCol = 1, IM - ipseed(iCol) = icseed(iCol) - enddo - endif - - ! Compute volume mixing-ratios for ozone (mmr) and specific-humidity. - vmr_h2o = merge((qlyr/(1-qlyr))*amdw, 0., qlyr .ne. 1.) - vmr_o3 = merge(olyr*amdo3, 0., olyr .gt. 0.) - - ! Compute ice/liquid cloud masks, needed by rrtmgp_cloud_optics - liqmask = (clouds1 .gt. 0 .and. clouds2 .gt. 0) - icemask = (clouds1 .gt. 0 .and. clouds4 .gt. 0) - - ! RRTMGP cloud_optics expects particle size to be in a certain range. bound here - cld_ref_ice2 = clouds5 - where(cld_ref_ice2 .gt. kdist_cldy_lw%get_max_radius_ice()) cld_ref_ice2=kdist_cldy_lw%get_max_radius_ice() - where(cld_ref_ice2 .lt. kdist_cldy_lw%get_min_radius_ice()) cld_ref_ice2=kdist_cldy_lw%get_min_radius_ice() - cld_ref_liq2 = clouds5 - where(cld_ref_liq2 .gt. kdist_cldy_lw%get_max_radius_liq()) cld_ref_liq2=kdist_cldy_lw%get_max_radius_liq() - where(cld_ref_liq2 .lt. kdist_cldy_lw%get_min_radius_liq()) cld_ref_liq2=kdist_cldy_lw%get_min_radius_liq() - - ! Allocate space for gas optical properties [ncol,nlay,ngpts] - ! Cloud optics [nCol,nLay,nBands] - print*,'In GFS_rrtmgp_pre: ' - call check_error_msg(optical_props_cloudsByBand%init(kdist_lw%get_band_lims_wavenumber())) - call check_error_msg(optical_props_cloudsByBand%alloc_1scl(IM, LMK)) - ! Aerosol optics [Ccol,nLay,nBands] - call check_error_msg(optical_props_aerosol%init(kdist_lw%get_band_lims_wavenumber())) - call check_error_msg(optical_props_aerosol%alloc_1scl(IM, LMK)) - ! Cloud optics [nCol,nLay,nGpts] - call check_error_msg(optical_props_clouds%alloc_1scl(IM, LMK, kdist_lw)) - - ! Set gas concentrations - call gas_concentrations%reset() - call check_error_msg(gas_concentrations%set_vmr('o2', gasvmr_o2)) - call check_error_msg(gas_concentrations%set_vmr('co2', gasvmr_co2)) - call check_error_msg(gas_concentrations%set_vmr('ch4', gasvmr_ch4)) - call check_error_msg(gas_concentrations%set_vmr('n2o', gasvmr_n2o)) - call check_error_msg(gas_concentrations%set_vmr('h2o', vmr_h2o)) - call check_error_msg(gas_concentrations%set_vmr('o3', vmr_o3)) - - ! Copy aerosol to RRTMGP DDT - optical_props_aerosol%tau = faerlw1 * (1. - faerlw2) - - ! Compute cloud-optics for RTE. - call check_error_msg(kdist_cldy_lw%cloud_optics(IM, LMK, kdist_lw%get_nband(), nrghice, & - liqmask, icemask, clouds2, clouds4, clouds3, clouds5, optical_props_cloudsByBand)) - - ! Call McICA to generate subcolumns. - if (isubclw .gt. 0) then - - ! Call RNG. Mersennse Twister accepts 1D array, so loop over columns and collapse along G-points - ! and layers. ([nGpts,nLayer,nColumn]-> [nGpts*nLayer]*nColumn) - do iCol=1,IM - call random_setseed(ipseed(icol),rng_stat) - call random_number(rng1D,rng_stat) - rng3D(:,:,iCol) = reshape(source = rng1D,shape=[kdist_lw%get_ngpt(),LMK]) - enddo - - ! Call McICA - select case ( iovrlw ) - ! Maximumn-random - case(1) - call check_error_msg(sampled_mask_max_ran(rng3D,clouds1,cldfracMCICA)) - end select - - ! Map band optical depth to each g-point using McICA - call check_error_msg(draw_samples(cldfracMCICA,optical_props_cloudsByBand,optical_props_clouds)) - endif + ! Attention - the output arguments lm, im, lmk, lmp must not be set + ! in the CCPP version - they are defined in the interstitial_create routine + ! ######################################################################################### + subroutine GFS_rrtmgp_pre_run (Model, Grid, Sfcprop, Statein, Tbd, Coupling, & ! IN + Radtend, & ! INOUT + lm, im, lmk, lmp, kdist_lw, kdist_sw, kdist_cldy_lw, kdist_cldy_sw, & ! IN + kd, kt, kb, raddt, delp, dz, plvl, plyr, tlvl, tlyr, tsfg, tsfa, qlyr, olyr, icseed, & ! OUT + aerodp, cldsa, mtopa, mbota, de_lgth, alb1d, & ! OUT + optical_propsLW_clouds, optical_propsLW_aerosol, optical_propsSW_clouds, & ! OUT + optical_propsSW_aerosol, gas_concentrations_lw, gas_concentrations_sw, & + sfc_emiss_byband, sfc_alb_nir_dir, sfc_alb_nir_dif, sfc_alb_uvvis_dir, & + sfc_alb_uvvis_dif, nday, idxday, errmsg, errflg) + + ! Inputs + type(GFS_control_type), intent(in) :: Model + type(GFS_grid_type), intent(in) :: Grid + type(GFS_sfcprop_type), intent(in) :: Sfcprop + type(GFS_statein_type), intent(in) :: Statein + type(GFS_radtend_type), intent(inout) :: Radtend + type(GFS_tbd_type), intent(in) :: Tbd + type(GFS_coupling_type), intent(in) :: Coupling + integer, intent(in) :: im, lm, lmk, lmp + type(ty_gas_optics_rrtmgp),intent(in) :: & + kdist_lw, & ! RRTMGP DDT containing spectral information for LW calculation + kdist_sw ! RRTMGP DDT containing spectral information for SW calculation + type(ty_cloud_optics),intent(in) :: & + kdist_cldy_lw, & + kdist_cldy_sw + type(ty_gas_concs),intent(out) :: & + gas_concentrations_lw,gas_concentrations_sw + integer,intent(in),dimension(IM) :: & + icseed ! auxiliary special cloud related array when module + ! variable isubclw=2, it provides permutation seed + ! for each column profile that are used for generating + ! random numbers. when isubclw /=2, it will not be used. + + ! Outputs + integer, intent(out) :: kd, kt, kb + real(kind_phys), intent(out) :: raddt + real(kind_phys), dimension(size(Grid%xlon,1)), intent(out) :: & + tsfg, tsfa + real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP), intent(out) :: & + delp, dz, plyr, tlyr, qlyr, olyr + real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+1+LTP), intent(out) :: & + plvl, tlvl + real(kind_phys), dimension(size(Grid%xlon,1),NSPC1), intent(out) :: & + aerodp + + real(kind_phys), dimension(size(Grid%xlon,1),5), intent(out) :: cldsa + integer, dimension(size(Grid%xlon,1),3), intent(out) :: mbota,mtopa + real(kind_phys), dimension(size(Grid%xlon,1)), intent(out) :: de_lgth,alb1d + character(len=*), intent(out) :: errmsg + integer, intent(out) :: errflg + type(ty_optical_props_1scl),intent(out) :: & + optical_propsLW_clouds, & + optical_propsLW_aerosol + type(ty_optical_props_2str),intent(out) :: & + optical_propsSW_clouds, & + optical_propsSW_aerosol + real(kind_phys),dimension(kdist_sw%get_nband(),size(Grid%xlon,1)),intent(out) :: & + sfc_emiss_byband, & ! Longwave surface emissivity in each band + sfc_alb_nir_dir, & ! Shortwave surface albedo (nIR-direct) + sfc_alb_nir_dif, & ! Shortwave surface albedo (nIR-diffuse) + sfc_alb_uvvis_dir, & ! Shortwave surface albedo (uvvis-direct) + sfc_alb_uvvis_dif ! Shortwave surface albedo (uvvis-diffuse) + + integer, intent(out) :: nday + integer, dimension(size(Grid%xlon,1)), intent(out) :: idxday + + ! Local variables + integer :: me, nfxr, ntrac, ntcw, ntiw, ncld, ntrw, ntsw, ntgl,i, j, k, k1, k2, lsk, & + LP1, lla, llb, lya, lyb, iCol, iBand + integer,dimension(IM) :: ipseed_lw,ipseed_sw + logical,dimension(IM,Model%levr+LTP) :: & + liqmask,icemask + real(kind_phys),dimension(IM,Model%levr+LTP) :: & + vmr_o3, vmr_h2o + real(kind_phys) :: es, qs, tem0d + real(kind_phys), dimension(size(Grid%xlon,1)) :: tem1d, tskn + real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP) :: & + rhly, tvly, qstl, prslk1, & + tem2da, cldcov, deltaq, cnvc, cnvw, effrl, effri, effrr, effrs + real (kind_phys) :: clwmin, clwm, clwt, onemrh, value, tem1, tem2 + real (kind_phys), parameter :: xrc3 = 100. + real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP+1) :: tem2db + real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,Model%ncnd) :: ccnd + real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,2:Model%ntrac) :: tracer1 + real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,NF_CLDS) :: clouds + real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,NF_VGAS) :: gasvmr + real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,kdist_sw%get_nband(),NF_AESW)::faersw,faersw2 + real(kind_phys), dimension(size(Grid%xlon,1),Model%levr+LTP,kdist_lw%get_nband(),NF_AELW)::faerlw + type(ty_optical_props_1scl) :: optical_propsLW_cloudsByBand + type(ty_optical_props_2str) :: optical_propsSW_cloudsByBand + real(kind_phys), dimension(kdist_lw%get_ngpt(),Model%levr+LTP,IM) :: & + rng3D_lw + real(kind_phys), dimension(kdist_lw%get_ngpt()*(Model%levr+LTP)) :: & + rng1D_lw + logical, dimension(IM,Model%levr+LTP,kdist_lw%get_ngpt()) :: & + cldfracMCICA_lw + real(kind_phys), dimension(kdist_sw%get_ngpt(),Model%levr+LTP,IM) :: & + rng3D_sw + real(kind_phys), dimension(kdist_sw%get_ngpt()*(Model%levr+LTP)) :: & + rng1D_sw + logical, dimension(IM,Model%levr+LTP,kdist_sw%get_ngpt()) :: & + cldfracMCICA_sw + type(random_stat) :: rng_stat + real(kind_phys), dimension(size(Grid%xlon,1),NF_ALBD) :: sfcalb + + + ! Initialize CCPP error handling variables + errmsg = '' + errflg = 0 + + if (.not. (Model%lsswr .or. Model%lslwr)) return + + ! Define some commonly used integers + me = Model%me ! MPI rank designator + NFXR = Model%nfxr ! second dimension for fluxr diagnostic variable (radiation) + NTRAC = Model%ntrac ! Number of tracers + ntcw = Model%ntcw ! Tracer index for cloud condensate (or liquid water) + ntiw = Model%ntiw ! Tracer index for ice + ncld = Model%ncld ! Cloud scheme + ntrw = Model%ntrw ! Tracer index for rain + ntsw = Model%ntsw ! Tracer index for snow + ntgl = Model%ntgl ! Tracer index for groupel + LP1 = LM + 1 ! num of in/out levels + + ! Set local /level/layer indexes corresponding to in/out variables + if ( lextop ) then + if ( ivflip == 1 ) then ! vertical from sfc upward + kd = 0 ! index diff between in/out and local + kt = 1 ! index diff between lyr and upper bound + kb = 0 ! index diff between lyr and lower bound + lla = LMK ! local index at the 2nd level from top + llb = LMP ! local index at toa level + lya = LM ! local index for the 2nd layer from top + lyb = LP1 ! local index for the top layer + else ! vertical from toa downward + kd = 1 ! index diff between in/out and local + kt = 0 ! index diff between lyr and upper bound + kb = 1 ! index diff between lyr and lower bound + lla = 2 ! local index at the 2nd level from top + llb = 1 ! local index at toa level + lya = 2 ! local index for the 2nd layer from top + lyb = 1 ! local index for the top layer + endif ! end if_ivflip_block + else + kd = 0 + if ( ivflip == 1 ) then ! vertical from sfc upward + kt = 1 ! index diff between lyr and upper bound + kb = 0 ! index diff between lyr and lower bound + else ! vertical from toa downward + kt = 0 ! index diff between lyr and upper bound + kb = 1 ! index diff between lyr and lower bound + endif ! end if_ivflip_block + endif ! end if_lextop_block + + ! Radiation time step (output) + raddt = min(Model%fhswr, Model%fhlwr) + + ! Setup surface ground temperature and ground/air skin temperature if required. + if ( itsfc == 0 ) then ! use same sfc skin-air/ground temp + tskn(1:IM) = Sfcprop%tsfc(1:IM) + tsfg(1:IM) = Sfcprop%tsfc(1:IM) + else ! use diff sfc skin-air/ground temp + tskn(1:IM) = Sfcprop%tsfc(1:IM) + tsfg(1:IM) = Sfcprop%tsfc(1:IM) + endif + + ! Prepare atmospheric profiles for radiation input. + lsk = 0 + if (ivflip == 0 .and. lm < Model%levs) lsk = Model%levs - lm + + ! Copy over state fields into fields, compute some needed quantities. + do k = 1, LM + k1 = k + kd + k2 = k + lsk + do i = 1, IM + plvl(i,k1+kb) = Statein%prsi(i,k2+kb) + plyr(i,k1) = Statein%prsl(i,k2) + tlyr(i,k1) = Statein%tgrs(i,k2) + prslk1(i,k1) = Statein%prslk(i,k2) + + ! Compute relative humidity. + es = min( Statein%prsl(i,k2), fpvs( Statein%tgrs(i,k2) ) ) ! fpvs and prsl in pa + qs = max( QMIN, eps * es / (Statein%prsl(i,k2) + epsm1*es) ) + rhly(i,k1) = max( 0.0, min( 1.0, max(QMIN, Statein%qgrs(i,k2,1))/qs ) ) + qstl(i,k1) = qs + enddo + plvl(i,k1+kb) = kdist_lw%get_press_min() + enddo + + ! Recast remaining all tracers (except sphum) forcing them all to be positive + do j = 2, NTRAC + do k = 1, LM + k1 = k + kd + k2 = k + lsk + tracer1(:,k1,j) = max(0.0, Statein%qgrs(:,k2,j)) + enddo + enddo + + ! Input data from toa to sfc + if (ivflip == 0) then + plvl(1:IM-1,1+kd) = Statein%prsi(1:IM-1,1) + plvl(IM,1+kd) = kdist_lw%get_press_min() + if (lsk /= 0) then + plvl(1:IM-1,1+kd) = 0.5 * (plvl(1:IM-1,2+kd) + plvl(1:IM-1,1+kd)) + plvl(IM,1+kd) = kdist_lw%get_press_min() + endif + ! Input data from sfc to top + else + plvl(1:IM-1,LP1+kd) = Statein%prsi(1:IM-1,LP1+lsk) + plvl(IM,LP1+kd) = kdist_lw%get_press_min() + if (lsk /= 0) then + plvl(1:IM-1,LM+kd) = 0.5 * (plvl(1:IM-1,LP1+kd) + plvl(1:IM-1,LM+kd)) + plvl(IM,LM+kd) = kdist_lw%get_press_min() + endif + endif + + if ( lextop ) then ! values for extra top layer + do i = 1, IM + plvl(i,llb) = kdist_lw%get_press_min() + if ( plvl(i,lla) <= kdist_lw%get_press_min() ) plvl(i,lla) = 2.0* kdist_lw%get_press_min() + plyr(i,lyb) = 0.5 * plvl(i,lla) + tlyr(i,lyb) = tlyr(i,lya) + prslk1(i,lyb) = (plyr(i,lyb)*0.001) ** rocp ! plyr in Pa + rhly(i,lyb) = rhly(i,lya) + qstl(i,lyb) = qstl(i,lya) + enddo - end subroutine GFS_rrtmgp_pre_run - + ! note: may need to take care the top layer amount + tracer1(:,lyb,:) = tracer1(:,lya,:) + endif + + ! Get layer ozone mass mixing ratio + if (Model%ntoz > 0) then + do k=1,lmk + do i=1,im + olyr(i,k) = max( QMIN, tracer1(i,k,Model%ntoz) ) + enddo + enddo + ! Use climatological ozone data + else + call getozn (prslk1, Grid%xlat, IM, LMK, olyr) + endif + + ! Compute cosine of zenith angle (only when SW is called) + if (Model%lsswr) then + call coszmn (Grid%xlon, Grid%sinlat, Grid%coslat, Model%solhr, IM, me, & + Radtend%coszen, Radtend%coszdg) + endif + + ! Call getgases(), to set up non-prognostic gas volume mixing ratios (gasvmr). + call getgases (plvl/100., Grid%xlon, Grid%xlat, IM, LMK, gasvmr) + + ! Get temperature at layer interface, and layer moisture. + tem2da(1:IM,2:LMK) = log( plyr(1:IM,2:LMK) ) + tem2db(1:IM,2:LMK) = log( plvl(1:IM,2:LMK) ) + + if (ivflip == 0) then ! input data from toa to sfc + do i = 1, IM + tem1d (i) = QME6 + tem2da(i,1) = log( plyr(i,1) ) + tem2db(i,1) = log( max( kdist_lw%get_press_min(), plvl(i,1)) ) + tem2db(i,LMP) = log( plvl(i,LMP) ) + tsfa (i) = tlyr(i,LMK) ! sfc layer air temp + tlvl(i,1) = tlyr(i,1) + tlvl(i,LMP) = tskn(i) + enddo + do k = 1, LM + k1 = k + kd + do i = 1, IM + qlyr(i,k1) = max( tem1d(i), Statein%qgrs(i,k,1) ) + tem1d(i) = min( QME5, qlyr(i,k1) ) + tvly(i,k1) = Statein%tgrs(i,k) * (1.0 + fvirt*qlyr(i,k1)) ! virtual T (K) + delp(i,k1) = plvl(i,k1+1) - plvl(i,k1) + enddo + enddo + + if ( lextop ) then + do i = 1, IM + qlyr(i,lyb) = qlyr(i,lya) + tvly(i,lyb) = tvly(i,lya) + delp(i,lyb) = plvl(i,lla) - plvl(i,llb) + enddo + endif + + do k = 2, LMK + do i = 1, IM + tlvl(i,k) = tlyr(i,k) + (tlyr(i,k-1) - tlyr(i,k)) * (tem2db(i,k) - tem2da(i,k)) / & + (tem2da(i,k-1) - tem2da(i,k)) + enddo + enddo + + ! Comput lvel height and layer thickness (km) + tem0d = 0.001 * rog + do i = 1, IM + do k = 1, LMK + dz(i,k) = tem0d * (tem2db(i,k+1) - tem2db(i,k)) * tvly(i,k) + enddo + enddo + + else ! input data from sfc to toa + do i = 1, IM + tem1d (i) = QME6 + tem2da(i,1) = log( plyr(i,1) ) + tem2db(i,1) = log( plvl(i,1) ) + tem2db(i,LMP) = log( max( kdist_lw%get_press_min(), plvl(i,LMP)) ) + tsfa (i) = tlyr(i,1) ! sfc layer air temp + tlvl(i,1) = tskn(i) + tlvl(i,LMP) = tlyr(i,LMK) + enddo + + do k = LM, 1, -1 + do i = 1, IM + qlyr(i,k) = max( tem1d(i), Statein%qgrs(i,k,1) ) + tem1d(i) = min( QME5, qlyr(i,k) ) + tvly(i,k) = Statein%tgrs(i,k) * (1.0 + fvirt*qlyr(i,k)) ! virtual T (K) + delp(i,k) = plvl(i,k) - plvl(i,k+1) + enddo + enddo + + if ( lextop ) then + do i = 1, IM + qlyr(i,lyb) = qlyr(i,lya) + tvly(i,lyb) = tvly(i,lya) + delp(i,lyb) = plvl(i,lla) - plvl(i,llb) + enddo + endif + + do k = 1, LMK-1 + do i = 1, IM + tlvl(i,k+1) = tlyr(i,k) + (tlyr(i,k+1) - tlyr(i,k)) * (tem2db(i,k+1) - tem2da(i,k)) / & + (tem2da(i,k+1) - tem2da(i,k)) + enddo + enddo + + ! Compute level height and layer thickness (km) + tem0d = 0.001 * rog + do i = 1, IM + do k = LMK, 1, -1 + dz(i,k) = tem0d * (tem2db(i,k) - tem2db(i,k+1)) * tvly(i,k) + enddo + enddo + endif ! end_if_ivflip + + ! Obtain cloud information for radiation calculations + ! (clouds,cldsa,mtopa,mbota) + ! for prognostic cloud: + ! - For Zhao/Moorthi's prognostic cloud scheme, + ! call module_radiation_clouds::progcld1() + ! - For Zhao/Moorthi's prognostic cloud+pdfcld, + ! call module_radiation_clouds::progcld3() + ! call module_radiation_clouds::progclduni() for unified cloud and ncld=2 + ccnd = 0.0_kind_phys + if (Model%ncnd == 1) then ! Zhao_Carr_Sundqvist + ccnd(1:IM,1:LMK,1) = tracer1(1:IM,1:LMK,ntcw) ! -liquid water/ice + elseif (Model%ncnd == 2) then ! MG + ccnd(1:IM,1:LMK,1) = tracer1(1:IM,1:LMK,ntcw) ! -liquid water + ccnd(1:IM,1:LMK,2) = tracer1(1:IM,1:LMK,ntiw) ! -ice water + elseif (Model%ncnd == 4) then ! MG2 + ccnd(1:IM,1:LMK,1) = tracer1(1:IM,1:LMK,ntcw) ! -liquid water + ccnd(1:IM,1:LMK,2) = tracer1(1:IM,1:LMK,ntiw) ! -ice water + ccnd(1:IM,1:LMK,3) = tracer1(1:IM,1:LMK,ntrw) ! -rain water + ccnd(1:IM,1:LMK,4) = tracer1(1:IM,1:LMK,ntsw) ! -snow water + elseif (Model%ncnd == 5) then ! GFDL MP, Thompson, MG3 + ccnd(1:IM,1:LMK,1) = tracer1(1:IM,1:LMK,ntcw) ! -liquid water + ccnd(1:IM,1:LMK,2) = tracer1(1:IM,1:LMK,ntiw) ! -ice water + ccnd(1:IM,1:LMK,3) = tracer1(1:IM,1:LMK,ntrw) ! -rain water + ccnd(1:IM,1:LMK,4) = tracer1(1:IM,1:LMK,ntsw) + & ! -snow + grapuel + tracer1(1:IM,1:LMK,ntgl) + endif + where(ccnd < epsq) ccnd = 0.0 + + if (Model%imp_physics == 11 ) then + if (.not. Model%lgfdlmprad) then + ccnd(:,:,1) = tracer1(:,1:LMK,ntcw) + ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:LMK,ntrw) + ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:LMK,ntiw) + ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:LMK,ntsw) + ccnd(:,:,1) = ccnd(:,:,1) + tracer1(:,1:LMK,ntgl) + endif + do k=1,LMK + do i=1,IM + if (ccnd(i,k,1) < EPSQ ) ccnd(i,k,1) = 0.0 + enddo + enddo + endif + + ! Add suspended convective cloud water to grid-scale cloud water + ! only for cloud fraction & radiation computation it is to enhance + ! cloudiness due to suspended convec cloud water for zhao/moorthi's + ! (imp_phys=99) & ferrier's (imp_phys=5) microphysics schemes + if ((Model%num_p3d == 4) .and. (Model%npdf3d == 3)) then ! same as Model%imp_physics = 99 + deltaq(1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,5) + cnvw (1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,6) + cnvc (1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,7) + elseif ((Model%npdf3d == 0) .and. (Model%ncnvcld3d == 1)) then ! same as MOdel%imp_physics=98 + deltaq(1:im,1+kd:lm+kd) = 0.0 + cnvw (1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,Model%num_p3d+1) + cnvc (1:im,1+kd:lm+kd) = 0.0 + else ! all the rest + deltaq(1:im,1:lmk) = 0.0 + cnvw (1:im,1:lmk) = 0.0 + cnvc (1:im,1:lmk) = 0.0 + endif + + if (lextop) then + cldcov(1:im,lyb) = cldcov(1:im,lya) + deltaq(1:im,lyb) = deltaq(1:im,lya) + cnvw (1:im,lyb) = cnvw (1:im,lya) + cnvc (1:im,lyb) = cnvc (1:im,lya) + if (Model%effr_in) then + effrl(1:im,lyb) = effrl(1:im,lya) + effri(1:im,lyb) = effri(1:im,lya) + effrr(1:im,lyb) = effrr(1:im,lya) + effrs(1:im,lyb) = effrs(1:im,lya) + endif + endif + + if (Model%imp_physics == 99) then + ccnd(1:IM,1:LMK,1) = ccnd(1:IM,1:LMK,1) + cnvw(1:IM,1:LMK) + endif + + if (Model%imp_physics == 10) then + ccnd(1:IM,1:LMK,1) = ccnd(1:IM,1:LMK,1) + cnvw(1:IM,1:LMK) + ccnd(1:IM,1:LMK,2) + endif + + ! DJS2019: START + ! Compute layer cloud fraction. + clwmin = 0.0 + cldcov(:,:) = 0.0 + if (.not. Model%lmfshal) then + do k = 1, LMK + do i = 1, IM + clwt = 1.0e-6 * (plyr(i,k)*0.1) + if (ccnd(i,k,1) > 0.) then + onemrh= max( 1.e-10, 1.0-rhly(i,k) ) + clwm = clwmin / max( 0.01, plyr(i,k)*0.1 ) + tem1 = min(max(sqrt(sqrt(onemrh*qstl(i,k))),0.0001),1.0) + tem1 = 2000.0 / tem1 + value = max( min( tem1*(ccnd(i,k,1)-clwm), 50.0 ), 0.0 ) + tem2 = sqrt( sqrt(rhly(i,k)) ) + cldcov(i,k) = max( tem2*(1.0-exp(-value)), 0.0 ) + endif + enddo + enddo + else + do k = 1, LMK + do i = 1, IM + clwt = 1.0e-6 * (plyr(i,k)*0.1) + if (ccnd(i,k,1) .gt. 0) then + onemrh= max( 1.e-10, 1.0-rhly(i,k) ) + clwm = clwmin / max( 0.01, plyr(i,k)*0.1 ) + tem1 = min(max((onemrh*qstl(i,k))**0.49,0.0001),1.0) !jhan + if (Model%lmfdeep2) then + tem1 = xrc3 / tem1 + else + tem1 = 100.0 / tem1 + endif + value = max( min( tem1*(ccnd(i,k,1)-clwm), 50.0 ), 0.0 ) + tem2 = sqrt( sqrt(rhly(i,k)) ) + cldcov(i,k) = max( tem2*(1.0-exp(-value)), 0.0 ) + endif + enddo + enddo + endif + ! DJS2019: END + + if (Model%uni_cld) then + if (Model%effr_in) then + cldcov(1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,Model%indcld) + effrl(1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,2) + effri(1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,3) + effrr(1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,4) + effrs(1:im,1+kd:lm+kd) = Tbd%phy_f3d(1:im,1:lm,5) + else + do k=1,lm + k1 = k + kd + do i=1,im + !cldcov(i,k1) = Tbd%phy_f3d(i,k,Model%indcld) + !if (tracer1(i,k,ntcw) .gt. 0 .or. tracer1(i,k,ntiw) .gt. 0) then + ! cldcov(i,k1) = 0.1 + !else + ! cldcov(i,k1) = 0.0 + !endif + enddo + enddo + endif + elseif (Model%imp_physics == Model%imp_physics_gfdl) then ! GFDL MP + cldcov(1:IM,1+kd:LM+kd) = tracer1(1:IM,1:LM,Model%ntclamt) + else ! neither of the other two cases + ! cldcov = 0.0 + endif + + ! MICROPHYSICS + ! *) zhao/moorthi's prognostic cloud scheme or unified cloud and/or with MG microphysics + if (Model%imp_physics == 99 .or. Model%imp_physics == 10) then + if (Model%uni_cld .and. Model%ncld >= 2) then + call progclduni (plyr/100., plvl/100., tlyr, tvly, ccnd, Model%ncnd, & ! IN + Grid%xlat, Grid%xlon, Sfcprop%slmsk, dz, delp/100.,IM, & ! IN + LMK, LMP, cldcov, effrl, effri, effrr, effrs, Model%effr_in, & ! IN + clouds, cldsa, mtopa, mbota, de_lgth) ! OUT + else + call progcld1 (plyr/100. ,plvl/100., tlyr, tvly, qlyr, qstl, rhly, & ! IN + ccnd(1:IM,1:LMK,1), Grid%xlat,Grid%xlon,Sfcprop%slmsk, dz, & ! IN + delp/100., IM, LMK, LMP, Model%uni_cld, Model%lmfshal, & ! IN + Model%lmfdeep2, cldcov, effrl, effri, effrr, effrs, & ! IN + Model%effr_in, & ! IN + clouds, cldsa, mtopa, mbota, de_lgth) ! OUT + endif + ! *) zhao/moorthi's prognostic cloud+pdfcld + elseif(Model%imp_physics == 98) then + call progcld3 (plyr/100., plvl/100., tlyr, tvly, qlyr, qstl, rhly, & ! IN + ccnd(1:IM,1:LMK,1), cnvw, cnvc, Grid%xlat, Grid%xlon, & ! IN + Sfcprop%slmsk, dz, delp/100., im, lmk, lmp, deltaq, Model%sup, & ! IN + Model%kdt, me, & ! IN + clouds, cldsa, mtopa, mbota, de_lgth) ! OUT + ! *) GFDL cloud scheme + elseif (Model%imp_physics == 11) then + if (.not.Model%lgfdlmprad) then + call progcld4 (plyr/100., plvl/100., tlyr, tvly, qlyr, qstl, rhly, & ! IN + ccnd(1:IM,1:LMK,1), cnvw, cnvc,Grid%xlat, Grid%xlon, & ! IN + Sfcprop%slmsk, cldcov, dz, delp/100., im, lmk, lmp, & ! IN + clouds, cldsa, mtopa, mbota, de_lgth) ! OUT + else + call progclduni (plyr/100., plvl/100., tlyr, tvly, ccnd, Model%ncnd, & ! IN + Grid%xlat, Grid%xlon, Sfcprop%slmsk, dz,delp/100., IM, LMK, & ! IN + LMP, cldcov, effrl, effri, effrr, effrs, Model%effr_in, & ! IN + clouds, cldsa, mtopa, mbota, de_lgth) ! OUT + endif + ! *) Thompson / WSM6 cloud micrphysics scheme + elseif(Model%imp_physics == 8 .or. Model%imp_physics == 6) then + if (Model%kdt == 1) then + Tbd%phy_f3d(:,:,1) = 10. + Tbd%phy_f3d(:,:,2) = 50. + Tbd%phy_f3d(:,:,3) = 250. + endif + + call progcld5 (plyr/100., plvl/100., tlyr, qlyr, qstl, rhly, tracer1, Grid%xlat, & ! IN + Grid%xlon,Sfcprop%slmsk,dz,delp/100., ntrac-1, ntcw-1, ntiw-1, & ! IN + ntrw-1, ntsw-1, ntgl-1, im, lmk, lmp, Model%uni_cld, Model%lmfshal,& ! IN + Model%lmfdeep2, cldcov(:,1:LMK),Tbd%phy_f3d(:,:,1), & ! IN + Tbd%phy_f3d(:,:,2), Tbd%phy_f3d(:,:,3), & ! IN + clouds,cldsa,mtopa,mbota, de_lgth) ! OUT + endif ! end if_imp_physics + + + ! mg, sfc-perts + ! --- scale random patterns for surface perturbations with + ! perturbation size + ! --- turn vegetation fraction pattern into percentile pattern + alb1d(:) = 0. + if (Model%do_sfcperts) then + if (Model%pertalb(1) > 0.) then + do i=1,im + call cdfnor(Coupling%sfc_wts(i,5),alb1d(i)) + enddo + endif + endif + ! mg, sfc-perts + + + ! ####################################################################################### + ! Call module_radiation_aerosols::setaer(),to setup aerosols property profile for both + ! LW and SW radiation. + ! ####################################################################################### + call setaer (plvl, plyr, prslk1, tvly, rhly, Sfcprop%slmsk, tracer1, Grid%xlon, & + Grid%xlat, IM, LMK, LMP, Model%lsswr, Model%lslwr, faersw, faerlw, aerodp) + + ! Store aerosol optical properties + ! SW. + ! For RRTMGP SW the bands are now ordered from [IR(band) -> nIR -> UV], in RRTMG the + ! band ordering was [nIR -> UV -> IR(band)] + faersw2(1:IM,1:LMK,1,1) = faersw(1:IM,1:LMK,kdist_sw%get_nband(),1) + faersw2(1:IM,1:LMK,1,2) = faersw(1:IM,1:LMK,kdist_sw%get_nband(),2) + faersw2(1:IM,1:LMK,1,3) = faersw(1:IM,1:LMK,kdist_sw%get_nband(),3) + faersw2(1:IM,1:LMK,2:kdist_sw%get_nband(),1) = faersw(1:IM,1:LMK,1:kdist_sw%get_nband()-1,1) + faersw2(1:IM,1:LMK,2:kdist_sw%get_nband(),2) = faersw(1:IM,1:LMK,1:kdist_sw%get_nband()-1,2) + faersw2(1:IM,1:LMK,2:kdist_sw%get_nband(),3) = faersw(1:IM,1:LMK,1:kdist_sw%get_nband()-1,3) + + ! ####################################################################################### + ! Call module_radiation_surface::setemis(),to setup surface emissivity for LW radiation. + ! ####################################################################################### + if (Model%lslwr) then + call setemis (Grid%xlon, Grid%xlat, Sfcprop%slmsk, Sfcprop%snowd, Sfcprop%sncovr, & + Sfcprop%zorl, tsfg, tsfa, Sfcprop%hprim, IM, Radtend%semis) + do iBand=1,kdist_lw%get_nband() + sfc_emiss_byband(iBand,1:IM) = Radtend%semis(1:IM) + enddo + endif + + ! ####################################################################################### + ! Check for daytime points for SW radiation. + ! ####################################################################################### + if (Model%lsswr) then + nday = 0 + idxday = 0 + do iCol = 1, IM + if (Radtend%coszen(iCol) >= 0.0001) then + nday = nday + 1 + idxday(nday) = iCol + endif + enddo + else + nday = 0 + idxday = 0 + endif + + ! ####################################################################################### + ! Call module_radiation_surface::setalb() to setup surface albedo for SW radiation. + ! ####################################################################################### + if (Model%lsswr) then + call setalb (Sfcprop%slmsk, Sfcprop%snowd, Sfcprop%sncovr, Sfcprop%snoalb, & + Sfcprop%zorl, Radtend%coszen, tsfg, tsfa, Sfcprop%hprim, Sfcprop%alvsf, & + Sfcprop%alnsf, Sfcprop%alvwf, Sfcprop%alnwf, Sfcprop%facsf, & + Sfcprop%facwf, Sfcprop%fice, Sfcprop%tisfc, IM, alb1d, Model%pertalb, & + sfcalb) + + ! Approximate mean surface albedo from vis- and nir- diffuse values. + Radtend%sfalb(:) = max(0.01, 0.5 * (sfcalb(:,2) + sfcalb(:,4))) + + ! Spread across all SW bands + do iBand=1,kdist_sw%get_nband() + sfc_alb_nir_dir(iBand,1:IM) = sfcalb(:,1) + sfc_alb_nir_dif(iBand,1:IM) = sfcalb(:,2) + sfc_alb_uvvis_dir(iBand,1:IM) = sfcalb(:,3) + sfc_alb_uvvis_dif(iBand,1:IM) = sfcalb(:,4) + enddo + else + sfc_alb_nir_dir(:,:) = 0._kind_phys + sfc_alb_nir_dif(:,:) = 0._kind_phys + sfc_alb_uvvis_dir(:,:) = 0._kind_phys + sfc_alb_uvvis_dif(:,:) = 0._kind_phys + endif + + ! ####################################################################################### + ! Compute radiative properties needed for RRTMGP + ! ####################################################################################### + + ! Change random number seed value for each radiation invocation (isubclw =1 or 2). + if(isubclw == 1) then ! advance prescribed permutation seed + do iCol = 1, IM + ipseed_lw(iCol) = ipsdlw0 + iCol + enddo + elseif (isubclw == 2) then ! use input array of permutaion seeds + do iCol = 1, IM + ipseed_lw(iCol) = icseed(iCol) + enddo + endif + ! Change random number seed value for each radiation invocation (isubcsw =1 or 2). + if(isubcsw == 1) then ! advance prescribed permutation seed + do iCol = 1, ncol + ipseed_sw(iCol) = ipsdsw0 + iCol + enddo + elseif (isubcsw == 2) then ! use input array of permutaion seeds + do iCol = 1, ncol + ipseed_sw(iCol) = icseed(iCol) + enddo + endif + + ! Compute volume mixing-ratios for ozone (mmr) and specific-humidity. + vmr_h2o = merge((qlyr/(1-qlyr))*amdw, 0., qlyr .ne. 1.) + vmr_o3 = merge(olyr*amdo3, 0., olyr .gt. 0.) + + ! Compute ice/liquid cloud masks, needed by rrtmgp_cloud_optics + liqmask = (clouds(:,:,1) .gt. 0 .and. clouds(:,:,2) .gt. 0) + icemask = (clouds(:,:,1) .gt. 0 .and. clouds(:,:,4) .gt. 0) + + ! ####################################################################################### + ! Allocate space for gas optical properties [ncol,nlay,ngpts] + ! ####################################################################################### + ! Longwave + if (Model%lslwr) then + ! Cloud optics [nCol,nLay,nBands] + call check_error_msg('GFS_rrtmgp_pre_run',optical_propsLW_cloudsByBand%init(kdist_lw%get_band_lims_wavenumber())) + call check_error_msg('GFS_rrtmgp_pre_run',optical_propsLW_cloudsByBand%alloc_1scl(IM, LMK)) + ! Aerosol optics [Ccol,nLay,nBands] + call check_error_msg('GFS_rrtmgp_pre_run',optical_propsLW_aerosol%init(kdist_lw%get_band_lims_wavenumber())) + call check_error_msg('GFS_rrtmgp_pre_run',optical_propsLW_aerosol%alloc_1scl(IM, LMK)) + ! Cloud optics [nCol,nLay,nGpts] + call check_error_msg('GFS_rrtmgp_pre_run',optical_propsLW_clouds%alloc_1scl(IM, LMK, kdist_lw)) + endif + ! Shortwave + if (Model%lsswr .and. nday .gt. 0) then + ! Cloud optics [nCol,nLay,nBands] + call check_error_msg('GFS_rrtmgp_pre_run',optical_propsSW_cloudsByBand%init(kdist_sw%get_band_lims_wavenumber())) + call check_error_msg('GFS_rrtmgp_pre_run',optical_propsSW_cloudsByBand%alloc_2str(nDay, LMK)) + ! Aerosol optics [Ccol,nLay,nBands] + call check_error_msg('GFS_rrtmgp_pre_run',optical_propsSW_aerosol%init(kdist_sw%get_band_lims_wavenumber())) + call check_error_msg('GFS_rrtmgp_pre_run',optical_propsSW_aerosol%alloc_2str(nDay, LMK)) + ! Cloud optics [nCol,nLay,nGpts] + call check_error_msg('GFS_rrtmgp_pre_run',optical_propsSW_clouds%alloc_2str(nDay, LMK, kdist_sw)) + endif + + ! ####################################################################################### + ! Set gas concentrations + ! ####################################################################################### + !if (Model%lslwr) then + call gas_concentrations_lw%reset() + call check_error_msg('GFS_rrtmgp_pre_run',gas_concentrations_lw%set_vmr('o2', gasvmr(:,:,4))) + call check_error_msg('GFS_rrtmgp_pre_run',gas_concentrations_lw%set_vmr('co2', gasvmr(:,:,1))) + call check_error_msg('GFS_rrtmgp_pre_run',gas_concentrations_lw%set_vmr('ch4', gasvmr(:,:,3))) + call check_error_msg('GFS_rrtmgp_pre_run',gas_concentrations_lw%set_vmr('n2o', gasvmr(:,:,2))) + call check_error_msg('GFS_rrtmgp_pre_run',gas_concentrations_lw%set_vmr('h2o', vmr_h2o)) + call check_error_msg('GFS_rrtmgp_pre_run',gas_concentrations_lw%set_vmr('o3', vmr_o3)) + !endif + !if (Model%lsswr .and. nday .gt. 0) then + call gas_concentrations_sw%reset() + call check_error_msg('GFS_rrtmgp_pre_run',gas_concentrations_sw%set_vmr('o2', gasvmr(idxday,:,4))) + call check_error_msg('GFS_rrtmgp_pre_run',gas_concentrations_sw%set_vmr('co2', gasvmr(idxday,:,1))) + call check_error_msg('GFS_rrtmgp_pre_run',gas_concentrations_sw%set_vmr('ch4', gasvmr(idxday,:,3))) + call check_error_msg('GFS_rrtmgp_pre_run',gas_concentrations_sw%set_vmr('n2o', gasvmr(idxday,:,2))) + call check_error_msg('GFS_rrtmgp_pre_run',gas_concentrations_sw%set_vmr('h2o', vmr_h2o(idxday,:))) + call check_error_msg('GFS_rrtmgp_pre_run',gas_concentrations_sw%set_vmr('o3', vmr_o3(idxday,:))) + !endif + + ! ####################################################################################### + ! Copy aerosol to RRTMGP DDT + ! ####################################################################################### + ! LW + if (Model%lslwr) then + optical_propsLW_aerosol%tau = faerlw(:,:,:,1) * (1. - faerlw(:,:,:,2)) + endif + ! SW + if (Model%lsswr .and. nday .gt. 0) then + optical_propsSW_aerosol%tau = faersw2(idxday,:,:,1) + optical_propsSW_aerosol%ssa = faersw2(idxday,:,:,2) + optical_propsSW_aerosol%g = faersw2(idxday,:,:,3) + endif + + ! ####################################################################################### + ! Compute cloud-optics for RTE. + ! ####################################################################################### + ! Longwave + if (Model%lslwr) then + call check_error_msg('GFS_rrtmgp_pre_run',kdist_cldy_lw%cloud_optics(IM, LMK, kdist_lw%get_nband(), & + nrghice_lw, liqmask, icemask, clouds(:,:,2), clouds(:,:,4), clouds(:,:,3), & + clouds(:,:,5), optical_propsLW_cloudsByBand)) + endif + ! Shortwave + if (Model%lsswr .and. nday .gt. 0) then + call check_error_msg('GFS_rrtmgp_pre_run',kdist_cldy_sw%cloud_optics(nDay, LMK, kdist_sw%get_nband(), & + nrghice_sw, liqmask(idxday,:), icemask(idxday,:), clouds(idxday,:,2), & + clouds(idxday,:,4), clouds(idxday,:,3), clouds(idxday,:,5), & + optical_propsSW_cloudsByBand)) + endif + + ! ####################################################################################### + ! Call McICA to generate subcolumns. + ! ####################################################################################### + ! Longwave + if (Model%lslwr .and. isubclw .gt. 0) then + + ! Call RNG. Mersennse Twister accepts 1D array, so loop over columns and collapse along G-points + ! and layers. ([nGpts,nLayer,nColumn]-> [nGpts*nLayer]*nColumn) + do iCol=1,IM + call random_setseed(ipseed_lw(icol),rng_stat) + call random_number(rng1D_lw,rng_stat) + rng3D_lw(:,:,iCol) = reshape(source = rng1D_lw,shape=[kdist_lw%get_ngpt(),LMK]) + enddo + + ! Call McICA + select case ( iovrlw ) + ! Maximumn-random + case(1) + call check_error_msg('GFS_rrtmgp_pre_run',sampled_mask_max_ran(rng3D_lw,clouds(:,:,1),cldfracMCICA_lw)) + end select + + ! Map band optical depth to each g-point using McICA + call check_error_msg('GFS_rrtmgp_pre_run',draw_samples(cldfracMCICA_lw,optical_propsLW_cloudsByBand,optical_propsLW_clouds)) + endif + + ! Shortwave + if (Model%lsswr .and. nday .gt. 0 .and. isubcsw .gt. 0) then + + ! Call RNG. Mersennse Twister accepts 1D array, so loop over columns and collapse along G-points + ! and layers. ([nGpts,nLayer,nColumn]-> [nGpts*nLayer]*nColumn) + do iCol=1,IM + call random_setseed(ipseed_sw(icol),rng_stat) + call random_number(rng1D_sw,rng_stat) + rng3D_sw(:,:,iCol) = reshape(source = rng1D_sw,shape=[kdist_sw%get_ngpt(),LMK]) + enddo + + ! Call McICA + select case ( iovrsw ) + ! Maximumn-random + case(1) + call check_error_msg('GFS_rrtmgp_pre_run',sampled_mask_max_ran(rng3D_sw,clouds(:,:,1),cldfracMCICA_sw)) + end select + + ! Map band optical depth to each g-point using McICA + call check_error_msg('GFS_rrtmgp_pre_run',draw_samples(cldfracMCICA_sw,optical_propsSW_cloudsByBand,optical_propsSW_clouds)) + endif + + end subroutine GFS_rrtmgp_pre_run + !> \section arg_table_GFS_rrtmgp_pre_finalize Argument Table !! - subroutine GFS_rrtmgp_pre_finalize () - close(58) - end subroutine GFS_rrtmgp_pre_finalize - + subroutine GFS_rrtmgp_pre_finalize () + end subroutine GFS_rrtmgp_pre_finalize + !! @} - end module GFS_rrtmgp_pre +end module GFS_rrtmgp_pre diff --git a/physics/rrtmgp_lw.F90 b/physics/rrtmgp_lw.F90 index 183bc8db6..3007146d7 100644 --- a/physics/rrtmgp_lw.F90 +++ b/physics/rrtmgp_lw.F90 @@ -10,6 +10,8 @@ module rrtmgp_lw use mo_rrtmgp_clr_all_sky, only: rte_lw use mo_gas_concentrations, only: ty_gas_concs use mo_fluxes_byband, only: ty_fluxes_byband + use rrtmgp_sw, only: check_error_msg + ! Parameters integer,parameter :: nGases = 6 @@ -441,9 +443,9 @@ subroutine rrtmgp_lw_init(Model, mpicomm, mpirank, mpiroot, kdist_lw, kdist_cldy ! Initialize gas concentrations and gas optics class with data do iGas=1,nGases - call check_error_msg(gas_concentrations%set_vmr(active_gases(iGas), 0._kind_phys)) + call check_error_msg('rrtmgp_lw_init',gas_concentrations%set_vmr(active_gases(iGas), 0._kind_phys)) enddo - call check_error_msg(kdist_lw%load(gas_concentrations, gas_names, key_species, band2gpt, & + call check_error_msg('rrtmgp_lw_init',kdist_lw%load(gas_concentrations, gas_names, key_species, band2gpt, & band_lims, press_ref, press_ref_trop, temp_ref, temp_ref_p, temp_ref_t, & vmr_ref, kmajor, kminor_lower, kminor_upper, gas_minor,identifier_minor, & minor_gases_lower, minor_gases_upper, minor_limits_gpt_lower, & @@ -642,14 +644,14 @@ subroutine rrtmgp_lw_init(Model, mpicomm, mpirank, mpiroot, kdist_lw, kdist_cldy ! Load tables data for RRTGMP cloud-optics if (rrtmgp_lw_cld_phys .eq. 1) then - call check_error_msg(kdist_cldy_lw%set_ice_roughness(nrghice)) - call check_error_msg(kdist_cldy_lw%load(band_lims_cldy, radliq_lwr, radliq_upr, & + call check_error_msg('rrtmgp_lw_init',kdist_cldy_lw%set_ice_roughness(nrghice)) + call check_error_msg('rrtmgp_lw_init',kdist_cldy_lw%load(band_lims_cldy, radliq_lwr, radliq_upr, & radliq_fac, radice_lwr, radice_upr, radice_fac, lut_extliq, lut_ssaliq, & lut_asyliq, lut_extice, lut_ssaice, lut_asyice)) endif if (rrtmgp_lw_cld_phys .eq. 2) then - call check_error_msg(kdist_cldy_lw%set_ice_roughness(nrghice)) - call check_error_msg(kdist_cldy_lw%load(band_lims_cldy, pade_extliq, & + call check_error_msg('rrtmgp_lw_init',kdist_cldy_lw%set_ice_roughness(nrghice)) + call check_error_msg('rrtmgp_lw_init',kdist_cldy_lw%load(band_lims_cldy, pade_extliq, & pade_ssaliq, pade_asyliq, pade_extice, pade_ssaice, pade_asyice, & pade_sizereg_extliq, pade_sizereg_ssaliq, pade_sizereg_asyliq, & pade_sizereg_extice, pade_sizereg_ssaice, pade_sizereg_asyice)) @@ -660,58 +662,57 @@ end subroutine rrtmgp_lw_init ! ######################################################################################### ! ######################################################################################### !! \section arg_table_rrtmgp_lw_run Argument Table -!! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | -!! |-----------------------|-----------------------------------------------------------------------------------------------|--------------------------------------------------------------------|-------|------|-----------------------|-----------|--------|----------| -!! | ncol | horizontal_loop_extent | horizontal dimension | count | 0 | integer | | in | F | -!! | nlay | adjusted_vertical_layer_dimension_for_radiation | number of vertical layers for radiation | count | 0 | integer | | in | F | -!! | p_lay | air_pressure_at_layer_for_radiation_in_hPa | air pressure layer | hPa | 2 | real | kind_phys | in | F | -!! | p_lev | air_pressure_at_interface_for_radiation_in_hPa | air pressure level | hPa | 2 | real | kind_phys | in | F | -!! | t_lay | air_temperature_at_layer_for_radiation | air temperature layer | K | 2 | real | kind_phys | in | F | -!! | skt | surface_ground_temperature_for_radiation | surface ground temperature for radiation | K | 1 | real | kind_phys | in | F | -!! | sfc_emiss | surface_longwave_emissivity_in_each_band | surface lw emissivity in fraction in each LW band | frac | 2 | real | kind_phys | in | F | -!! | kdist_lw | K_distribution_file_for_RRTMGP_LW_scheme | DDT containing spectral information for RRTMGP LW radiation scheme | DDT | 0 | ty_gas_optics_rrtmgp | | in | F | -!! | optical_props_clds | optical_properties_for_cloudy_atmosphere | Fortran DDT containing RRTMGP optical properties | DDT | 0 | ty_optical_props_1scl | | in | F | -!! | optical_props_aerosol | optical_properties_for_aerosols | Fortran DDT containing RRTMGP optical properties | DDT | 0 | ty_optical_props_1scl | | in | F | -!! | gas_concentrations | Gas_concentrations_for_RRTMGP_suite | DDT containing gas concentrations for RRTMGP radiation scheme | DDT | 0 | ty_gas_concs | | in | F | -!! | lslwr | flag_to_calc_lw | flag to calculate LW irradiances | flag | 0 | logical | | in | F | -!! | fluxLW_allsky | lw_flux_profiles_byband_allsky | Fortran DDT containing RRTMGP 3D fluxes | DDT | 0 | ty_fluxes_byband | | out | F | -!! | fluxLW_clrsky | lw_flux_profiles_byband_clrsky | Fortran DDT containing RRTMGP 3D fluxes | DDT | 0 | ty_fluxes_byband | | out | F | -!! | hlw0 | tendency_of_air_temperature_due_to_longwave_heating_assuming_clear_sky_on_radiation_time_step | longwave clear sky heating rate | K s-1 | 2 | real | kind_phys | inout | T | -!! | hlwb | lw_heating_rate_spectral | longwave total sky heating rate (spectral) | K s-1 | 3 | real | kind_phys | inout | T | -!! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | -!! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | +!! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | +!! |-------------------------|-----------------------------------------------------------------------------------------------|--------------------------------------------------------------------|-------|------|-----------------------|-----------|--------|----------| +!! | ncol | horizontal_loop_extent | horizontal dimension | count | 0 | integer | | in | F | +!! | nlay | adjusted_vertical_layer_dimension_for_radiation | number of vertical layers for radiation | count | 0 | integer | | in | F | +!! | p_lay | air_pressure_at_layer_for_radiation_in_hPa | air pressure layer | hPa | 2 | real | kind_phys | in | F | +!! | p_lev | air_pressure_at_interface_for_radiation_in_hPa | air pressure level | hPa | 2 | real | kind_phys | in | F | +!! | t_lay | air_temperature_at_layer_for_radiation | air temperature layer | K | 2 | real | kind_phys | in | F | +!! | skt | surface_ground_temperature_for_radiation | surface ground temperature for radiation | K | 1 | real | kind_phys | in | F | +!! | sfc_emiss | surface_longwave_emissivity_in_each_band | surface lw emissivity in fraction in each LW band | frac | 2 | real | kind_phys | in | F | +!! | kdist_lw | K_distribution_file_for_RRTMGP_LW_scheme | DDT containing spectral information for RRTMGP LW radiation scheme | DDT | 0 | ty_gas_optics_rrtmgp | | in | F | +!! | optical_propsLW_clds | longwave_optical_properties_for_cloudy_atmosphere | Fortran DDT containing RRTMGP optical properties | DDT | 0 | ty_optical_props_1scl | | in | F | +!! | optical_propsLW_aerosol | longwave_optical_properties_for_aerosols | Fortran DDT containing RRTMGP optical properties | DDT | 0 | ty_optical_props_1scl | | in | F | +!! | gas_concentrations | Gas_concentrations_for_RRTMGP_suite_lw | DDT containing gas concentrations for RRTMGP radiation scheme | DDT | 0 | ty_gas_concs | | in | F | +!! | lslwr | flag_to_calc_lw | flag to calculate LW irradiances | flag | 0 | logical | | in | F | +!! | fluxLW_allsky | lw_flux_profiles_byband_allsky | Fortran DDT containing RRTMGP 3D fluxes | DDT | 0 | ty_fluxes_byband | | out | F | +!! | fluxLW_clrsky | lw_flux_profiles_byband_clrsky | Fortran DDT containing RRTMGP 3D fluxes | DDT | 0 | ty_fluxes_byband | | out | F | +!! | hlw0 | tendency_of_air_temperature_due_to_longwave_heating_assuming_clear_sky_on_radiation_time_step | longwave clear sky heating rate | K s-1 | 2 | real | kind_phys | in | T | +!! | hlwb | lw_heating_rate_spectral | longwave total sky heating rate (spectral) | K s-1 | 3 | real | kind_phys | in | T | +!! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | +!! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | !! subroutine rrtmgp_lw_run(ncol, nlay, kdist_lw, p_lay, t_lay, p_lev, skt, & - sfc_emiss, gas_concentrations, optical_props_clds, optical_props_aerosol,& + sfc_emiss, gas_concentrations, optical_propsLW_clds, optical_propsLW_aerosol,& lslwr, fluxLW_allsky, fluxLW_clrsky, hlw0, hlwb, errmsg, errflg) ! Inputs integer, intent(in) :: & - ncol, & ! Number of horizontal gridpoints - nlay ! Number of vertical layers + ncol, & ! Number of horizontal gridpoints + nlay ! Number of vertical layers real(kind_phys), dimension(ncol,nlay), intent(in) :: & - p_lay, & ! Pressure @ model layer-centers (hPa) - t_lay ! Temperature (K) + p_lay, & ! Pressure @ model layer-centers (hPa) + t_lay ! Temperature (K) real(kind_phys), dimension(ncol,nlay+1), intent(in) :: & - p_lev ! Pressure @ model layer-interfaces (hPa) + p_lev ! Pressure @ model layer-interfaces (hPa) real(kind_phys), dimension(ncol), intent(in) :: & - skt ! Surface(skin) temperature (K) + skt ! Surface(skin) temperature (K) type(ty_gas_optics_rrtmgp),intent(in) :: & - kdist_lw ! DDT containing LW spectral information + kdist_lw ! DDT containing LW spectral information real(kind_phys), dimension(kdist_lw%get_nband(),ncol) :: & - sfc_emiss ! Surface emissivity (1) + sfc_emiss ! Surface emissivity (1) type(ty_optical_props_1scl),intent(in) :: & - optical_props_clds, & ! RRTMGP DDT: cloud radiative properties - optical_props_aerosol ! RRTMGP DDT: aerosol radiative properties + optical_propsLW_clds, & ! RRTMGP DDT: longwave cloud radiative properties + optical_propsLW_aerosol ! RRTMGP DDT: longwave aerosol radiative properties type(ty_gas_concs),intent(in) :: & - gas_concentrations ! RRTMGP DDT: trace gas concentrations (vmr) + gas_concentrations ! RRTMGP DDT: trace gas concentrations (vmr) logical, intent(in) :: & - lslwr ! Flag to calculate LW irradiances + lslwr ! Flag to calculate LW irradiances ! Outputs character(len=*), intent(out) :: errmsg integer, intent(out) :: errflg - type(ty_fluxes_byband),intent(out) :: & fluxLW_allsky, & ! All-sky flux (W/m2) fluxLW_clrsky ! Clear-sky flux (W/m2) @@ -724,9 +725,9 @@ subroutine rrtmgp_lw_run(ncol, nlay, kdist_lw, p_lay, t_lay, p_lev, skt, & ! Local variables real(kind_phys), dimension(ncol,nlay+1),target :: & - flux_up_allsky, flux_up_clrsky, flux_dn_allsky, flux_dn_clrsky + fluxLW_up_allsky, fluxLW_up_clrsky, fluxLW_dn_allsky, fluxLW_dn_clrsky real(kind_phys), dimension(ncol,nlay+1,kdist_lw%get_nband()),target :: & - fluxBB_up_allsky, fluxBB_dn_allsky + fluxLWBB_up_allsky, fluxLWBB_dn_allsky logical :: l_ClrSky_HR, l_AllSky_HR_byband ! Initialize CCPP error handling variables @@ -739,43 +740,33 @@ subroutine rrtmgp_lw_run(ncol, nlay, kdist_lw, p_lay, t_lay, p_lev, skt, & l_AllSky_HR_byband = present(hlwb) ! Initialize RRTMGP DDT containing 2D(3D) fluxes - fluxLW_allsky%flux_up => flux_up_allsky - fluxLW_allsky%flux_dn => flux_dn_allsky - fluxLW_clrsky%flux_up => flux_up_clrsky - fluxLW_clrsky%flux_dn => flux_dn_clrsky + fluxLW_allsky%flux_up => fluxLW_up_allsky + fluxLW_allsky%flux_dn => fluxLW_dn_allsky + fluxLW_clrsky%flux_up => fluxLW_up_clrsky + fluxLW_clrsky%flux_dn => fluxLW_dn_clrsky ! Only calculate fluxes by-band, only when heating-rate profiles by band are requested. if (l_AllSky_HR_byband) then - fluxLW_allsky%bnd_flux_up => fluxBB_up_allsky - fluxLW_allsky%bnd_flux_dn => fluxBB_dn_allsky + fluxLW_allsky%bnd_flux_up => fluxLWBB_up_allsky + fluxLW_allsky%bnd_flux_dn => fluxLWBB_dn_allsky endif ! Call RRTMGP LW scheme - call check_error_msg(rte_lw( & - kdist_lw, & ! IN - spectral information - gas_concentrations, & ! IN - gas concentrations (vmr) - p_lay, & ! IN - pressure at layer interfaces (Pa) - t_lay, & ! IN - temperature at layer interfaes (K) - p_lev, & ! IN - pressure at layer centers (Pa) - skt, & ! IN - skin temperature (K) - sfc_emiss, & ! IN - surface emissivity in each LW band - optical_props_clds, & ! IN - DDT containing cloud optical information - fluxLW_allsky, & ! OUT - Fluxes, all-sky, 3D (nCol,nLay,nBand) - fluxLW_clrsky, & ! OUT - Fluxes, clear-sky, 3D (nCol,nLay,nBand) - aer_props = optical_props_aerosol)) ! IN(optional) - DDT containing aerosol optical information + call check_error_msg('rrtmgp_lw_run',rte_lw( & + kdist_lw, & ! IN - spectral information + gas_concentrations, & ! IN - gas concentrations (vmr) + p_lay, & ! IN - pressure at layer interfaces (Pa) + t_lay, & ! IN - temperature at layer interfaes (K) + p_lev, & ! IN - pressure at layer centers (Pa) + skt, & ! IN - skin temperature (K) + sfc_emiss, & ! IN - surface emissivity in each LW band + optical_propsLW_clds, & ! IN - DDT containing cloud optical information + fluxLW_allsky, & ! OUT - Fluxes, all-sky, 3D (nCol,nLay,nBand) + fluxLW_clrsky, & ! OUT - Fluxes, clear-sky, 3D (nCol,nLay,nBand) + aer_props = optical_propsLW_aerosol)) ! IN(optional) - DDT containing aerosol optical information end subroutine rrtmgp_lw_run subroutine rrtmgp_lw_finalize() end subroutine rrtmgp_lw_finalize - subroutine check_error_msg(error_msg) - character(len=*), intent(in) :: error_msg - - if(error_msg /= "") then - print*,"ERROR(rrtmgp_lw.F90): " - print*,trim(error_msg) - return - end if - end subroutine check_error_msg - end module rrtmgp_lw diff --git a/physics/rrtmgp_sw.F90 b/physics/rrtmgp_sw.F90 index b1eec1fe1..0961a0086 100644 --- a/physics/rrtmgp_sw.F90 +++ b/physics/rrtmgp_sw.F90 @@ -1,62 +1,28 @@ ! ########################################################################################### ! ########################################################################################### module rrtmgp_sw - use GFS_typedefs, only: GFS_control_type - use physparam, only: iovrsw, icldflg, iswcliq, isubcsw - use machine, only: kind_phys - use mo_rte_kind, only: wl - use mo_gas_optics_rrtmgp, only: ty_gas_optics_rrtmgp - use mo_cloud_optics, only: ty_cloud_optics - use mo_gas_concentrations, only: ty_gas_concs - use mo_fluxes, only: ty_fluxes_broadband - use mo_fluxes_byband, only: ty_fluxes_byband - use mo_optical_props, only: ty_optical_props_1scl,ty_optical_props_2str - use mo_heating_rates, only: compute_heating_rate - use mo_rrtmgp_constants, only: grav, avogad - use module_radsw_parameters, only: topfsw_type, sfcfsw_type, profsw_type, cmpfsw_type - use mo_rrtmgp_sw_cloud_optics, only: rrtmgp_sw_cloud_optics,mcica_subcol_sw - use mo_cloud_sampling, only: sampled_mask_max_ran, sampled_mask_exp_ran, draw_samples - use mersenne_twister, only: random_setseed, random_number, random_stat - use mo_rrtmgp_clr_all_sky, only: rte_sw - - implicit none + use machine, only: kind_phys + use GFS_typedefs, only: GFS_control_type + use mo_rte_kind, only: wl + use mo_gas_optics_rrtmgp, only: ty_gas_optics_rrtmgp + use mo_cloud_optics, only: ty_cloud_optics + use mo_optical_props, only: ty_optical_props_2str + use mo_rrtmgp_clr_all_sky, only: rte_sw + use mo_gas_concentrations, only: ty_gas_concs + use mo_fluxes_byband, only: ty_fluxes_byband + use module_radsw_parameters, only: cmpfsw_type ! Parameters integer,parameter :: nGases = 6 real(kind_phys),parameter :: epsilon=1.0e-6 - real (kind=kind_phys), parameter :: ftiny = 1.0e-12 character(len=3),parameter, dimension(nGases) :: & active_gases = (/ 'h2o', 'co2', 'o3 ', 'n2o', 'ch4', 'o2 '/) - - ! Molecular weight ratios (for converting mmr to vmr) - real(kind_phys), parameter :: & - amd = 28.9644_kind_phys, & ! Molecular weight of dry-air (g/mol) - amw = 18.0154_kind_phys, & ! Molecular weight of water vapor (g/mol) - amo3 = 47.9982_kind_phys, & ! Modelular weight of ozone (g/mol) - amdw = amd/amw, & ! Molecular weight of dry air / water vapor - amdo3 = amd/amo3 ! Molecular weight of dry air / ozone - ! - real (kind_phys), parameter :: & - s0 = 1368.22 ! Solar constant (W/m2) - - ! Logical flags for optional output fields in rrtmgp_sw_run(), default=.false. - logical :: & - l_AllSky_HR_byband = .false., & ! 2D [ncol,nlay] all-sky heating rates, in each band [ncol,nlay,nBandsSW]? - l_ClrSky_HR = .false., & ! 2D [ncol,nlay] clear-sky heating rate? - l_fluxes2D = .false., & ! 2D [ncol,nlay] radiative fluxes *Note* fluxes is a DDT w/ 4 fields. - l_sfcFluxes1D = .false. ! 1D [ncol] surface fluxes *Note* fluxes is a DDT w/ 6 fields. - - ! Module parameters (set during rrtmgp_sw_init()) - integer :: & - rrtmgp_sw_cld_phys, & ! RRTMGP cloud-physics (0-RRTMG, 1-RRTGMP(LUT), 2-RRTMGP(Pade)) - nGptsSW, & ! Number of SW spectral g-points - nBandsSW, & ! Number of SW bands - nrghice, & ! Number of ice roughness categories - ipsdsw0 ! Initial seed for McICA + integer :: nrghice, ipsdsw0 public rrtmgp_sw_init, rrtmgp_sw_run, rrtmgp_sw_finalize contains + !! \section arg_table_rrtmgp_sw_init Argument Table !! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | !! |--------------------|-------------------------------------------------|---------------------------------------------------------------------------|-------|------|----------------------|-----------|--------|----------| @@ -202,26 +168,6 @@ subroutine rrtmgp_sw_init(Model,mpicomm, mpirank, mpiroot, kdist_sw, kdist_cldy_ errmsg = '' errflg = 0 - ! Ensure that requested cloud overlap is reasonable. - if ( iovrsw<0 .or. iovrsw>3 ) then - print *,' *** Error in specification of cloud overlap flag', & - ' IOVRSW=',iovrsw,' in RSWINIT !!' - stop - endif - - ! Check cloud flags for consistency. - if ((icldflg == 0 .and. iswcliq /= 0) .or. & - (icldflg == 1 .and. iswcliq == 0)) then - print *,' *** Model cloud scheme inconsistent with SW', & - ' radiation cloud radiative property setup !!' - stop - endif - if ( isubcsw==0 .and. iovrsw>2 ) then - print *,' *** IOVRSW=',iovrsw,' is not available for ISUBCSW=0 setting!!' - print *,' The program will use maximum/random overlap instead.' - iovrsw = 1 - endif - ! How are we handling cloud-optics? rrtmgp_sw_cld_phys = Model%rrtmgp_cld_phys @@ -497,9 +443,9 @@ subroutine rrtmgp_sw_init(Model,mpicomm, mpirank, mpiroot, kdist_sw, kdist_cldy_ ! Initialize gas concentrations and gas optics class with data do iGas=1,nGases - call check_error_msg(gas_concentrations%set_vmr(active_gases(iGas), 0._kind_phys)) + call check_error_msg('rrtmgp_sw_init',gas_concentrations%set_vmr(active_gases(iGas), 0._kind_phys)) enddo - call check_error_msg(kdist_sw%load(gas_concentrations, gas_names_sw, key_species_sw, band2gpt_sw, & + call check_error_msg('rrtmgp_sw_init',kdist_sw%load(gas_concentrations, gas_names_sw, key_species_sw, band2gpt_sw, & band_lims_sw, press_ref_sw, press_ref_trop_sw, temp_ref_sw, temp_ref_p_sw, temp_ref_t_sw, & vmr_ref_sw, kmajor_sw, kminor_lower_sw, kminor_upper_sw, gas_minor_sw,identifier_minor_sw, & minor_gases_lower_sw, minor_gases_upper_sw, minor_limits_gpt_lower_sw, & @@ -702,578 +648,167 @@ subroutine rrtmgp_sw_init(Model,mpicomm, mpirank, mpiroot, kdist_sw, kdist_cldy_ ! Load tables data for RRTGMP cloud-optics if (rrtmgp_sw_cld_phys .eq. 1) then - call check_error_msg(kdist_cldy_sw%set_ice_roughness(nrghice)) - call check_error_msg(kdist_cldy_sw%load(band_lims_cldy_sw, radliq_lwr_sw, & + call check_error_msg('rrtmgp_sw_init',kdist_cldy_sw%set_ice_roughness(nrghice)) + call check_error_msg('rrtmgp_sw_init',kdist_cldy_sw%load(band_lims_cldy_sw, radliq_lwr_sw, & radliq_upr_sw, radliq_fac_sw, radice_lwr_sw, radice_upr_sw, radice_fac_sw, & lut_extliq_sw, lut_ssaliq_sw, lut_asyliq_sw, lut_extice_sw, lut_ssaice_sw, & lut_asyice_sw)) endif if (rrtmgp_sw_cld_phys .eq. 2) then - call check_error_msg(kdist_cldy_sw%set_ice_roughness(nrghice)) - call check_error_msg(kdist_cldy_sw%load(band_lims_cldy_sw, pade_extliq_sw, & + call check_error_msg('rrtmgp_sw_init',kdist_cldy_sw%set_ice_roughness(nrghice)) + call check_error_msg('rrtmgp_sw_init', kdist_cldy_sw%load(band_lims_cldy_sw, pade_extliq_sw, & pade_ssaliq_sw, pade_asyliq_sw, pade_extice_sw, pade_ssaice_sw, pade_asyice_sw, & pade_sizereg_extliq_sw, pade_sizereg_ssaliq_sw, pade_sizereg_asyliq_sw, & pade_sizereg_extice_sw, pade_sizereg_ssaice_sw, pade_sizereg_asyice_sw)) endif end subroutine rrtmgp_sw_init + ! ######################################################################################### - ! GFS_RRTMGP_SW_RUN ! ######################################################################################### !! \section arg_table_rrtmgp_sw_run Argument Table -!! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | -!! |-----------------|------------------------------------------------------------------------------------------------|---------------------------------------------------------------------------|---------|------|---------------------------|-----------|--------|----------| -!! | p_lay | air_pressure_at_layer_for_radiation_in_hPa | air pressure layer | hPa | 2 | real | kind_phys | in | F | -!! | p_lev | air_pressure_at_interface_for_radiation_in_hPa | air pressure level | hPa | 2 | real | kind_phys | in | F | -!! | t_lay | air_temperature_at_layer_for_radiation | air temperature layer | K | 2 | real | kind_phys | in | F | -!! | t_lev | air_temperature_at_interface_for_radiation | air temperature level | K | 2 | real | kind_phys | in | F | -!! | q_lay | water_vapor_specific_humidity_at_layer_for_radiation | specific humidity layer | kg kg-1 | 2 | real | kind_phys | in | F | -!! | o3_lay | ozone_concentration_at_layer_for_radiation | ozone concentration layer | kg kg-1 | 2 | real | kind_phys | in | F | -!! | vmr_co2 | volume_mixing_ratio_co2 | volume mixing ratio co2 | kg kg-1 | 2 | real | kind_phys | in | F | -!! | vmr_n2o | volume_mixing_ratio_n2o | volume mixing ratio no2 | kg kg-1 | 2 | real | kind_phys | in | F | -!! | vmr_ch4 | volume_mixing_ratio_ch4 | volume mixing ratio ch4 | kg kg-1 | 2 | real | kind_phys | in | F | -!! | vmr_o2 | volume_mixing_ratio_o2 | volume mixing ratio o2 | kg kg-1 | 2 | real | kind_phys | in | F | -!! | vmr_co | volume_mixing_ratio_co | volume mixing ratio co | kg kg-1 | 2 | real | kind_phys | in | F | -!! | vmr_cfc11 | volume_mixing_ratio_cfc11 | volume mixing ratio cfc11 | kg kg-1 | 2 | real | kind_phys | in | F | -!! | vmr_cfc12 | volume_mixing_ratio_cfc12 | volume mixing ratio cfc12 | kg kg-1 | 2 | real | kind_phys | in | F | -!! | vmr_cfc22 | volume_mixing_ratio_cfc22 | volume mixing ratio cfc22 | kg kg-1 | 2 | real | kind_phys | in | F | -!! | vmr_ccl4 | volume_mixing_ratio_ccl4 | volume mixing ratio ccl4 | kg kg-1 | 2 | real | kind_phys | in | F | -!! | icseed | seed_random_numbers_sw | seed for random number generation for shortwave radiation | none | 1 | integer | | in | F | -!! | tau_aer | aerosol_optical_depth_for_longwave_bands_01-16 | aerosol optical depth for shortwave bands 01-16 | none | 3 | real | kind_phys | in | F | -!! | ssa_aer | aerosol_single_scattering_albedo_for_longwave_bands_01-16 | aerosol single scattering albedo for shortwave bands 01-16 | frac | 3 | real | kind_phys | in | F | -!! | asy_aer | aerosol_asymmetry_parameter_for_shortwave_bands_01-16 | aerosol asymmetry paramter for shortwave bands 01-16 | none | 3 | real | kind_phys | in | F | -!! | sfcalb_nir_dir | surface_albedo_due_to_near_IR_direct | surface albedo due to near IR direct beam | frac | 1 | real | kind_phys | in | F | -!! | sfcalb_nir_dif | surface_albedo_due_to_near_IR_diffused | surface albedo due to near IR diffused beam | frac | 1 | real | kind_phys | in | F | -!! | sfcalb_uvis_dir | surface_albedo_due_to_UV_and_VIS_direct | surface albedo due to UV+VIS direct beam | frac | 1 | real | kind_phys | in | F | -!! | sfcalb_uvis_dif | surface_albedo_due_to_UV_and_VIS_diffused | surface albedo due to UV+VIS diffused beam | frac | 1 | real | kind_phys | in | F | -!! | dzlyr | layer_thickness_for_radiation | layer thickness | km | 2 | real | kind_phys | in | F | -!! | delpin | layer_pressure_thickness_for_radiation | layer pressure thickness | hPa | 2 | real | kind_phys | in | F | -!! | de_lgth | cloud_decorrelation_length | cloud decorrelation length | km | 1 | real | kind_phys | in | F | -!! | cossza | cosine_of_zenith_angle | cosine of the solar zenit angle | none | 1 | real | kind_phys | in | F | -!! | solcon | solar_constant | solar constant | W m-2 | 0 | real | kind_phys | in | F | -!! | nday | daytime_points_dimension | daytime points dimension | count | 0 | integer | | in | F | -!! | idxday | daytime_points | daytime points | index | 1 | integer | | in | F | -!! | ncol | horizontal_loop_extent | horizontal dimension | count | 0 | integer | | in | F | -!! | nlay | adjusted_vertical_layer_dimension_for_radiation | number of vertical layers for radiation | count | 0 | integer | | in | F | -!! | lprint | flag_print | flag to print | flag | 0 | logical | | in | F | -!! | cldfrac | total_cloud_fraction | total cloud fraction | frac | 2 | real | kind_phys | in | F | -!! | lsswr | flag_to_calc_sw | flag to calculate SW irradiances | flag | 0 | logical | | in | F | -!! | hswc | tendency_of_air_temperature_due_to_shortwave_heating_on_radiation_time_step | shortwave total sky heating rate | K s-1 | 2 | real | kind_phys | inout | F | -!! | topflx | sw_fluxes_top_atmosphere | shortwave total sky fluxes at the top of the atm | W m-2 | 1 | topfsw_type | | inout | F | -!! | sfcflx | sw_fluxes_sfc | shortwave total sky fluxes at the Earth surface | W m-2 | 1 | sfcfsw_type | | inout | F | -!! | cldtau | cloud_optical_depth_layers_at_0.55mu_band | approx .55mu band layer cloud optical depth | none | 2 | real | kind_phys | inout | F | -!! | hsw0 | tendency_of_air_temperature_due_to_shortwave_heating_assuming_clear_sky_on_radiation_time_step | shortwave clear sky heating rate | K s-1 | 2 | real | kind_phys | inout | T | -!! | hswb | sw_heating_rate_spectral | shortwave total sky heating rate (spectral) | K s-1 | 3 | real | kind_phys | inout | T | -!! | flxprf | sw_fluxes | sw fluxes total sky / csk and up / down at levels | W m-2 | 2 | profsw_type | | inout | T | -!! | fdncmp | components_of_surface_downward_shortwave_fluxes | derived type for special components of surface downward shortwave fluxes | W m-2 | 1 | cmpfsw_type | | inout | T | -!! | cld_lwp | cloud_liquid_water_path | cloud liquid water path | g m-2 | 2 | real | kind_phys | in | T | -!! | cld_ref_liq | mean_effective_radius_for_liquid_cloud | mean effective radius for liquid cloud | micron | 2 | real | kind_phys | in | T | -!! | cld_iwp | cloud_ice_water_path | cloud ice water path | g m-2 | 2 | real | kind_phys | in | T | -!! | cld_ref_ice | mean_effective_radius_for_ice_cloud | mean effective radius for ice cloud | micron | 2 | real | kind_phys | in | T | -!! | cld_rwp | cloud_rain_water_path | cloud rain water path | g m-2 | 2 | real | kind_phys | in | T | -!! | cld_ref_rain | mean_effective_radius_for_rain_drop | mean effective radius for rain drop | micron | 2 | real | kind_phys | in | T | -!! | cld_swp | cloud_snow_water_path | cloud snow water path | g m-2 | 2 | real | kind_phys | in | T | -!! | cld_ref_snow | mean_effective_radius_for_snow_flake | mean effective radius for snow flake | micron | 2 | real | kind_phys | in | T | -!! | cld_od | cloud_optical_depth | cloud optical depth | none | 2 | real | kind_phys | in | T | -!! | cld_ssa | cloud_single_scattering_albedo | cloud single scattering albedo | frac | 2 | real | kind_phys | in | T | -!! | cld_asy | cloud_asymmetry_parameter | cloud asymmetry parameter | none | 2 | real | kind_phys | in | T | -!! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | -!! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | -!! | kdist_sw | K_distribution_file_for_RRTMGP_SW_scheme | DDT containing spectral information for RRTMGP SW radiation scheme | DDT | 0 | ty_gas_optics_rrtmgp | | in | F | -!! | kdist_cldy_sw | K_distribution_file_for_cloudy_RRTMGP_SW_scheme | DDT containing spectral information for cloudy RRTMGP SW radiation scheme | DDT | 0 | ty_cloud_optics | | in | F | +!! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | +!! |-------------------------|------------------------------------------------------------------------------------------------|--------------------------------------------------------------------------|-------|------|-----------------------|-----------|--------|----------| +!! | ncol | horizontal_loop_extent | horizontal dimension | count | 0 | integer | | in | F | +!! | nlay | adjusted_vertical_layer_dimension_for_radiation | number of vertical layers for radiation | count | 0 | integer | | in | F | +!! | p_lay | air_pressure_at_layer_for_radiation_in_hPa | air pressure layer | hPa | 2 | real | kind_phys | in | F | +!! | p_lev | air_pressure_at_interface_for_radiation_in_hPa | air pressure level | hPa | 2 | real | kind_phys | in | F | +!! | t_lay | air_temperature_at_layer_for_radiation | air temperature layer | K | 2 | real | kind_phys | in | F | +!! | kdist_sw | K_distribution_file_for_RRTMGP_SW_scheme | DDT containing spectral information for RRTMGP SW radiation scheme | DDT | 0 | ty_gas_optics_rrtmgp | | in | F | +!! | optical_propsSW_clds | shortwave_optical_properties_for_cloudy_atmosphere | Fortran DDT containing RRTMGP optical properties | DDT | 0 | ty_optical_props_2str | | in | F | +!! | optical_propsSW_aerosol | shortwave_optical_properties_for_aerosols | Fortran DDT containing RRTMGP optical properties | DDT | 0 | ty_optical_props_2str | | in | F | +!! | gas_concentrations | Gas_concentrations_for_RRTMGP_suite_sw | DDT containing gas concentrations for RRTMGP radiation scheme | DDT | 0 | ty_gas_concs | | in | F | +!! | lsswr | flag_to_calc_sw | flag to calculate SW irradiances | flag | 0 | logical | | in | F | +!! | sfcalb_nir_dir | surface_shortwave_albedo_near_infrared_direct_in_each_band | surface sw near-infrared direct albedo in each SW band | frac | 2 | real | kind_phys | out | F | +!! | sfcalb_nir_dif | surface_shortwave_albedo_near_infrared_diffuse_in_each_band | surface sw near-infrared diffuse albedo in each SW band | frac | 2 | real | kind_phys | out | F | +!! | cossza | cosine_of_zenith_angle | cosine of the solar zenit angle | none | 1 | real | kind_phys | in | F | +!! | nday | daytime_points_dimension | daytime points dimension | count | 0 | integer | | in | F | +!! | idxday | daytime_points | daytime points | index | 1 | integer | | in | F | +!! | fluxSW_allsky | sw_flux_profiles_byband_allsky | Fortran DDT containing RRTMGP 3D fluxes | DDT | 0 | ty_fluxes_byband | | out | F | +!! | fluxSW_clrsky | sw_flux_profiles_byband_clrsky | Fortran DDT containing RRTMGP 3D fluxes | DDT | 0 | ty_fluxes_byband | | out | F | +!! | hsw0 | tendency_of_air_temperature_due_to_shortwave_heating_assuming_clear_sky_on_radiation_time_step | shortwave clear sky heating rate | K s-1 | 2 | real | kind_phys | in | T | +!! | hswb | sw_heating_rate_spectral | shortwave total sky heating rate (spectral) | K s-1 | 3 | real | kind_phys | in | T | +!! | scmpsw | components_of_surface_downward_shortwave_fluxes | derived type for special components of surface downward shortwave fluxes | W m-2 | 1 | cmpfsw_type | | inout | F | +!! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | +!! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | !! - subroutine rrtmgp_sw_run(p_lay, p_lev, t_lay, t_lev, q_lay, o3_lay, vmr_co2, vmr_n2o, & ! IN - vmr_ch4, vmr_o2, vmr_co, vmr_cfc11, vmr_cfc12, vmr_cfc22, vmr_ccl4, icseed, tau_aer, & ! IN - ssa_aer, asy_aer, sfcalb_nir_dir, sfcalb_nir_dif, sfcalb_uvis_dir, sfcalb_uvis_dif, & ! IN - dzlyr, delpin, de_lgth, cossza, solcon, nday, idxday, ncol, nlay, lprint, cldfrac, & ! IN - lsswr, kdist_sw, kdist_cldy_sw, & ! IN - hswc, topflx, sfcflx, cldtau, & ! OUT - hsw0, hswB, flxprf, fdncmp, cld_lwp, cld_ref_liq, cld_iwp, cld_ref_ice, cld_rwp, & ! OUT(optional) - cld_ref_rain, cld_swp, cld_ref_snow, cld_od, cld_ssa, cld_asy, & ! OUT(optional) - errmsg, errflg) + subroutine rrtmgp_sw_run(ncol, nlay, kdist_sw, p_lay, t_lay, p_lev, gas_concentrations, & + optical_propsSW_clds, optical_propsSW_aerosol,& + lsswr, sfcalb_nir_dir, sfcalb_nir_dif, cossza, nday, idxday, fluxSW_allsky, fluxSW_clrsky, hsw0, hswb, scmpsw, errmsg, errflg) ! Inputs integer, intent(in) :: & - ncol, & ! Number of horizontal grid-points - nlay, & ! Number of vertical layers - nday ! Number of daytime points - integer, intent(in), dimension(ncol) :: & - icseed ! auxiliary special cloud related array when module - ! variable isubcsw=2, it provides permutation seed - ! for each column profile that are used for generating - ! random numbers. when isubcsw /=2, it will not be used. + ncol, & ! Number of horizontal gridpoints + nlay, & ! Number of vertical layers + nday ! Number of daytime points integer, intent(in), dimension(nday) :: & - idxday ! Index array for daytime points - logical, intent(in) :: & - lprint, & ! Control flag for diagnostics - lsswr ! Flag to calculate RRTMGP SW? + idxday ! Index array for daytime points + real(kind_phys), dimension(ncol,nlay), intent(in) :: & + p_lay, & ! Pressure @ model layer-centers (hPa) + t_lay ! Temperature (K) + real(kind_phys), dimension(ncol,nlay+1), intent(in) :: & + p_lev ! Pressure @ model layer-interfaces (hPa) type(ty_gas_optics_rrtmgp),intent(in) :: & - kdist_sw ! DDT containing LW spectral information - type(ty_cloud_optics),intent(in) :: & - kdist_cldy_sw - !type(ty_gas_concs),intent(inout) :: & - ! gas_concentrations - real(kind_phys), intent(in) :: & - solcon ! Solar constant (W/m2) - real(kind_phys), dimension(ncol), intent(in) :: & + kdist_sw ! DDT containing SW spectral information + real(kind_phys), dimension(kdist_sw%get_nband(),ncol), intent(in) :: & sfcalb_nir_dir, & ! Surface albedo direct (near-IR) (1) - sfcalb_nir_dif, & ! Surface albedo diffuse (near-IR) (1) - sfcalb_uvis_dir, & ! Surface albedo direct (UV and Visible) (1) - sfcalb_uvis_dif, & ! Surface albedo diffuse (UV and Visible) (1) - de_lgth, & ! Cloud decorrelation length (km) + sfcalb_nir_dif ! Surface albedo diffuse (near-IR) (1) + real(kind_phys), dimension(ncol), intent(in) :: & cossza ! Cosine of solar zenith angle (1) - real(kind_phys), dimension(ncol,nlay), intent(in) :: & - dzlyr, & ! layer thinkness (km) - delpin, & ! layer thickness (mb) - cldfrac, & ! Cloud-fraction (1) - p_lay, & ! Pressure @ model layer-centers (mb) - t_lay, & ! Temperature (K) - q_lay, & ! Specific humidity (kg/kg) - o3_lay, & ! O3 mass mixing-ratio (kg/kg) - vmr_co2, & ! Co2 volume-mixing ratio (kg/kg) - vmr_n2o, & ! N2o volume-mixing ratio (kg/kg) - vmr_ch4, & ! Ch4 volume-mixing ratio (kg/kg) - vmr_o2, & ! O2 volume-mixing ratio (kg/kg) - vmr_co, & ! Co volume-mixing ratio (kg/kg) - vmr_cfc11, & ! CFC11 volume-mixing ratio (kg/kg) - vmr_cfc12, & ! CFC12 volume-mixing ratio (kg/kg) - vmr_cfc22, & ! CFC22 volume-mixing ratio (kg/kg) - vmr_ccl4 ! CCl4 volume-mixing ratio (kg/kg) - real(kind_phys), dimension(ncol,nlay+1), intent(in) :: & - p_lev, & ! Pressure @ model layer-interfaces (mb) - t_lev ! Temperature (K) - real(kind_phys), dimension(ncol,nlay,nbandsSW), intent(in) :: & - tau_aer, & ! Aerosol optical depth (1) - ssa_aer, & ! Aerosol single-scattering albedo (1) - asy_aer ! Aerosol asymmetry parameter (1) - ! Inputs (optional) - real(kind_phys), dimension(ncol,nlay), intent(in), optional:: & - cld_lwp, & ! Cloud liquid water path (g/m2) - cld_ref_liq, & ! Effective radius (liquid) (micron) - cld_iwp, & ! Cloud ice water path (g/m2) - cld_ref_ice, & ! Effective radius (ice) (micron) - cld_rwp, & ! Cloud rain water path (g/m2) - cld_ref_rain, & ! Effective radius (rain-drop) (micron) - cld_swp, & ! Cloud snow-water path (g/m2) - cld_ref_snow, & ! Effective radius (snow-flake) (micron) - cld_od, & ! Cloud optical-depth (1) - cld_ssa, & ! Cloud single-scattering albedo (1) - cld_asy ! Cloud asymmetry parameter (1) + type(ty_optical_props_2str),intent(in) :: & + optical_propsSW_clds, & ! RRTMGP DDT: longwave cloud radiative properties + optical_propsSW_aerosol ! RRTMGP DDT: longwave aerosol radiative properties - ! Outputs (mandatory) - character(len=*), intent(out) :: & - errmsg ! Error message - integer, intent(out) :: & - errflg ! Error code - real(kind_phys),dimension(ncol,nlay), intent(inout) :: & - hswc, & ! All-sky heating-rate (K/sec) - cldtau ! ~0.55mu band layer tau (1) - type(topfsw_type), dimension(ncol), intent(inout) :: & - topflx ! radiation fluxes at top, components: - ! upfxc - total sky upward flux at top (w/m2) - ! upfx0 - clear sky upward flux at top (w/m2) - type(sfcfsw_type), dimension(ncol), intent(inout) :: & - sfcflx ! radiation fluxes at sfc, components: - ! upfxc - total sky upward flux at sfc (w/m2) - ! upfx0 - clear sky upward flux at sfc (w/m2) - ! dnfxc - total sky downward flux at sfc (w/m2) - ! dnfx0 - clear sky downward flux at sfc (w/m2) - ! Outputs (optional) - real(kind_phys), dimension(ncol,nlay), intent(inout), optional :: & - hsw0 ! Clear-sky heating-rate (K/sec) - real(kind_phys), dimension(ncol,nlay,nBandsSW), intent(inout), optional :: & - hswb ! All-sky heating rate, in each band (K/sec) - type(profsw_type), dimension(ncol,nlay+1), intent(inout), optional :: & - flxprf ! 2D radiative fluxes, components: - ! upfxc - total sky upward flux (W/m2) - ! dnfxc - total sky dnward flux (W/m2) - ! upfx0 - clear sky upward flux (W/m2) - ! dnfx0 - clear sky dnward flux (W/m2) - type(cmpfsw_type), dimension(ncol), intent(inout), optional :: & - fdncmp ! 2D surface fluxes, components: - ! uvbfc - total sky downward uv-b flux at (W/m2) - ! uvbf0 - clear sky downward uv-b flux at (W/m2) - ! nirbm - downward nir direct beam flux (W/m2) - ! nirdf - downward nir diffused flux (W/m2) - ! visbm - downward uv+vis direct beam flux(W/m2) - ! visdf - downward uv+vis diffused flux (W/m2) + type(ty_gas_concs),intent(in) :: & + gas_concentrations ! RRTMGP DDT: trace gas concentrations (vmr) + logical, intent(in) :: & + lsswr ! Flag to calculate SW irradiances - ! RTE+RRTMGP classes - type(ty_gas_concs) :: & - gas_concentrations - type(ty_optical_props_2str) :: & - optical_props_clr, & ! Optical properties for gaseous atmosphere - optical_props_aer, & ! Optical properties for aerosols - optical_props_mcica, & ! Optical properties for clouds (sampled) - optical_props_cldy ! Optical properties for clouds (by-band) - type(ty_fluxes_byband) :: & - fluxAllSky, & ! All-sky flux (W/m2) - fluxClrSky ! Clear-sky flux (W/m2) + ! Outputs + character(len=*), intent(out) :: errmsg + integer, intent(out) :: errflg + type(ty_fluxes_byband),intent(out) :: & + fluxSW_allsky, & ! All-sky flux (W/m2) + fluxSW_clrsky ! Clear-sky flux (W/m2) + + ! Inputs (optional) (NOTE. We only need the optional arguments to know what fluxes to output, HR's are computed later) + real(kind_phys), dimension(ncol,nlay), optional, intent(in) :: & + hsw0 ! Clear-sky heating rate (K/sec) + real(kind_phys), dimension(ncol,nlay,kdist_sw%get_nband()), intent(in), optional :: & + hswb ! All-sky heating rate, by band (K/sec) + ! Outputs (optional) + type(cmpfsw_type), dimension(ncol), intent(out),optional :: & + scmpsw ! 2D surface fluxes, components: + ! uvbfc - total sky downward uv-b flux at (W/m2) + ! uvbf0 - clear sky downward uv-b flux at (W/m2) + ! nirbm - downward nir direct beam flux (W/m2) + ! nirdf - downward nir diffused flux (W/m2) + ! visbm - downward uv+vis direct beam flux (W/m2) + ! visdf - downward uv+vis diffused flux (W/m2) - ! Types used by Random Number Generator - type(random_stat) :: rng_stat ! Local variables - integer :: iCol, iBand, iGpt, iDay, iLay, iTOA, iSFC - integer,dimension(ncol) :: ipseed - real(kind_phys) :: cfrac, asyw, ssaw, za1, za2 - logical :: top_at_1=.false. - real(kind_phys), dimension(ncol) :: clrfracSFC, cldfracSFC - real(kind_phys), dimension(ncol,nlay) :: vmr_o3, vmr_h2o, coldry, tem0, & - cld_ref_liq2,cld_ref_ice2 - real(kind_phys), dimension(ncol,nlay,nBandsSW) :: thetaTendByBandAllSky - real(kind_phys), dimension(nday,nlay) :: cld_lwp2,thetaTendClrSky, & - thetaTendAllSky real(kind_phys), dimension(nday,nlay+1),target :: & - flux_up_allSky, flux_up_clrSky, flux_dn_allSky, flux_dn_clrSky - real(kind_phys), dimension(nday,nlay+1,nBandsSW),target :: & - fluxBB_up_allSky, fluxBB_dn_allSky - real(kind_phys), dimension(nday,nGptsSW) :: toa_flux - real(kind_phys), dimension(nday,nlay,nBandsSW) :: tau_cld, asy_cld, ssa_cld - real(kind_phys), dimension(nGptsSW,nlay,ncol) :: & - rng3D - real(kind_phys), dimension(nGptsSW*nLay) :: & - rng1D - logical,dimension(ncol,nlay) :: & - liqmask,icemask - logical, dimension(ncol,nlay,nGptsSW) :: & - cldfracMCICA + fluxSW_up_allsky, fluxSW_up_clrsky, fluxSW_dn_allsky, fluxSW_dn_clrsky + real(kind_phys), dimension(nday,nlay+1,kdist_sw%get_nband()),target :: & + fluxSWBB_up_allsky, fluxSWBB_dn_allsky + logical :: l_ClrSky_HR=.false., l_AllSky_HR_byband=.false., l_scmpsw=.false. - ! Initialize + ! Initialize CCPP error handling variables errmsg = '' errflg = 0 if (.not. lsswr) return - if (nday <= 0) return - ! Are any optional outputs requested? + ! Are any optional outputs requested? Need to know now to compute correct fluxes. l_ClrSky_HR = present(hsw0) l_AllSky_HR_byband = present(hswb) - l_fluxes2D = present(flxprf) - l_sfcFluxes1D = present(fdncmp) - - ! Check for optional input arguments, this depends on cloud method - if (iswcliq > 0) then ! use prognostic cloud method - if ( .not.present(cld_lwp) .or. .not.present(cld_ref_liq) .or. & - .not.present(cld_iwp) .or. .not.present(cld_ref_ice) .or. & - .not.present(cld_rwp) .or. .not.present(cld_ref_rain) .or. & - .not.present(cld_swp) .or. .not.present(cld_ref_snow) )then - write(errmsg,'(*(a))') & - 'Logic error: iswcliq>0 requires the following', & - ' optional arguments to be present:', & - ' cld_lwp, cld_ref_liq, cld_iwp, cld_ref_ice,', & - ' cld_rwp, cld_ref_rain, cld_swp, cld_ref_snow' - errflg = 1 - return - end if - else ! use diagnostic cloud method - if ( .not.present(cld_od) .or. .not.present(cld_ssa) .or. & - .not.present(cld_asy)) then - write(errmsg,'(*(a))') & - 'Logic error: iswcliq<=0 requires the following', & - ' optional arguments to be present:', & - ' cld_od, cld_ssa, cld_asy' - errflg = 1 - return - end if - endif - - ! What is vertical ordering? - top_at_1 = (p_lay(1,1) .lt. p_lay(1,nlay)) - if (top_at_1) then - iSFC = nlay+1 - iTOA = 1 - else - iSFC = 1 - iTOA = nlay+1 - endif - - ! Change random number seed value for each radiation invocation (isubcsw =1 or 2). - if(isubcsw == 1) then ! advance prescribed permutation seed - do iCol = 1, ncol - ipseed(iCol) = ipsdsw0 + iCol - enddo - elseif (isubcsw == 2) then ! use input array of permutaion seeds - do iCol = 1, ncol - ipseed(iCol) = icseed(iCol) - enddo + l_scmpsw = present(scmpsw) + if ( l_scmpsw ) then + scmpsw = cmpfsw_type (0., 0., 0., 0., 0., 0.) endif - ! Compute volume mixing-ratios for ozone (mmr) and specific-humidity. - vmr_h2o = merge((q_lay/(1-q_lay))*amdw, 0., q_lay .ne. 1.) - vmr_o3 = merge(o3_lay*amdo3, 0., o3_lay .gt. 0.) - - ! Compute ice/liquid cloud masks, needed by rrtmgp_cloud_optics - liqmask = (cldfrac .gt. 0 .and. cld_lwp .gt. 0) - icemask = (cldfrac .gt. 0 .and. cld_iwp .gt. 0) - - ! RRTMGP cloud_optics expects particle size to be in a certain range. bound here - if (rrtmgp_sw_cld_phys .gt. 0) then - cld_ref_ice2 = cld_ref_ice - where(cld_ref_ice2 .gt. kdist_cldy_sw%get_max_radius_ice()) cld_ref_ice2=kdist_cldy_sw%get_max_radius_ice() - where(cld_ref_ice2 .lt. kdist_cldy_sw%get_min_radius_ice()) cld_ref_ice2=kdist_cldy_sw%get_min_radius_ice() - cld_ref_liq2 = cld_ref_liq - where(cld_ref_liq2 .gt. kdist_cldy_sw%get_max_radius_liq()) cld_ref_liq2=kdist_cldy_sw%get_max_radius_liq() - where(cld_ref_liq2 .lt. kdist_cldy_sw%get_min_radius_liq()) cld_ref_liq2=kdist_cldy_sw%get_min_radius_liq() - endif - - ! Compute dry air column amount - tem0 = (1. - vmr_h2o)*amd + vmr_h2o*amw - coldry = ( 1.0e-20 * 1.0e3 *avogad)*delpin / (100.*grav*tem0*(1. + vmr_h2o)) - - ! Compute fractions of clear sky view at surface. *NOTE* This is only used if cloud radiative - ! properties are provided directly. - clrfracSFC = 1._kind_phys - if (iovrsw == 0) then ! random overlapping - do iCol=1,nCol - do iLay = 1, nlay - clrfracSFC(iCol) = clrfracSFC(iCol) * (1._kind_phys - cldfrac(iCol,iLay)) - enddo - enddo - else if (iovrsw == 1) then ! max/ran overlapping - do iLay = 1, nlay - if (cldfrac(iCol,iLay) > ftiny) then ! cloudy layer - cldfracSFC(iCol) = min ( cldfracSFC(iCol), 1._kind_phys-cldfrac(iCol,iLay) ) - elseif (cldfracSFC(iCol) < 1._kind_phys) then ! clear layer - clrfracSFC(iCol) = clrfracSFC(iCol) * cldfracSFC(iCol) - cldfracSFC(iCol) = 1._kind_phys - endif - enddo - clrfracSFC(iCol) = clrfracSFC(iCol) * cldfracSFC(iCol) - else if (iovrsw >= 2) then - do iLay = 1, nlay - clrfracSFC(iCol) = min ( clrfracSFC(iCol), 1._kind_phys-cldfrac(iCol,iLay) ) ! used only as clear/cloudy indicator - enddo - endif - if (clrfracSFC(iCol) <= ftiny) clrfracSFC(iCol) = 0._kind_phys - if (clrfracSFC(iCol) > (1._kind_phys-epsilon)) clrfracSFC(iCol) = 1._kind_phys - cldfracSFC(iCol) = 1._kind_phys - clrfracSFC(iCol) - - ! Initialize outputs - hswc(:,:) = 0. - cldtau(:,:) = 0. - topflx = topfsw_type ( 0., 0., 0. ) - sfcflx = sfcfsw_type ( 0., 0., 0., 0. ) - if (l_ClrSky_HR) then - hsw0(:,:) = 0. - endif - if(l_AllSky_HR_byband) then - hswb(:,:,:) = 0. - endif - if (l_fluxes2D) then - flxprf = profsw_type ( 0., 0., 0., 0. ) - endif - if (l_sfcFluxes1D) then - fdncmp = cmpfsw_type (0.,0.,0.,0.,0.,0.) - endif - - ! ####################################################################################### - ! CALL RRTMGP (Only for daylit (idxday) points) - ! ####################################################################################### if (nDay .gt. 0) then - - ! Allocate space for gas optical properties - ! Cloud optics [nCol,nLay,nBands] - call check_error_msg(optical_props_cldy%init(kdist_sw%get_band_lims_wavenumber())) - call check_error_msg(optical_props_cldy%alloc_2str(ncol,nlay)) - ! Aerosol optics [Ccol,nLay,nBands] - call check_error_msg(optical_props_aer%init(kdist_sw%get_band_lims_wavenumber())) - call check_error_msg(optical_props_aer%alloc_2str(ncol,nlay)) - ! Cloud optics sampled [nCol,nLay,nGpts] - call check_error_msg(optical_props_mcica%alloc_2str(ncol, nlay, kdist_sw)) - - ! Initialize RRTMGP DDT containing 2D(3D) fluxes - fluxAllSky%flux_up => flux_up_allSky - fluxAllsky%flux_dn => flux_dn_allSky - fluxClrSky%flux_up => flux_up_clrSky - fluxClrsky%flux_dn => flux_dn_clrSky + fluxSW_allsky%flux_up => fluxSW_up_allsky + fluxSW_allsky%flux_dn => fluxSW_dn_allsky + fluxSW_clrsky%flux_up => fluxSW_up_clrsky + fluxSW_clrsky%flux_dn => fluxSW_dn_clrsky ! Only calculate fluxes by-band, only when heating-rate profiles by band are requested. if (l_AllSky_HR_byband) then - fluxAllSky%bnd_flux_up => fluxBB_up_allSky - fluxAllsky%bnd_flux_dn => fluxBB_dn_allSky - endif - - ! ####################################################################################### - ! Set gas concentrations - ! ####################################################################################### - call gas_concentrations%reset() - call check_error_msg(gas_concentrations%set_vmr('o2', vmr_o2(idxday,1:nlay))) - call check_error_msg(gas_concentrations%set_vmr('co2', vmr_co2(idxday,1:nlay))) - call check_error_msg(gas_concentrations%set_vmr('ch4', vmr_ch4(idxday,1:nlay))) - call check_error_msg(gas_concentrations%set_vmr('n2o', vmr_n2o(idxday,1:nlay))) - call check_error_msg(gas_concentrations%set_vmr('h2o', vmr_h2o(idxday,1:nlay))) - call check_error_msg(gas_concentrations%set_vmr('o3', vmr_o3(idxday,1:nlay))) - - ! ####################################################################################### - ! Copy aerosol to RRTMG DDT - ! ####################################################################################### - optical_props_aer%tau = tau_aer(idxday,:,:) - optical_props_aer%ssa = ssa_aer(idxday,:,:) - optical_props_aer%g = asy_aer(idxday,:,:) - - ! #################################################################################### - ! Compute cloud-optics for RTE. - ! #################################################################################### - - ! Compute in-cloud radiative properties. - if (any(cldfrac(idxday,:) .gt. 0)) then - ! i) RRTMG cloud optics. - ! Cloud-optical properties by type provided. Compute optical-depth, single- - ! scattering albedo, and asymmetry parameter - if (rrtmgp_sw_cld_phys .eq. 0) then - if (.not. present(cld_od)) then - call rrtmgp_sw_cloud_optics(nday, nlay, nBandsSW, & - cld_lwp(idxday,1:nLay), cld_ref_liq(idxday,1:nLay), & - cld_iwp(idxday,1:nLay), cld_ref_ice(idxday,1:nLay), & - cld_rwp(idxday,1:nLay), cld_ref_rain(idxday,1:nLay), & - cld_swp(idxday,1:nLay), cld_ref_snow(idxday,1:nLay), & - cldfrac(idxday,1:nLay), & - tau_cld, ssa_cld, asy_cld) - optical_props_cldy%tau = tau_cld - optical_props_cldy%ssa = ssa_cld - optical_props_cldy%g = asy_cld - else - ! Cloud-optical depth, single scattering albedo, and asymmetry parameter provided. - do iDay=1,nDay - do iLay=1,nLay - if (cldfrac(iCol,iLay) .gt. 1e-20_kind_phys) then - optical_props_cldy%tau(iDay,iLay,:) = cld_od(idxday(iDay),iLay) - optical_props_cldy%ssa(iDay,iLay,:) = cld_ssa(idxday(iDay),iLay) - optical_props_cldy%g(iDay,iLay,:) = cld_asy(idxday(iDay),iLay) - else - optical_props_cldy%tau(iDay,iLay,:) = 0. - optical_props_cldy%ssa(iDay,iLay,:) = 1. - optical_props_cldy%g(iDay,iLay,:) = 0. - endif - end do - end do - endif - endif - - ! ii) Use RRTMGP cloud-optics. - if (rrtmgp_sw_cld_phys .gt. 0) then - call check_error_msg(kdist_cldy_sw%cloud_optics(nday, nlay, nBandsSW, nrghice, & - liqmask(idxday,1:nLay), icemask(idxday,1:nLay), cld_lwp(idxday,1:nLay), & - cld_iwp(idxday,1:nLay), cld_ref_liq2(idxday,1:nLay), & - cld_ref_ice2(idxday,1:nLay), optical_props_cldy)) - end if - endif - - ! ####################################################################################### - ! Call McICA to sample clouds. - ! ####################################################################################### - if (isubcsw .gt. 0) then - ! Call RNG. Mersennse Twister accepts 1D array, so loop over columns and collapse along G-points - ! and layers. ([nGpts,nLayer,nColumn]-> [nGpts*nLayer]*nColumn) - do iCol=1,nCol - call random_setseed(ipseed(icol),rng_stat) - call random_number(rng1D,rng_stat) - rng3D(:,:,iCol) = reshape(source = rng1D,shape=[nGptsSW,nLay]) - enddo - - ! Call McICA - select case ( iovrsw ) - ! Maximumn-random - case(1) - call check_error_msg(sampled_mask_max_ran(rng3D,cldfrac,cldfracMCICA)) - end select - - ! Map band optical depth to each g-point using McICA - call check_error_msg(draw_samples(cldfracMCICA,optical_props_cldy,optical_props_mcica)) - endif - - ! ####################################################################################### - ! Compute fluxes - ! ####################################################################################### - print*,'In rrmtgp_sw(): ' - print*,' shape(optical_props_aerosol%tau): ',shape(optical_props_aer%tau) - print*,' shape(optical_props_clds%tau): ',shape(optical_props_mcica%tau) - call check_error_msg(rte_sw( & - kdist_sw, & - gas_concentrations, & - p_lay(idxday,1:nlay), & - t_lay(idxday,1:nlay), & - p_lev(idxday,1:nlay+1), & - cossza(idxday), & - spread(sfcalb_nir_dir(idxday),1, ncopies = nBandsSW), & - spread(sfcalb_nir_dif(idxday),1, ncopies = nBandsSW), & - optical_props_mcica, & - fluxAllSky, & - fluxClrSky, & - aer_props = optical_props_aer)) - - ! ####################################################################################### - ! Compute heating rates - ! ####################################################################################### - if (l_ClrSky_HR) then - call check_error_msg(compute_heating_rate( & - fluxClrSky%flux_up, & - fluxClrSky%flux_dn, & - p_lev(idxday,1:nlay+1), & - thetaTendClrSky)) - endif - if (l_AllSky_HR_byband) then - do iBand=1,nBandsSW - call check_error_msg(compute_heating_rate( & - fluxAllSky%bnd_flux_up(:,:,iBand), & - fluxAllSky%bnd_flux_dn(:,:,iBand), & - p_lev(idxday,1:nlay+1), & - thetaTendByBandAllSky(:,:,iBand))) - enddo - else - call check_error_msg(compute_heating_rate( & - fluxAllSky%flux_up, & - fluxAllSky%flux_dn, & - p_lev(idxday,1:nlay+1), & - thetaTendAllSky)) + fluxSW_allsky%bnd_flux_up => fluxSWBB_up_allsky + fluxSW_allsky%bnd_flux_dn => fluxSWBB_dn_allsky endif - - end if ! Daylit days - - ! ####################################################################################### - ! Copy fluxes from RRTGMP types into model radiation types. - ! ####################################################################################### - ! Mandatory outputs - topflx(idxday)%upfxc = fluxAllSky%flux_up(:,iTOA) - topflx(idxday)%upfx0 = fluxClrSky%flux_up(:,iTOA) - sfcflx(idxday)%upfxc = fluxAllSky%flux_up(:,iSFC) - sfcflx(idxday)%upfx0 = fluxClrSky%flux_up(:,iSFC) - sfcflx(idxday)%dnfxc = fluxAllSky%flux_dn(:,iSFC) - sfcflx(idxday)%dnfx0 = fluxClrSky%flux_dn(:,iSFC) - cldtau(idxday,:) = optical_props_cldy%tau(:,:,10) - hswc(idxday,:) = thetaTendAllSky - - ! Optional output - if(l_fluxes2D) then - flxprf(idxday,:)%upfxc = fluxAllSky%flux_up - flxprf(idxday,:)%dnfxc = fluxAllSky%flux_dn - flxprf(idxday,:)%upfx0 = fluxClrSky%flux_up - flxprf(idxday,:)%dnfx0 = fluxClrSky%flux_dn - endif - if (l_AllSky_HR_byband) then - hswb(idxday,:,:) = thetaTendByBandAllSky - endif - if (l_ClrSky_HR) then - hsw0(idxday,:) = thetaTendClrSky + + ! Call RRTMGP SW scheme + call check_error_msg('rrtmgp_sw_run',rte_sw( & + kdist_sw, & ! IN - spectral information + gas_concentrations, & ! IN - gas concentrations (vmr) + p_lay(idxday,1:nlay), & ! IN - pressure at layer interfaces (Pa) + t_lay(idxday,1:nlay), & ! IN - temperature at layer interfaes (K) + p_lev(idxday,1:nlay+1), & ! IN - pressure at layer centers (Pa) + cossza(idxday), & ! IN - Cosine of solar zenith angle + sfcalb_nir_dir(:,idxday), & ! IN - Shortwave surface albedo (direct) + sfcalb_nir_dif(:,idxday), & ! IN - Shortwave surface albedo (diffuse) + optical_propsSW_clds, & ! IN - DDT containing cloud optical information + fluxSW_allsky, & ! OUT - Fluxes, all-sky, 3D (nCol,nLay,nBand) + fluxSW_clrsky, & ! OUT - Fluxes, clear-sky, 3D (nCol,nLay,nBand) + aer_props = optical_propsSW_aerosol)) ! IN(optional) - DDT containing aerosol optical information endif + end subroutine rrtmgp_sw_run - ! ######################################################################################### - ! ######################################################################################### + subroutine rrtmgp_sw_finalize() end subroutine rrtmgp_sw_finalize - - ! ######################################################################################### - ! Ancillary functions - ! ######################################################################################### - subroutine check_error_msg(error_msg) - character(len=*), intent(in) :: error_msg + subroutine check_error_msg(routine_name, error_msg) + character(len=*), intent(in) :: & + error_msg, routine_name if(error_msg /= "") then - print*,"ERROR(rrtmgp_sw.F90): " + print*,"ERROR("//trim(routine_name)//"): " print*,trim(error_msg) return end if - end subroutine check_error_msg - ! ######################################################################################### - ! ######################################################################################### + end subroutine check_error_msg + + end module rrtmgp_sw