From 949ca840e63edd9dc43055d90f83c527b4d44e19 Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Sat, 9 Oct 2021 20:07:58 -0400 Subject: [PATCH 01/15] reverting to my versions of GFS_rrtmgp_*post.F90 files --- physics/GFS_rrtmgp_lw_post.F90 | 72 +++++++++++++------------- physics/GFS_rrtmgp_sw_post.F90 | 92 ++++++++++++++++++---------------- physics/rte-rrtmgp | 2 +- 3 files changed, 85 insertions(+), 81 deletions(-) diff --git a/physics/GFS_rrtmgp_lw_post.F90 b/physics/GFS_rrtmgp_lw_post.F90 index ff0346fe4..2e30fdbd0 100644 --- a/physics/GFS_rrtmgp_lw_post.F90 +++ b/physics/GFS_rrtmgp_lw_post.F90 @@ -1,4 +1,4 @@ -module GFS_rrtmgp_lw_post +module GFS_rrtmgp_lw_post use machine, only: kind_phys use module_radiation_aerosols, only: NSPC1 use module_radlw_parameters, only: topflw_type, sfcflw_type, proflw_type @@ -6,9 +6,9 @@ module GFS_rrtmgp_lw_post 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 radiation_tools, only: check_error_msg + use radiation_tools, only: check_error_msg implicit none - + public GFS_rrtmgp_lw_post_init,GFS_rrtmgp_lw_post_run,GFS_rrtmgp_lw_post_finalize contains @@ -29,16 +29,16 @@ subroutine GFS_rrtmgp_lw_post_run (nCol, nLev, lslwr, do_lw_clrsky_hr, save_diag fluxlwDOWN_clrsky, raddt, aerodp, cldsa, mtopa, mbota, cld_frac, cldtaulw, fluxr, & sfcdlw, sfculw, sfcflw, tsflw, htrlw, topflw, flxprf_lw, htrlwc, errmsg, errflg) - ! Inputs - integer, intent(in) :: & + ! Inputs + integer, intent(in) :: & nCol, & ! Horizontal loop extent nLev ! Number of vertical layers - logical, intent(in) :: & + logical, intent(in) :: & lslwr, & ! Logical flags for lw radiation calls - do_lw_clrsky_hr, & ! Output clear-sky SW heating-rate? - save_diag ! Output radiation diagnostics? + do_lw_clrsky_hr, & ! Output clear-sky SW heating-rate? + save_diag ! Output radiation diagnostics? real(kind_phys), intent(in) :: & - fhlwr ! Frequency for SW radiation + fhlwr ! Frequency for SW radiation real(kind_phys), dimension(nCol), intent(in) :: & tsfa ! Lowest model layer air temperature for radiation (K) real(kind_phys), dimension(nCol, nLev), intent(in) :: & @@ -52,25 +52,25 @@ subroutine GFS_rrtmgp_lw_post_run (nCol, nLev, lslwr, do_lw_clrsky_hr, save_diag real(kind_phys), intent(in) :: & raddt ! Radiation time step real(kind_phys), dimension(nCol,NSPC1), intent(in) :: & - aerodp ! Vertical integrated optical depth for various aerosol species + aerodp ! Vertical integrated optical depth for various aerosol species real(kind_phys), dimension(nCol,5), intent(in) :: & - cldsa ! Fraction of clouds for low, middle, high, total and BL + cldsa ! Fraction of clouds for low, middle, high, total and BL integer, dimension(nCol,3), intent(in) ::& mbota, & ! vertical indices for low, middle and high cloud tops mtopa ! vertical indices for low, middle and high cloud bases real(kind_phys), dimension(nCol,nLev), intent(in) :: & cld_frac, & ! Total cloud fraction in each layer - cldtaulw ! approx 10.mu band layer cloud optical depth - + cldtaulw ! approx 10.mu band layer cloud optical depth + real(kind=kind_phys), dimension(:,:), intent(inout) :: fluxr - + ! Outputs (mandatory) real(kind_phys), dimension(nCol), intent(inout) :: & sfcdlw, & ! Total sky sfc downward lw flux (W/m2) sfculw, & ! Total sky sfc upward lw flux (W/m2) tsflw ! surface air temp during lw calculation (K) type(sfcflw_type), dimension(nCol), intent(inout) :: & - sfcflw ! LW radiation fluxes at sfc + sfcflw ! LW radiation fluxes at sfc real(kind_phys), dimension(nCol,nLev), intent(inout) :: & htrlw ! LW all-sky heating rate type(topflw_type), dimension(nCol), intent(out) :: & @@ -79,7 +79,7 @@ subroutine GFS_rrtmgp_lw_post_run (nCol, nLev, lslwr, do_lw_clrsky_hr, save_diag errmsg integer, intent(out) :: & errflg - + ! Outputs (optional) type(proflw_type), dimension(nCol, nLev+1), optional, intent(inout) :: & flxprf_lw ! 2D radiative fluxes, components: @@ -89,7 +89,7 @@ subroutine GFS_rrtmgp_lw_post_run (nCol, nLev, lslwr, do_lw_clrsky_hr, save_diag ! dnfx0 - clear sky dnward flux (W/m2) real(kind_phys),dimension(nCol, nLev),intent(inout),optional :: & htrlwc ! Longwave clear-sky heating-rate (K/sec) - + ! Local variables integer :: i, j, k, iSFC, iTOA, itop, ibtc logical :: l_fluxeslw2d, top_at_1 @@ -118,7 +118,7 @@ subroutine GFS_rrtmgp_lw_post_run (nCol, nLev, lslwr, do_lw_clrsky_hr, save_diag endif ! ####################################################################################### - ! Compute LW heating-rates. + ! Compute LW heating-rates. ! ####################################################################################### ! Clear-sky heating-rate (optional) if (do_lw_clrsky_hr) then @@ -128,7 +128,7 @@ subroutine GFS_rrtmgp_lw_post_run (nCol, nLev, lslwr, do_lw_clrsky_hr, save_diag p_lev, & ! IN - Pressure @ layer-interfaces (Pa) htrlwc)) ! OUT - Longwave clear-sky heating rate (K/sec) endif - + ! All-sky heating-rate (mandatory) call check_error_msg('GFS_rrtmgp_post',compute_heating_rate( & fluxlwUP_allsky, & ! IN - RRTMGP upward longwave all-sky flux profiles (W/m2) @@ -147,7 +147,7 @@ subroutine GFS_rrtmgp_lw_post_run (nCol, nLev, lslwr, do_lw_clrsky_hr, save_diag sfcflw(:)%upfx0 = fluxlwUP_clrsky(:,iSFC) sfcflw(:)%dnfxc = fluxlwDOWN_allsky(:,iSFC) sfcflw(:)%dnfx0 = fluxlwDOWN_clrsky(:,iSFC) - + ! Optional outputs if(l_fluxeslw2d) then flxprf_lw%upfxc = fluxlwUP_allsky @@ -155,7 +155,7 @@ subroutine GFS_rrtmgp_lw_post_run (nCol, nLev, lslwr, do_lw_clrsky_hr, save_diag flxprf_lw%upfx0 = fluxlwUP_clrsky flxprf_lw%dnfx0 = fluxlwDOWN_clrsky endif - + ! Save surface air temp for diurnal adjustment at model t-steps tsflw (:) = tsfa(:) @@ -165,8 +165,8 @@ subroutine GFS_rrtmgp_lw_post_run (nCol, nLev, lslwr, do_lw_clrsky_hr, save_diag ! ####################################################################################### ! Save LW diagnostics - ! - 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 + ! - 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 @@ -182,24 +182,24 @@ subroutine GFS_rrtmgp_lw_post_run (nCol, nLev, lslwr, do_lw_clrsky_hr, save_diag fluxr(i,30) = fluxr(i,30) + fhlwr * fluxlwDOWN_clrsky(i,iSFC) ! clear sky sfc lw dn fluxr(i,33) = fluxr(i,33) + fhlwr * fluxlwUP_clrsky( i,iSFC) ! clear sky sfc lw up enddo - - do i=1,nCol - fluxr(i,17) = fluxr(i,17) + raddt * cldsa(i,4) - fluxr(i,18) = fluxr(i,18) + raddt * cldsa(i,5) - enddo + +! do i=1,nCol +! fluxr(i,17) = fluxr(i,17) + raddt * cldsa(i,4) +! fluxr(i,18) = fluxr(i,18) + raddt * cldsa(i,5) +! enddo ! Save cld frac,toplyr,botlyr and top temp, note that the order of h,m,l cloud is reversed for ! the fluxr output. save interface pressure (pa) of top/bot do j = 1, 3 do i = 1, nCol - tem0d = raddt * cldsa(i,j) - itop = mtopa(i,j) - ibtc = mbota(i,j) - fluxr(i, 8-j) = fluxr(i, 8-j) + tem0d - fluxr(i,11-j) = fluxr(i,11-j) + tem0d * p_lev(i,itop) - fluxr(i,14-j) = fluxr(i,14-j) + tem0d * p_lev(i,ibtc) - fluxr(i,17-j) = fluxr(i,17-j) + tem0d * t_lay(i,itop) - +! tem0d = raddt * cldsa(i,j) +! itop = mtopa(i,j) +! ibtc = mbota(i,j) +! fluxr(i, 8-j) = fluxr(i, 8-j) + tem0d +! fluxr(i,11-j) = fluxr(i,11-j) + tem0d * p_lev(i,itop) +! fluxr(i,14-j) = fluxr(i,14-j) + tem0d * p_lev(i,ibtc) +! fluxr(i,17-j) = fluxr(i,17-j) + tem0d * t_lay(i,itop) + ! Add optical depth and emissivity output tem2 = 0. do k=ibtc,itop diff --git a/physics/GFS_rrtmgp_sw_post.F90 b/physics/GFS_rrtmgp_sw_post.F90 index 23a681826..d8514650a 100644 --- a/physics/GFS_rrtmgp_sw_post.F90 +++ b/physics/GFS_rrtmgp_sw_post.F90 @@ -1,14 +1,14 @@ -module GFS_rrtmgp_sw_post +module GFS_rrtmgp_sw_post use machine, only: kind_phys use module_radiation_aerosols, only: NSPC1 use module_radsw_parameters, only: topfsw_type, sfcfsw_type, profsw_type, cmpfsw_type 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 radiation_tools, only: check_error_msg + use radiation_tools, only: check_error_msg use rrtmgp_sw_gas_optics, only: sw_gas_props implicit none - + public GFS_rrtmgp_sw_post_init,GFS_rrtmgp_sw_post_run,GFS_rrtmgp_sw_post_finalize contains @@ -33,23 +33,23 @@ subroutine GFS_rrtmgp_sw_post_run (nCol, nLev, nDay, idxday, lsswr, do_sw_clrsky nirbmdi, nirdfdi, visbmdi, visdfdi, nirbmui, nirdfui, visbmui, visdfui, sfcnsw, & sfcdsw, htrsw, sfcfsw, topfsw, htrswc, flxprf_sw, scmpsw, errmsg, errflg) - ! Inputs - integer, intent(in) :: & - nCol, & ! Horizontal loop extent + ! Inputs + integer, intent(in) :: & + nCol, & ! Horizontal loop extent nLev, & ! Number of vertical layers nDay ! Number of daylit columns integer, intent(in), dimension(nday) :: & idxday ! Index array for daytime points - logical, intent(in) :: & - lsswr, & ! Call SW radiation? - do_sw_clrsky_hr, & ! Output clear-sky SW heating-rate? - save_diag ! Output radiation diagnostics? + logical, intent(in) :: & + lsswr, & ! Call SW radiation? + do_sw_clrsky_hr, & ! Output clear-sky SW heating-rate? + save_diag ! Output radiation diagnostics? real(kind_phys), intent(in) :: & fhswr ! Frequency for SW radiation real(kind_phys), dimension(nCol), intent(in) :: & t_lay, & ! Temperature at model layer centers (K) - coszen, & ! Cosine(SZA) - coszdg ! Cosine(SZA), daytime + coszen, & ! Cosine(SZA) + coszdg ! Cosine(SZA), daytime real(kind_phys), dimension(nCol, nLev+1), intent(in) :: & p_lev ! Pressure @ model layer-interfaces (Pa) real(kind_phys), dimension(sw_gas_props%get_nband(),ncol), intent(in) :: & @@ -65,17 +65,17 @@ subroutine GFS_rrtmgp_sw_post_run (nCol, nLev, nDay, idxday, lsswr, do_sw_clrsky real(kind_phys), intent(in) :: & raddt ! Radiation time step real(kind_phys), dimension(nCol,NSPC1), intent(in) :: & - aerodp ! Vertical integrated optical depth for various aerosol species + aerodp ! Vertical integrated optical depth for various aerosol species real(kind_phys), dimension(nCol,5), intent(in) :: & - cldsa ! Fraction of clouds for low, middle, high, total and BL + cldsa ! Fraction of clouds for low, middle, high, total and BL integer, dimension(nCol,3), intent(in) ::& mbota, & ! vertical indices for low, middle and high cloud tops mtopa ! vertical indices for low, middle and high cloud bases real(kind_phys), dimension(nCol,nLev), intent(in) :: & cld_frac, & ! Total cloud fraction in each layer cldtausw ! approx .55mu band layer cloud optical depth - - ! Inputs (optional) + + ! Inputs (optional) type(cmpfsw_type), dimension(nCol), intent(inout), optional :: & scmpsw ! 2D surface fluxes, components: ! uvbfc - total sky downward uv-b flux at (W/m2) @@ -83,10 +83,10 @@ subroutine GFS_rrtmgp_sw_post_run (nCol, nLev, nDay, idxday, lsswr, do_sw_clrsky ! 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) - + ! visdf - downward uv+vis diffused flux (W/m2) + real(kind=kind_phys), dimension(:,:), intent(inout) :: fluxr - + ! Outputs (mandatory) real(kind_phys), dimension(nCol), intent(inout) :: & nirbmdi, & ! sfc nir beam sw downward flux (W/m2) @@ -96,7 +96,7 @@ subroutine GFS_rrtmgp_sw_post_run (nCol, nLev, nDay, idxday, lsswr, do_sw_clrsky nirbmui, & ! sfc nir beam sw upward flux (W/m2) nirdfui, & ! sfc nir diff sw upward flux (W/m2) visbmui, & ! sfc uv+vis beam sw upward flux (W/m2) - visdfui, & ! sfc uv+vis diff sw upward flux (W/m2) + visdfui, & ! sfc uv+vis diff sw upward flux (W/m2) sfcnsw, & ! total sky sfc netsw flx into ground sfcdsw ! real(kind_phys), dimension(nCol,nLev), intent(inout) :: & @@ -119,7 +119,7 @@ subroutine GFS_rrtmgp_sw_post_run (nCol, nLev, nDay, idxday, lsswr, do_sw_clrsky ! dnfx0 - clear sky dnward flux (W/m2) real(kind_phys),dimension(nCol, nLev),intent(inout),optional :: & htrswc ! Clear-sky heating rate (K/s) - + ! Local variables integer :: i, j, k, iSFC, iTOA, itop, ibtc real(kind_phys) :: tem0d, tem1, tem2 @@ -135,7 +135,7 @@ subroutine GFS_rrtmgp_sw_post_run (nCol, nLev, nDay, idxday, lsswr, do_sw_clrsky ! Are any optional outputs requested? l_fluxessw2d = present(flxprf_sw) - + ! Are the components of the surface fluxes provided? l_scmpsw = present(scmpsw) @@ -150,7 +150,7 @@ subroutine GFS_rrtmgp_sw_post_run (nCol, nLev, nDay, idxday, lsswr, do_sw_clrsky iSFC = 1 iTOA = nLev+1 endif - + ! ####################################################################################### ! Compute SW heating-rates ! ####################################################################################### @@ -194,10 +194,10 @@ subroutine GFS_rrtmgp_sw_post_run (nCol, nLev, nDay, idxday, lsswr, do_sw_clrsky flxprf_sw(:,:)%upfx0 = fluxswUP_clrsky(:,:) flxprf_sw(:,:)%dnfx0 = fluxswDOWN_clrsky(:,:) endif - + ! Surface down and up spectral component fluxes ! - Save two spectral bands' surface downward and upward fluxes for output. - if (l_scmpsw) then + if (l_scmpsw) then do i=1,nCol nirbmdi(i) = scmpsw(i)%nirbm nirdfdi(i) = scmpsw(i)%nirdf @@ -209,15 +209,17 @@ subroutine GFS_rrtmgp_sw_post_run (nCol, nLev, nDay, idxday, lsswr, do_sw_clrsky visdfui(i) = scmpsw(i)%visdf * sfc_alb_uvvis_dif(1,i) enddo else - nirbmdi(:) = 0.0 - nirdfdi(:) = 0.0 - visbmdi(:) = 0.0 - visdfdi(:) = 0.0 - nirbmui(:) = 0.0 - nirdfui(:) = 0.0 - visbmui(:) = 0.0 - visdfui(:) = 0.0 - endif + do i=1,nCol + nirbmdi(i) = 0.0 + nirdfdi(i) = 0.0 + visbmdi(i) = 0.0 + visdfdi(i) = 0.0 + nirbmui(i) = 0.0 + nirdfui(i) = 0.0 + visbmui(i) = 0.0 + visdfui(i) = 0.0 + enddo + endif else ! if_nday_block ! ####################################################################################### ! Dark everywhere @@ -225,15 +227,17 @@ subroutine GFS_rrtmgp_sw_post_run (nCol, nLev, nDay, idxday, lsswr, do_sw_clrsky htrsw(:,:) = 0.0 sfcfsw = sfcfsw_type( 0.0, 0.0, 0.0, 0.0 ) topfsw = topfsw_type( 0.0, 0.0, 0.0 ) - nirbmdi(:) = 0.0 - nirdfdi(:) = 0.0 - visbmdi(:) = 0.0 - visdfdi(:) = 0.0 - nirbmui(:) = 0.0 - nirdfui(:) = 0.0 - visbmui(:) = 0.0 - visdfui(:) = 0.0 - + do i=1,nCol + nirbmdi(i) = 0.0 + nirdfdi(i) = 0.0 + visbmdi(i) = 0.0 + visdfdi(i) = 0.0 + nirbmui(i) = 0.0 + nirdfui(i) = 0.0 + visbmui(i) = 0.0 + visdfui(i) = 0.0 + enddo + if (do_sw_clrsky_hr) then htrswc(:,:) = 0 endif @@ -279,7 +283,7 @@ subroutine GFS_rrtmgp_sw_post_run (nCol, nLev, nDay, idxday, lsswr, do_sw_clrsky fluxr(i,27) = fluxr(i,27) + nirdfdi(i) * tem0d ! nir diff sw dn ! SW clear-sky fluxes fluxr(i,29) = fluxr(i,29) + topfsw(i)%upfx0 * tem0d - fluxr(i,31) = fluxr(i,31) + sfcfsw(i)%upfx0 * tem0d + fluxr(i,31) = fluxr(i,31) + sfcfsw(i)%upfx0 * tem0d fluxr(i,32) = fluxr(i,32) + sfcfsw(i)%dnfx0 * tem0d endif enddo diff --git a/physics/rte-rrtmgp b/physics/rte-rrtmgp index 9588c7bd8..d9594c46c 160000 --- a/physics/rte-rrtmgp +++ b/physics/rte-rrtmgp @@ -1 +1 @@ -Subproject commit 9588c7bd89e4f51a924f766e313bc42830fb4479 +Subproject commit d9594c46c877a2ab8001f5cd37961efdcf08ad8e From 215c80e4a8caa4032eb50f4df99a0414a8bbce09 Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Fri, 15 Oct 2021 14:05:54 +0000 Subject: [PATCH 02/15] additional updates to emissivity calculation etc --- physics/GFS_debug.F90 | 2 +- physics/GFS_radiation_surface.F90 | 24 +++---- physics/GFS_radiation_surface.meta | 24 ++++++- physics/GFS_surface_composites.F90 | 27 +++---- physics/GFS_surface_composites.meta | 36 ---------- physics/radiation_surface.f | 107 +++++++++++----------------- 6 files changed, 89 insertions(+), 131 deletions(-) diff --git a/physics/GFS_debug.F90 b/physics/GFS_debug.F90 index deb88458b..33e5beec0 100644 --- a/physics/GFS_debug.F90 +++ b/physics/GFS_debug.F90 @@ -1307,7 +1307,7 @@ subroutine GFS_interstitialtoscreen_run (Model, Statein, Stateout, Sfcprop, Coup call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Interstitial%scmpsw%visdf ', Interstitial%scmpsw%visdf ) call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Interstitial%semis_ice ', Interstitial%semis_ice ) call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Interstitial%semis_land ', Interstitial%semis_land ) - call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Interstitial%semis_water ', Interstitial%semis_water ) +! call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Interstitial%semis_water ', Interstitial%semis_water ) call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Interstitial%sfcalb ', Interstitial%sfcalb ) call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Interstitial%sigma ', Interstitial%sigma ) call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Interstitial%sigmaf ', Interstitial%sigmaf ) diff --git a/physics/GFS_radiation_surface.F90 b/physics/GFS_radiation_surface.F90 index 02d0f1c57..d48bce332 100644 --- a/physics/GFS_radiation_surface.F90 +++ b/physics/GFS_radiation_surface.F90 @@ -58,11 +58,11 @@ end subroutine GFS_radiation_surface_init subroutine GFS_radiation_surface_run ( & im, frac_grid, lslwr, lsswr, lsm, lsm_noahmp, lsm_ruc, & xlat, xlon, slmsk, lndp_type, n_var_lndp, sfc_alb_pert, & - lndp_var_list, lndp_prt_list, landfrac, snowd, sncovr, & + lndp_var_list, lndp_prt_list, landfrac, snodl, snodi, sncovr, & sncovr_ice, fice, zorl, hprime, tsfg, tsfa, tisfc, coszen, & cplice, min_seaice, min_lakeice, lakefrac, use_flake, & alvsf, alnsf, alvwf, alnwf, facsf, facwf, & - semis_lnd, semis_ice, snoalb, use_cice_alb, & + semis_lnd, semis_ice, semis_wat, snoalb, use_cice_alb, & albdvis_lnd, albdnir_lnd, albivis_lnd, albinir_lnd, & albdvis_ice, albdnir_ice, albivis_ice, albinir_ice, & semisbase, semis, sfcalb, sfc_alb_dif, errmsg, errflg) @@ -82,7 +82,7 @@ subroutine GFS_radiation_surface_run ( & real(kind=kind_phys), dimension(:), intent(in) :: xlat, xlon, slmsk, & sfc_alb_pert, lndp_prt_list, & landfrac, lakefrac, & - snowd, sncovr, & + snodl, snodi, sncovr, & sncovr_ice, fice, zorl, & hprime, tsfg, tsfa, tisfc, & coszen, alvsf, alnsf, alvwf, & @@ -93,7 +93,7 @@ subroutine GFS_radiation_surface_run ( & real(kind=kind_phys), dimension(:), intent(inout) :: albdvis_lnd, albdnir_lnd, & albivis_lnd, albinir_lnd, & - semis_lnd, semis_ice + semis_lnd, semis_ice, semis_wat real(kind=kind_phys), dimension(:), intent(inout) :: semisbase, semis real(kind=kind_phys), dimension(:,:), intent(inout) :: sfcalb real(kind=kind_phys), dimension(:), intent(inout) :: sfc_alb_dif @@ -161,13 +161,13 @@ subroutine GFS_radiation_surface_run ( & if (lslwr) then !> - Call module_radiation_surface::setemis(),to set up surface !! emissivity for LW radiation. - call setemis (lsm, lsm_noahmp, lsm_ruc, frac_grid, cplice, & - use_flake, lakefrac, xlon, xlat, slmsk, & -! frac_grid, min_seaice, xlon, xlat, slmsk, & - snowd, sncovr, sncovr_ice, zorl, tsfg, tsfa, & - hprime, semis_lnd, semis_ice, im, & - fracl, fraco, fraci, icy, & ! --- inputs - semisbase, semis) ! --- outputs + call setemis (lsm, lsm_noahmp, lsm_ruc, frac_grid, cplice, & + use_flake, lakefrac, xlon, xlat, slmsk, & +! frac_grid, min_seaice, xlon, xlat, slmsk, & + snodl, snodi, sncovr, sncovr_ice, zorl, tsfg, & + tsfa, hprime, semis_lnd, semis_ice, semis_wat,& + im, fracl, fraco, fraci, icy, & ! --- inputs + semisbase, semis) ! --- outputs endif if (lsswr) then @@ -184,7 +184,7 @@ subroutine GFS_radiation_surface_run ( & !> - Call module_radiation_surface::setalb(),to set up surface !! albedor for SW radiation. - call setalb (slmsk, lsm, lsm_noahmp, lsm_ruc, use_cice_alb, snowd, sncovr, sncovr_ice, & + call setalb (slmsk, lsm, lsm_noahmp, lsm_ruc, use_cice_alb, snodl, sncovr, sncovr_ice, & snoalb, zorl, coszen, tsfg, tsfa, hprime, frac_grid, lakefrac, & ! snoalb, zorl, coszen, tsfg, tsfa, hprime, frac_grid, min_seaice, & alvsf, alnsf, alvwf, alnwf, facsf, facwf, fice, tisfc, & diff --git a/physics/GFS_radiation_surface.meta b/physics/GFS_radiation_surface.meta index d360b37d8..5aa40ff1f 100644 --- a/physics/GFS_radiation_surface.meta +++ b/physics/GFS_radiation_surface.meta @@ -197,9 +197,18 @@ kind = kind_phys intent = in optional = F -[snowd] - standard_name = lwe_surface_snow - long_name = water equivalent snow depth +[snodl] + standard_name = surface_snow_thickness_water_equivalent_over_land + long_name = water equivalent snow depth over land + units = mm + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in + optional = F +[snodi] + standard_name = surface_snow_thickness_water_equivalent_over_ice + long_name = water equivalent snow depth over ice units = mm dimensions = (horizontal_loop_extent) type = real @@ -402,6 +411,15 @@ kind = kind_phys intent = in optional = F +[semis_wat] + standard_name = surface_longwave_emissivity_over_water + long_name = surface lw emissivity in fraction over water + units = frac + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = inout + optional = F [snoalb] standard_name = upper_bound_of_max_albedo_assuming_deep_snow long_name = maximum snow albedo diff --git a/physics/GFS_surface_composites.F90 b/physics/GFS_surface_composites.F90 index ca5ea2765..2ad6ef3d8 100644 --- a/physics/GFS_surface_composites.F90 +++ b/physics/GFS_surface_composites.F90 @@ -34,8 +34,7 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra tprcp_lnd, tprcp_ice, uustar, uustar_wat, uustar_lnd, uustar_ice, & weasd, weasd_lnd, weasd_ice, ep1d_ice, tsfc, tsfco, tsfcl, tsfc_wat, & tisfc, tsurf_wat, tsurf_lnd, tsurf_ice, & - gflx_ice, tgice, islmsk, islmsk_cice, slmsk, semis_wat, semis_lnd, semis_ice, & - qss, qss_wat, qss_lnd, qss_ice, & + gflx_ice, tgice, islmsk, islmsk_cice, slmsk, qss, qss_wat, qss_lnd, qss_ice, & min_lakeice, min_seaice, kdt, huge, errmsg, errflg) implicit none @@ -57,7 +56,7 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra qss_wat, qss_lnd, qss_ice, ep1d_ice, gflx_ice real(kind=kind_phys), intent(in ) :: tgice integer, dimension(:), intent(inout) :: islmsk, islmsk_cice - real(kind=kind_phys), dimension(:), intent(inout) :: semis_wat, semis_lnd, semis_ice, slmsk + real(kind=kind_phys), dimension(:), intent(inout) :: slmsk real(kind=kind_phys), intent(in ) :: min_lakeice, min_seaice, huge ! real(kind=kind_phys), dimension(:), intent(inout) :: zorlo, zorll, zorli @@ -212,11 +211,6 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra uustar_wat(i) = uustar(i) tsfc_wat(i) = tsfco(i) tsurf_wat(i) = tsfco(i) - !-- reference emiss value for surface emissivity in setemis - ! 1-open water, 2-grass/shrub land, 3-bare soil, tundra, - ! 4-sandy desert, 5-rocky desert, 6-forest, 7-ice, 8-snow - !data emsref / 0.97, 0.95, 0.94, 0.90, 0.93, 0.96, 0.96, 0.99 / - semis_wat(i) = 0.97_kind_phys ! consistent with setemis ! DH* else zorlo(i) = huge @@ -325,8 +319,9 @@ end subroutine GFS_surface_composites_inter_finalize !> \section arg_table_GFS_surface_composites_inter_run Argument Table !! \htmlinclude GFS_surface_composites_inter_run.html !! - subroutine GFS_surface_composites_inter_run (im, dry, icy, wet, semis_wat, semis_lnd, semis_ice, adjsfcdlw, & - gabsbdlw_lnd, gabsbdlw_ice, gabsbdlw_wat, & +! subroutine GFS_surface_composites_inter_run (im, dry, icy, wet, semis_wat, semis_lnd, semis_ice, adjsfcdlw, & + subroutine GFS_surface_composites_inter_run (im, dry, icy, wet, semis_lnd, semis_ice, adjsfcdlw, & + gabsbdlw_lnd, gabsbdlw_ice, gabsbdlw_wat, & adjsfcusw, adjsfcdsw, adjsfcnsw, errmsg, errflg) implicit none @@ -334,7 +329,8 @@ subroutine GFS_surface_composites_inter_run (im, dry, icy, wet, semis_wat, semis ! Interface variables integer, intent(in ) :: im logical, dimension(:), intent(in ) :: dry, icy, wet - real(kind=kind_phys), dimension(:), intent(in ) :: semis_wat, semis_lnd, semis_ice, adjsfcdlw, & +! real(kind=kind_phys), dimension(:), intent(in ) :: semis_wat, semis_lnd, semis_ice, adjsfcdlw, & + real(kind=kind_phys), dimension(:), intent(in ) :: semis_lnd, semis_ice, adjsfcdlw, & adjsfcdsw, adjsfcnsw real(kind=kind_phys), dimension(:), intent(inout) :: gabsbdlw_lnd, gabsbdlw_ice, gabsbdlw_wat real(kind=kind_phys), dimension(:), intent(out) :: adjsfcusw @@ -343,6 +339,13 @@ subroutine GFS_surface_composites_inter_run (im, dry, icy, wet, semis_wat, semis character(len=*), intent(out) :: errmsg integer, intent(out) :: errflg +! + !-- reference emiss value for surface emissivity in setemis + ! 1-open water, 2-grass/shrub land, 3-bare soil, tundra, + ! 4-sandy desert, 5-rocky desert, 6-forest, 7-ice, 8-snow + !data emsref / 0.97, 0.95, 0.94, 0.90, 0.93, 0.96, 0.96, 0.99 / + real(kind=kind_phys), parameter :: semis_wat = 0.97_kind_phys ! consistent with setemis + ! Local variables integer :: i @@ -371,7 +374,7 @@ subroutine GFS_surface_composites_inter_run (im, dry, icy, wet, semis_wat, semis do i=1,im if (dry(i)) gabsbdlw_lnd(i) = semis_lnd(i) * adjsfcdlw(i) if (icy(i)) gabsbdlw_ice(i) = semis_ice(i) * adjsfcdlw(i) - if (wet(i)) gabsbdlw_wat(i) = semis_wat(i) * adjsfcdlw(i) + if (wet(i)) gabsbdlw_wat(i) = semis_wat * adjsfcdlw(i) adjsfcusw(i) = adjsfcdsw(i) - adjsfcnsw(i) enddo diff --git a/physics/GFS_surface_composites.meta b/physics/GFS_surface_composites.meta index a8f76e2ed..7d60d2b82 100644 --- a/physics/GFS_surface_composites.meta +++ b/physics/GFS_surface_composites.meta @@ -459,33 +459,6 @@ kind = kind_phys intent = inout optional = F -[semis_wat] - standard_name = surface_longwave_emissivity_over_water - long_name = surface lw emissivity in fraction over water - units = frac - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = inout - optional = F -[semis_lnd] - standard_name = surface_longwave_emissivity_over_land - long_name = surface lw emissivity in fraction over land - units = frac - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = inout - optional = F -[semis_ice] - standard_name = surface_longwave_emissivity_over_ice - long_name = surface lw emissivity in fraction over ice - units = frac - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = inout - optional = F [qss] standard_name = surface_specific_humidity long_name = surface air saturation specific humidity @@ -617,15 +590,6 @@ type = logical intent = in optional = F -[semis_wat] - standard_name = surface_longwave_emissivity_over_water - long_name = surface lw emissivity in fraction over water - units = frac - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = in - optional = F [semis_lnd] standard_name = surface_longwave_emissivity_over_land long_name = surface lw emissivity in fraction over land diff --git a/physics/radiation_surface.f b/physics/radiation_surface.f index 3bb53be73..29cae3992 100644 --- a/physics/radiation_surface.f +++ b/physics/radiation_surface.f @@ -298,7 +298,7 @@ end subroutine sfc_init !! \n 1) climatological surface albedo scheme (\cite briegleb_1992) !! \n 2) MODIS retrieval based scheme from Boston univ. !!\param slmsk (IMAX), sea(0),land(1),ice(2) mask on fcst model grid -!!\param snowf (IMAX), snow depth water equivalent in mm +!!\param snowf (IMAX), snow depth water equivalent in mm over land !!\param sncovr (IMAX), snow cover over land !!\param snoalb (IMAX), maximum snow albedo over land (for deep snow) !!\param zorlf (IMAX), surface roughness in cm @@ -712,7 +712,8 @@ end subroutine setalb !! or -pi -> +pi ranges !!\param xlat (IMAX), latitude in radiance, default to pi/2 -> !! -pi/2 range, otherwise see in-line comment -!!\param snowf (IMAX), snow depth water equivalent in mm +!!\param snodl (IMAX), snow depth water equivalent in mm land +!!\param snodi (IMAX), snow depth water equivalent in mm ice !!\param sncovr (IMAX), snow cover over land !!\param zorlf (IMAX), surface roughness in cm !!\param tsknf (IMAX), ground surface temperature in K @@ -725,9 +726,9 @@ end subroutine setalb !----------------------------------- subroutine setemis & & ( lsm,lsm_noahmp,lsm_ruc,frac_grid,cplice,use_flake, & ! --- inputs: - & lakefrac,xlon,xlat,slmsk,snowf,sncovr,sncovr_ice, & + & lakefrac,xlon,xlat,slmsk,snodl,snodi,sncovr,sncovr_ice, & & zorlf,tsknf,tairf,hprif, & - & semis_lnd,semis_ice,IMAX,fracl,fraco,fraci,icy, & + & semis_lnd,semis_ice,semis_wat,IMAX,fracl,fraco,fraci,icy, & & semisbase, sfcemis & ! --- outputs: & ) @@ -748,7 +749,8 @@ subroutine setemis & ! xlat (IMAX) - latitude in radiance, default to pi/2 -> -pi/2 ! ! range, otherwise see in-line comment ! ! slmsk (IMAX) - sea(0),land(1),ice(2) mask on fcst model grid ! -! snowf (IMAX) - snow depth water equivalent in mm ! +! snodl (IMAX) - snow depth water equivalent in mm over land ! +! snodi (IMAX) - snow depth water equivalent in mm over ice ! ! sncovr(IMAX) - ialbflg=1: snow cover over land in fraction ! ! sncovr_ice(IMAX) - snow cover over ice in fraction ! ! zorlf (IMAX) - surface roughness in cm ! @@ -757,6 +759,7 @@ subroutine setemis & ! hprif (IMAX) - topographic sdv in m ! ! semis_lnd (IMAX) - land emissivity ! ! semis_ice (IMAX) - ice emissivity ! +! semis_wat (IMAX) - water emissivity ! ! IMAX - array horizontal dimension ! ! ! ! outputs: ! @@ -787,12 +790,12 @@ subroutine setemis & real (kind=kind_phys), dimension(:), intent(in) :: lakefrac real (kind=kind_phys), dimension(:), intent(in) :: & - & xlon,xlat, slmsk, snowf,sncovr, sncovr_ice, & + & xlon,xlat, slmsk, snodl, snodi, sncovr, sncovr_ice, & & zorlf, tsknf, tairf, hprif real (kind=kind_phys), dimension(:), intent(in) :: & & fracl, fraco, fraci real (kind=kind_phys), dimension(:), intent(inout) :: & - & semis_lnd, semis_ice + & semis_lnd, semis_ice, semis_wat logical, dimension(:), intent(in) :: & & icy @@ -805,7 +808,7 @@ subroutine setemis & integer :: ivgtyp real (kind=kind_phys) :: dltg, hdlt, tmp1, tmp2, & - & asnow, argh, hrgh, fsno, fsnol, fsnoi, snowc + & asnow, argh, hrgh, fsno real (kind=kind_phys) :: sfcemis_land, sfcemis_ice ! --- reference emiss value for diff surface emiss index @@ -819,6 +822,8 @@ subroutine setemis & !===> ... begin here ! !> -# Set emissivity by surface type and conditions + + semis_wat = emsref(1) if ( iemslw == 1 ) then dltg = 360.0 / float(IMXEMS) @@ -830,20 +835,16 @@ subroutine setemis & lab_do_IMAX : do i = 1, IMAX - snowc = sncovr(i) * fracl(i) if (.not. cplice .or. lakefrac(i) > f_zero) then semis_ice(i) = emsref(7) - snowc = sncovr(i) + sncovr_ice(i)*fraci(i) endif if (fracl(i) < epsln) then ! no land if ( abs(fraco(i)-f_one) < epsln ) then ! open water point sfcemis(i) = emsref(1) elseif ( abs(fraci(i)-f_one) < epsln ) then ! complete sea/lake ice -! sfcemis(i) = emsref(7) sfcemis(i) = semis_ice(i) else !-- fractional sea ice -! sfcemis(i) = fraco(i)*emsref(1) + fraci(i)*emsref(7) sfcemis(i) = fraco(i)*emsref(1) + fraci(i)*semis_ice(i) endif @@ -889,65 +890,39 @@ subroutine setemis & semisbase(i) = sfcemis(i) semis_lnd(i) = emsref(idx) - endif ! end if_slmsk_block + endif !> - Check for snow covered area. +!> it is assume here that "sncovr" is the fraction of land covered by snow +!> and "sncovr_ice" is the fraction of ice coverd by snow - if (snowc > f_zero) then ! input land/ice area snow cover - -! it is assume here that "sncovr" is the fraction of land covered by snow -! and "sncovr_ice" is the fraction of ice coverd by snow - + if (fracl(i) > epsln) then if (sncovr(i) > f_zero) then semis_lnd(i) = semis_lnd(i) * (f_one - sncovr(i)) & & + emsref(8) * sncovr(i) + elseif (snodl(i) > f_zero) then + asnow = 0.02*snodl(i) + argh = min(0.50, max(.025, 0.01*zorlf(i))) + hrgh = min(f_one, max(0.20, 1.0577-1.1538e-3*hprif(i) ) ) + fsno = min(f_one, max(f_zero, asnow/(argh+asnow) * hrgh)) + semis_lnd(i) = semis_lnd(i)*(f_one-fsno) + emsref(8)*fsno endif - if (sncovr_ice(i) > f_zero .and. & - & (lakefrac(i) > f_zero .or. .not. cplice)) then + endif + if (fraci(i) > epsln .and. & + & (lakefrac(i) > f_zero .or. .not. cplice)) then + if (sncovr_ice(i) > f_zero) then semis_ice(i) = semis_ice(i) * (f_one - sncovr_ice(i)) & & + emsref(8) * sncovr_ice(i) - endif - sfcemis(i) = fracl(i)*semis_lnd(i) + fraco(i)*emsref(1) & - & + fraci(i)*semis_ice(i) - - else ! compute snow cover from snow depth - if (abs(fraco(i)-f_one) > epsln .and. & - & snowf(i) > f_zero) then - asnow = 0.02*snowf(i) + elseif (snodi(i) > f_zero) then + asnow = 0.02*snodi(i) argh = min(0.50, max(.025, 0.01*zorlf(i))) hrgh = min(f_one, max(0.20, 1.0577-1.1538e-3*hprif(i) ) ) - tmp1 = fracl(i) + fraci(i) - if (tmp1 > f_zero) then - fsno = min(tmp1, asnow / (argh + asnow) * hrgh) - tmp2 = fsno / tmp1 - fsnol = fracl(i) * tmp2 - fsnoi = fraci(i) * tmp2 - - if (fracl(i) > f_zero) then - if (fracl(i) <= fsnol) then - semis_lnd(i) = emsref(8) - else - tmp1 = (fracl(i)-fsnol) / fracl(i) - semis_lnd(i) = semis_lnd(i) * tmp1 & - & + emsref(8) * (f_one-tmp1) - endif - endif - if (fraci(i) > f_zero .and. & - & (lakefrac(i) > f_zero .or. .not. cplice)) then - if (fraci(i) <= fsnoi) then - semis_ice(i) = emsref(8) - else - tmp1 = (fraci(i)-fsnoi) / fraci(i) - semis_ice(i) = semis_ice(i) * tmp1 & - & + emsref(8) * (f_one-tmp1) - endif - endif - endif + fsno = min(f_one, max(f_zero, asnow/(argh+asnow) * hrgh)) + semis_ice(i) = semis_ice(i)*(f_one-fsno) + emsref(8)*fsno endif - sfcemis(i) = fracl(i)*semis_lnd(i) + fraco(i)*emsref(1) & - & + fraci(i)*semis_ice(i) - - endif ! end if_ialbflg + endif + sfcemis(i) = fracl(i)*semis_lnd(i) + fraco(i)*emsref(1) & + & + fraci(i)*semis_ice(i) enddo lab_do_IMAX @@ -964,13 +939,12 @@ subroutine setemis & if (sncovr_ice(i) > f_zero) then sfcemis_ice = emsref(7) * (f_one-sncovr_ice(i)) & & + emsref(8) * sncovr_ice(i) - elseif (snowf(i) > f_zero) then - asnow = 0.02*snowf(i) + elseif (snodi(i) > f_zero) then + asnow = 0.02*snodi(i) argh = min(0.50, max(.025,0.01*zorlf(i))) hrgh = min(f_one,max(0.20,1.0577-1.1538e-3*hprif(i))) fsno = asnow / (argh + asnow) * hrgh - fsnoi = min(f_one, fsno / (fraci(i)+fracl(i))) - sfcemis_ice = emsref(7)*(f_one-fsnoi)+emsref(8)*fsnoi + sfcemis_ice = emsref(7)*(f_one-fsno) + emsref(8)*fsno endif semis_ice(i) = sfcemis_ice else @@ -981,13 +955,12 @@ subroutine setemis & if (sncovr_ice(i) > f_zero) then sfcemis_ice = emsref(7) * (f_one-sncovr_ice(i)) & & + emsref(8) * sncovr_ice(i) - elseif (snowf(i) > f_zero) then - asnow = 0.02*snowf(i) + elseif (snodi(i) > f_zero) then + asnow = 0.02*snodi(i) argh = min(0.50, max(.025,0.01*zorlf(i))) hrgh = min(f_one,max(0.20,1.0577-1.1538e-3*hprif(i))) fsno = asnow / (argh + asnow) * hrgh - fsnoi = min(f_one, fsno / (fraci(i)+fracl(i))) - sfcemis_ice = emsref(7)*(f_one-fsnoi)+emsref(8)*fsnoi + sfcemis_ice = emsref(7)*(f_one-fsno) + emsref(8)*fsno endif semis_ice(i) = sfcemis_ice else From 3d9b79ce9a55c27dc3d88733828e142eb6d3c526 Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Wed, 20 Oct 2021 08:01:08 -0400 Subject: [PATCH 03/15] updating .gitmodules --- .gitmodules | 6 +++--- physics/rte-rrtmgp | 2 +- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/.gitmodules b/.gitmodules index 5bcc65869..ec9a09abb 100644 --- a/.gitmodules +++ b/.gitmodules @@ -1,4 +1,4 @@ [submodule "physics/rte-rrtmgp"] - path = physics/rte-rrtmgp - url = https://github.com/earth-system-radiation/rte-rrtmgp - branch = dtc/ccpp + path = physics/rte-rrtmgp + url = https://github.com/earth-system-radiation/rte-rrtmgp + branch = dtc/ccpp diff --git a/physics/rte-rrtmgp b/physics/rte-rrtmgp index d9594c46c..9588c7bd8 160000 --- a/physics/rte-rrtmgp +++ b/physics/rte-rrtmgp @@ -1 +1 @@ -Subproject commit d9594c46c877a2ab8001f5cd37961efdcf08ad8e +Subproject commit 9588c7bd89e4f51a924f766e313bc42830fb4479 From c4c8bad1e157672afeea5b5cb12b07ac2b06c4c5 Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Wed, 20 Oct 2021 19:25:20 +0000 Subject: [PATCH 04/15] changes consistent with branch SM_Sept21_PR --- physics/radiation_surface.f | 30 +++++++++++++++++------------- 1 file changed, 17 insertions(+), 13 deletions(-) diff --git a/physics/radiation_surface.f b/physics/radiation_surface.f index 05e62ae88..066bcfbef 100644 --- a/physics/radiation_surface.f +++ b/physics/radiation_surface.f @@ -20,18 +20,22 @@ ! ! ! 'setalb' -- set up four-component surface albedoes ! ! inputs: ! -! (slmsk,snowf,sncovr,snoalb,zorlf,coszf,tsknf,tairf,hprif, ! +! (slmsk,snodi,sncovr,snoalb,zorlf,coszf,tsknf,tairf,hprif, ! ! alvsf,alnsf,alvwf,alnwf,facsf,facwf,fice,tisfc ! ! IMAX) ! ! outputs: ! ! (sfcalb) ! ! ! ! 'setemis' -- set up surface emissivity for lw radiation ! -! inputs: ! -! (xlon,xlat,slmsk,snowf,sncovr,zorlf,tsknf,tairf,hprif, ! -! IMAX) ! -! outputs: ! -! (sfcemis) ! +! ( lsm,lsm_noahmp,lsm_ruc,frac_grid,cplice,use_flake, ! +! --- inputs: +! lakefrac,xlon,xlat,slmsk,snodl,snodi,sncovr,sncovr_ice, ! +! zorlf,tsknf,tairf,hprif, ! +! semis_lnd,semis_ice,semis_wat,IMAX,fracl,fraco,fraci,icy, ! +! +! --- outputs: +! semisbase, sfcemis ! +! ! ! ! external modules referenced: ! ! ! @@ -298,7 +302,7 @@ end subroutine sfc_init !! \n 1) climatological surface albedo scheme (\cite briegleb_1992) !! \n 2) MODIS retrieval based scheme from Boston univ. !!\param slmsk (IMAX), sea(0),land(1),ice(2) mask on fcst model grid -!!\param snowf (IMAX), snow depth water equivalent in mm over land +!!\param snodi (IMAX), snow depth water equivalent in mm over ice !!\param sncovr (IMAX), snow cover over land !!\param snoalb (IMAX), maximum snow albedo over land (for deep snow) !!\param zorlf (IMAX), surface roughness in cm @@ -332,7 +336,7 @@ end subroutine sfc_init !! @{ !----------------------------------- subroutine setalb & - & ( slmsk,lsm,lsm_noahmp,lsm_ruc,use_cice_alb,snowf, & ! --- inputs: + & ( slmsk,lsm,lsm_noahmp,lsm_ruc,use_cice_alb,snodi, & ! --- inputs: & sncovr,sncovr_ice,snoalb,zorlf,coszf, & & tsknf,tairf,hprif,frac_grid, lakefrac, & & alvsf,alnsf,alvwf,alnwf,facsf,facwf,fice,tisfc, & @@ -358,7 +362,7 @@ subroutine setalb & ! ! ! inputs: ! ! slmsk (IMAX) - sea(0),land(1),ice(2) mask on fcst model grid ! -! snowf (IMAX) - snow depth water equivalent in mm ! +! snodi (IMAX) - snow depth water equivalent in mm over ice ! ! sncovr(IMAX) - ialgflg=0: not used ! ! ialgflg=1: snow cover over land in fraction ! ! sncovr_ice(IMAX) - ialgflg=0: not used ! @@ -410,7 +414,7 @@ subroutine setalb & real (kind=kind_phys), dimension(:), intent(in) :: & & lakefrac, & - & slmsk, snowf, zorlf, coszf, tsknf, tairf, hprif, & + & slmsk, snodi, zorlf, coszf, tsknf, tairf, hprif, & & alvsf, alnsf, alvwf, alnwf, facsf, facwf, fice, tisfc, & & icealbdvis, icealbdnir, icealbivis, icealbinir, & & sncovr, sncovr_ice, snoalb, albPpert ! sfc-perts, mgehne @@ -478,7 +482,7 @@ subroutine setalb & asevb_ice = icealbdvis(i) asenb_ice = icealbdnir(i) else - asnow = 0.02*snowf(i) + asnow = 0.02*snodi(i) argh = min(0.50, max(.025, 0.01*zorlf(i))) hrgh = min(f_one,max(0.20,1.0577-1.1538e-3*hprif(i))) fsno0 = asnow / (argh + asnow) * hrgh ! snow fraction on ice @@ -614,7 +618,7 @@ subroutine setalb & asenb_ice = icealbdnir(i) else !-- Computation of ice albedo - asnow = 0.02*snowf(i) + asnow = 0.02*snodi(i) argh = min(0.50, max(.025, 0.01*zorlf(i))) hrgh = min(f_one,max(0.20,1.0577-1.1538e-3*hprif(i))) fsno0 = asnow / (argh + asnow) * hrgh @@ -750,7 +754,7 @@ subroutine setemis & ! range, otherwise see in-line comment ! ! slmsk (IMAX) - sea(0),land(1),ice(2) mask on fcst model grid ! ! snodl (IMAX) - snow depth water equivalent in mm over land ! -! snodi (IMAX) - snow depth water equivalent in mm over ice ! +! snodi (IMAX) - snow depth water equivalent in mm over ice ! ! sncovr(IMAX) - ialbflg=1: snow cover over land in fraction ! ! sncovr_ice(IMAX) - snow cover over ice in fraction ! ! zorlf (IMAX) - surface roughness in cm ! From 3233ebac0bd1ee938da11c6a4410c2117e2bbf9d Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Thu, 4 Nov 2021 16:23:05 +0000 Subject: [PATCH 05/15] after submodule sync --- physics/rte-rrtmgp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/physics/rte-rrtmgp b/physics/rte-rrtmgp index 9588c7bd8..9c51cb7c3 160000 --- a/physics/rte-rrtmgp +++ b/physics/rte-rrtmgp @@ -1 +1 @@ -Subproject commit 9588c7bd89e4f51a924f766e313bc42830fb4479 +Subproject commit 9c51cb7c3e227c9e84c2bff29ce4f438c7a54ae6 From 49bf5c4dd7fa0b4ef85f94611f02a417d1c3f8a9 Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Wed, 24 Nov 2021 00:44:34 +0000 Subject: [PATCH 06/15] some cosmetic change in mg3 --- physics/micro_mg3_0.F90 | 222 ++++++++++++++++++++-------------------- 1 file changed, 111 insertions(+), 111 deletions(-) diff --git a/physics/micro_mg3_0.F90 b/physics/micro_mg3_0.F90 index dde143c4d..bca005bc9 100644 --- a/physics/micro_mg3_0.F90 +++ b/physics/micro_mg3_0.F90 @@ -16,20 +16,20 @@ !! This version:https://svn-ccsm-models.cgd.ucar.edu/cam1/branch_tags/mg3_tags/mg3_33_cam5_4_153/ !! !! \version 2 history: Sep 2011: Development begun. -!!\n Feb 2013: Added of prognostic precipitation. -!!\n Aug 2015: Published and released version +!!\n Feb 2013: Added of prognostic precipitation. +!!\n Aug 2015: Published and released version !! !! Contributions from: Sean Santos, Peter Caldwell, Xiaohong Liu and Steve Ghan !! !! - Anning Cheng adopted mg2 for FV3GFS 9/29/2017 -!!\n add GMAO ice conversion and Liu et. al liquid water +!!\n add GMAO ice conversion and Liu et. al liquid water !!\n conversion in 10/12/2017 !! !! - Anning showed promising results for FV3GFS on 10/15/2017 !! - S. Moorthi - Oct/Nov 2017 - optimized the MG2 code !! - S. Moorthi - Nov 2017 - made the sedimentation quasi-implicit !! - S. Moorthi - Feb 2018 - updated to MG3 - modified graupel sedimentation -!! other modifications to eliminate blowup. +!! other modifications to eliminate blowup. !! - S. Moorthi - Mar 2018 - fixed a few bugs and added option to run as MG2 !! - S. Moorthi - Oct,29,2018 - change nlb from nlev/3 to levels with p/ps < 0.05 (nlball) !! @@ -1162,7 +1162,7 @@ subroutine micro_mg_tend ( & qsfm(i,k) = qsatfac(i,k) enddo enddo - end if + endif ! if (lprnt) write(0,*)' cldm=',cldm(1,nlev-20:nlev) ! if (lprnt) write(0,*)' liqcldf=',liqcldf(1,nlev-20:nlev) @@ -1234,8 +1234,8 @@ subroutine micro_mg_tend ( & ! esl(i,k) = qsfm(i,k) * esl(i,k) relhum(i,k) = max(zero, min(q(i,k)/max(qvl(i,k), qsmall), two)) - end do - end do + enddo + enddo !=============================================== @@ -1590,7 +1590,7 @@ subroutine micro_mg_tend ( & mnuccd = zero end where - end if + endif !============================================================================= @@ -1628,11 +1628,11 @@ subroutine micro_mg_tend ( & ns(i,k) = max(ns(i,k) - ninstsm(i,k), zero) qr(i,k) = max(qr(i,k) + minstsm(i,k), zero) nr(i,k) = max(nr(i,k) + ninstsm(i,k), zero) - end if - end if + endif + endif - end do - end do + enddo + enddo ! if (lprnt) write(0,*)' tlat1=',tlat(1,:)*deltat ! if (lprnt) write(0,*)' qg1=',qg(1,:) @@ -1670,11 +1670,11 @@ subroutine micro_mg_tend ( & ng(i,k) = max(ng(i,k) - ninstgm(i,k), zero) qr(i,k) = max(qr(i,k) + minstgm(i,k), zero) nr(i,k) = max(nr(i,k) + ninstgm(i,k), zero) - end if - end if + endif + endif - end do - end do + enddo + enddo endif ! if (lprnt) write(0,*)' tlat1g=',tlat(1,:)*deltat @@ -1719,10 +1719,10 @@ subroutine micro_mg_tend ( & end if !--ag - end if - end if - end do - end do + endif + endif + enddo + enddo ! if (lprnt) then ! write(0,*)' tlat2=',tlat(1,:)*deltat @@ -1750,11 +1750,11 @@ subroutine micro_mg_tend ( & ! specify droplet concentration if (nccons) then ncic(i,k) = ncnst * rhoinv(i,k) - end if + endif else qcic(i,k) = zero ncic(i,k) = zero - end if + endif ! if (qi(i,k) >= qsmall .and. icldm(i,k) > mincld) then if (qi(i,k) >= qsmall) then @@ -1766,14 +1766,14 @@ subroutine micro_mg_tend ( & ! switch for specification of cloud ice number if (nicons) then niic(i,k) = ninst * rhoinv(i,k) - end if + endif else qiic(i,k) = zero niic(i,k) = zero - end if + endif - end do - end do + enddo + enddo !======================================================================== @@ -1802,12 +1802,12 @@ subroutine micro_mg_tend ( & ! then leave precip_frac as cloud fraction at current level if (k /= 1) then !++ag -! where (qr(:,k-1) >= qsmall .or. qs(:,k-1) >= qsmall .or. qg(:,k-1) >= qsmall) +! where (qr(:,k-1) >= qsmall .or. qs(:,k-1) >= qsmall .or. qg(:,k-1) >= qsmall) !--ag where (qr(:,k-1) >= qsmall .or. qs(:,k-1) >= qsmall) precip_frac(:,k) = max(precip_frac(:,k-1), precip_frac(:,k)) end where - end if + endif endif @@ -1916,13 +1916,13 @@ subroutine micro_mg_tend ( & else call ice_autoconversion(t(:,k), qiic(:,k), lami(:,k), n0i(:,k), & dcs, ts_au_loc(:), prci(:,k), nprci(:,k), mgncol) - end if + endif !else ! Add in the particles that we have already converted to snow, and ! don't do any further autoconversion of ice. !prci(:,k) = tnd_qsnow(:,k) / cldm(:,k) !nprci(:,k) = tnd_nsnow(:,k) / cldm(:,k) - end if + endif ! note, currently we don't have this ! inside the do_cldice block, should be changed later @@ -2125,7 +2125,7 @@ subroutine micro_mg_tend ( & !mnudep(:,k) = zero !end where - end if + endif else do i=1,mgncol @@ -2136,7 +2136,7 @@ subroutine micro_mg_tend ( & mnudep(i,k) = zero nnudep(i,k) = zero enddo - end if + endif call snow_self_aggregation(t(:,k), rho(:,k), asn(:,k), rhosn, qsic(:,k), nsic(:,k), & nsagg(:,k), mgncol) @@ -2150,7 +2150,7 @@ subroutine micro_mg_tend ( & else nsacwi(:,k) = zero msacwi(:,k) = zero - end if + endif call accrete_rain_snow(t(:,k), rho(:,k), umr(:,k), ums(:,k), unr(:,k), uns(:,k), & qric(:,k), qsic(:,k), lamr(:,k), n0r(:,k), lams(:,k), n0s(:,k), & @@ -2175,7 +2175,7 @@ subroutine micro_mg_tend ( & else prai(:,k) = zero nprai(:,k) = zero - end if + endif !++ag Moved below graupel conditional, now two different versions ! if (.not. (do_hail .or. do_graupel)) then @@ -2223,9 +2223,9 @@ subroutine micro_mg_tend ( & ! all ql is removed (which is handled elsewhere) !in fact, nothing in this entire file makes nsubc nonzero. nsubc(i,k) = zero - end do + enddo - end if !do_cldice + endif !do_cldice !---PMC 12/3/12 !++ag Process rate calls for graupel here. @@ -2337,7 +2337,7 @@ subroutine micro_mg_tend ( & pre(:,k), prds(:,k), am_evp_st(:,k), mgncol) - end if ! end do_graupel/hail loop + endif ! end do_graupel/hail loop !--ag do i=1,mgncol @@ -2382,7 +2382,7 @@ subroutine micro_mg_tend ( & qcrat(i,k) = ratio else qcrat(i,k) = one - end if + endif ! if(lprnt) write(0,*)' bergs2=',bergs(1,k),' k=',k,' ratio=',ratio @@ -2391,9 +2391,9 @@ subroutine micro_mg_tend ( & !deposition for the remaining frac of the timestep. if (qc(i,k) >= qsmall) then vap_dep(i,k) = vap_dep(i,k) * (one-qcrat(i,k)) - end if + endif - end do + enddo do i=1,mgncol @@ -2419,10 +2419,10 @@ subroutine micro_mg_tend ( & dum1 = mnuccd(i,k) / (vap_dep(i,k)+mnuccd(i,k)) mnuccd(i,k) = dum*dum1 vap_dep(i,k) = dum - mnuccd(i,k) - end if - end if + endif + endif - end do + enddo do i=1,mgncol @@ -2492,9 +2492,9 @@ subroutine micro_mg_tend ( & pracs(i,k) = ratio * pracs(i,k) mnuccr(i,k) = ratio * mnuccr(i,k) mnuccri(i,k) = ratio * mnuccri(i,k) - end if + endif - end do + enddo do i=1,mgncol @@ -2507,9 +2507,9 @@ subroutine micro_mg_tend ( & nsubr(i,k) = dum*nr(i,k) * oneodt else nsubr(i,k) = zero - end if + endif - end do + enddo do i=1,mgncol @@ -2535,9 +2535,9 @@ subroutine micro_mg_tend ( & nnuccr(i,k) = ratio * nnuccr(i,k) nsubr(i,k) = ratio * nsubr(i,k) nnuccri(i,k) = ratio * nnuccri(i,k) - end if + endif - end do + enddo if (do_cldice) then @@ -2569,11 +2569,11 @@ subroutine micro_mg_tend ( & prci(i,k) = ratio * prci(i,k) prai(i,k) = ratio * prai(i,k) ice_sublim(i,k) = ratio * ice_sublim(i,k) - end if + endif - end do + enddo - end if + endif if (do_cldice) then @@ -2585,7 +2585,7 @@ subroutine micro_mg_tend ( & tmpfrz = nnuccc(i,k) else tmpfrz = zero - end if + endif !++ag dum1 = (nprci(i,k)+nprai(i,k)-nsubi(i,k))*icldm(i,k) ! dum2 = nnuccd(i,k)+(nnucct(i,k)+tmpfrz+nnudep(i,k)+nsacwi(i,k))*lcldm(i,k) & @@ -2601,11 +2601,11 @@ subroutine micro_mg_tend ( & nprci(i,k) = ratio * nprci(i,k) nprai(i,k) = ratio * nprai(i,k) nsubi(i,k) = ratio * nsubi(i,k) - end if + endif - end do + enddo - end if + endif do i=1,mgncol @@ -2648,9 +2648,9 @@ subroutine micro_mg_tend ( & ! ratio = (qs(i,k)*oneodt + dum2) / dum1 * omsm ! prds(i,k) = ratio * prds(i,k) -! end if +! endif - end do + enddo do i=1,mgncol @@ -2678,7 +2678,7 @@ subroutine micro_mg_tend ( & ratio = (ns(i,k)*oneodt + dum2) / dum1 * omsm nscng(i,k) = ratio * nscng(i,k) ngracs(i,k) = ratio * ngracs(i,k) - end if + endif else dum1 = precip_frac(i,k)* (-nsubs(i,k)-nsagg(i,k)) @@ -2692,7 +2692,7 @@ subroutine micro_mg_tend ( & nsubs(i,k) = ratio * nsubs(i,k) nsagg(i,k) = ratio * nsagg(i,k) - end do + enddo !++ag Graupel Conservation Checks !------------------------------------------------------------------- @@ -2714,13 +2714,13 @@ subroutine micro_mg_tend ( & prdg(i,k) = ratio * prdg(i,k) - end if + endif - end do + enddo ! conservation of graupel number: not needed, no sinks !------------------------------------------------------------------- - end if + endif !--ag @@ -2801,10 +2801,10 @@ subroutine micro_mg_tend ( & dum1 = one - dum1 - dum2 - dum3 !--ag ice_sublim(i,k) = dum*dum1*oneodt - end if - end if + endif + endif - end do + enddo ! Big "administration" loop enforces conservation, updates variables ! that accumulate over substeps, and sets output variables. @@ -2902,7 +2902,7 @@ subroutine micro_mg_tend ( & qstend(i,k) = qstend(i,k) + (prai(i,k)+prci(i,k))*icldm(i,k)+(psacws(i,k)+bergs(i,k))*lcldm(i,k) & + (prds(i,k)+pracs(i,k)+mnuccr(i,k))*precip_frac(i,k) - end if + endif !--ag @@ -2938,7 +2938,7 @@ subroutine micro_mg_tend ( & else prodsnow(i,k) = (prai(i,k)+prci(i,k))*icldm(i,k)+(psacws(i,k)+bergs(i,k))*lcldm(i,k)+ & (pracs(i,k)+mnuccr(i,k))*precip_frac(i,k) - end if + endif ! following are used to calculate 1st order conversion rate of cloud water ! to rain and snow (1/s), for later use in aerosol wet removal routine @@ -3013,7 +3013,7 @@ subroutine micro_mg_tend ( & + (nnucct(i,k)+tmpfrz+nnudep(i,k)+nsacwi(i,k)+nmultg(i,k))*lcldm(i,k) & + (nsubi(i,k)-nprci(i,k)-nprai(i,k))*icldm(i,k) & + (nnuccri(i,k)+nmultrg(i,k))*precip_frac(i,k) - end if + endif if(do_graupel.or.do_hail) then ! nstend(i,k) = nstend(i,k) + (nsubs(i,k)+nsagg(i,k)+nnuccr(i,k))*precip_frac(i,k) & @@ -3030,7 +3030,7 @@ subroutine micro_mg_tend ( & nstend(i,k) = nstend(i,k) + (nsubs(i,k)+nsagg(i,k)+nnuccr(i,k))*precip_frac(i,k) & + nprci(i,k)*icldm(i,k) - end if + endif ! nrtend(i,k) = nrtend(i,k) + nprc(i,k)*lcldm(i,k)+(nsubr(i,k)-npracs(i,k)-nnuccr(i,k) & ! - nnuccri(i,k)+nragg(i,k))*precip_frac(i,k) @@ -3047,9 +3047,9 @@ subroutine micro_mg_tend ( & if (do_cldice .and. nitend(i,k) > zero .and. ni(i,k)+nitend(i,k)*deltat > nimax(i,k)) then nitend(i,k) = max(zero, (nimax(i,k)-ni(i,k))*oneodt) - end if + endif - end do + enddo ! End of "administration" loop @@ -3298,7 +3298,7 @@ subroutine micro_mg_tend ( & else fi(i,k) = zero fni(i,k)= zero - end if + endif ! fallspeed for rain @@ -3318,7 +3318,7 @@ subroutine micro_mg_tend ( & else fr(i,k) = zero fnr(i,k) = zero - end if + endif ! fallspeed for snow @@ -3337,7 +3337,7 @@ subroutine micro_mg_tend ( & else fs(i,k) = zero fns(i,k) = zero - end if + endif if (do_graupel .or. do_hail) then !++ag @@ -3359,7 +3359,7 @@ subroutine micro_mg_tend ( & else fg(i,k) = zero fng(i,k) = zero - end if + endif endif ! redefine dummy variables - sedimentation is calculated over grid-scale @@ -3386,7 +3386,7 @@ subroutine micro_mg_tend ( & if (dumg(i,k) < qsmall) dumng(i,k) = zero enddo - end do !!! vertical loop + enddo !!! vertical loop do k=1,nlev do i=1,mgncol @@ -3488,8 +3488,8 @@ subroutine micro_mg_tend ( & prect(i) = prect(i) + falouti(nlev) * (tx3*0.001_r8) preci(i) = preci(i) + falouti(nlev) * (tx3*0.001_r8) - end do - end if + enddo + endif ! if (lprnt) write(0,*)' tlat4=',tlat(1,:)*deltat ! calculate number of split time steps to ensure courant stability criteria @@ -3559,7 +3559,7 @@ subroutine micro_mg_tend ( & tlat(i,k) = tlat(i,k) - dum2 * xxlv lflx(i,k+1) = lflx(i,k+1) + faloutc(k) * tx3 ! Liquid condensate flux here - end do + enddo prect(i) = prect(i) + faloutc(nlev) * (tx3*0.001_r8) @@ -3629,11 +3629,11 @@ subroutine micro_mg_tend ( & faloutnr(k) = fnr(i,k) * dumnr(i,k) rflx(i,k+1) = rflx(i,k+1) + faloutr(k) * tx3 ! Rain Flux - end do + enddo prect(i) = prect(i) + faloutr(nlev) * (tx3*0.001_r8) - end do + enddo ! if (lprnt) write(0,*)' prectaftrain=',prect(i),' preci=',preci(i) @@ -3698,7 +3698,7 @@ subroutine micro_mg_tend ( & prect(i) = prect(i) + falouts(nlev) * (tx3*0.001_r8) preci(i) = preci(i) + falouts(nlev) * (tx3*0.001_r8) - end do !! nstep loop + enddo !! nstep loop ! if (lprnt) write(0,*)' prectaftssno=',prect(i),' preci=',preci(i) ! if (lprnt) write(0,*)' qgtnd1=',qgtend(1,:) @@ -3761,7 +3761,7 @@ subroutine micro_mg_tend ( & faloutng(k) = fng(i,k) * dumng(i,k) gflx(i,k+1) = gflx(i,k+1) + faloutg(k) * tx3 ! Ice flux - end do + enddo ! units below are m/s ! sedimentation flux at surface is added to precip flux at surface @@ -3770,7 +3770,7 @@ subroutine micro_mg_tend ( & prect(i) = prect(i) + faloutg(nlev) * (tx3*0.001_r8) preci(i) = preci(i) + faloutg(nlev) * (tx3*0.001_r8) - end do !! nstep loop + enddo !! nstep loop endif ! if (lprnt) write(0,*)' qgtnds=',qgtend(1,:) !--ag @@ -3813,18 +3813,18 @@ subroutine micro_mg_tend ( & ! switch for specification of droplet and crystal number if (nccons) then dumnc(i,k) = ncnst*rhoinv(i,k)*lcldm(i,k) - end if + endif ! switch for specification of cloud ice number if (nicons) then dumni(i,k) = ninst*rhoinv(i,k)*icldm(i,k) - end if + endif !++ag ! switch for specification of graupel number if (ngcons) then dumng(i,k) = ngnst*rhoinv(i,k)*precip_frac(i,k) - end if + endif !--ag if (dumc(i,k) < qsmall) dumnc(i,k) = zero @@ -3868,8 +3868,8 @@ subroutine micro_mg_tend ( & dum1 = - xlf * tx2 * dums(i,k) tlat(i,k) = dum1 + tlat(i,k) meltsdttot(i,k) = dum1 + meltsdttot(i,k) - end if - end if + endif + endif enddo enddo @@ -3904,8 +3904,8 @@ subroutine micro_mg_tend ( & dum1 = - xlf*tx2*dumg(i,k) tlat(i,k) = dum1 + tlat(i,k) meltsdttot(i,k) = dum1 + meltsdttot(i,k) - end if - end if + endif + endif enddo enddo @@ -3957,8 +3957,8 @@ subroutine micro_mg_tend ( & frzrdttot(i,k) = dum1 + frzrdttot(i,k) tlat(i,k) = dum1 + tlat(i,k) - end if - end if + endif + endif enddo enddo @@ -3992,8 +3992,8 @@ subroutine micro_mg_tend ( & qitend(i,k) = ((one-dum)*dumi(i,k)-qi(i,k)) * oneodt nitend(i,k) = ((one-dum)*dumni(i,k)-ni(i,k)) * oneodt tlat(i,k) = tlat(i,k) - xlf*tx2*dumi(i,k) - end if - end if + endif + endif enddo enddo @@ -4029,8 +4029,8 @@ subroutine micro_mg_tend ( & qctend(i,k) = ((one-dum)*dumc(i,k)-qc(i,k)) * oneodt nctend(i,k) = ((one-dum)*dumnc(i,k)-nc(i,k)) * oneodt tlat(i,k) = tlat(i,k) + xlf*tx2 - end if - end if + endif + endif enddo enddo ! remove any excess over-saturation, which is possible due to non-linearity when adding @@ -4077,10 +4077,10 @@ subroutine micro_mg_tend ( & ! for output qvres(i,k) = -dum tlat(i,k) = tlat(i,k) + dum*tx1 - end if + endif enddo enddo - end if + endif ! if (lprnt) write(0,*)' tlat7=',tlat(1,:)*deltat @@ -4249,7 +4249,7 @@ subroutine micro_mg_tend ( & lamcrad(i,k) = zero pgamrad(i,k) = zero effc_fn(i,k) = ten - end if + endif enddo enddo ! recalculate 'final' rain size distribution parameters @@ -4266,9 +4266,9 @@ subroutine micro_mg_tend ( & if (dum /= dumnr(i,k)) then ! adjust number conc if needed to keep mean size in reasonable range nrtend(i,k) = (dumnr(i,k)*precip_frac(i,k)-nr(i,k)) *oneodt - end if + endif - end if + endif enddo enddo ! recalculate 'final' snow size distribution parameters @@ -4290,10 +4290,10 @@ subroutine micro_mg_tend ( & tx1 = (two*pi*1.e-2_r8) / (lams(i,k)*lams(i,k)*lams(i,k)) sadsnow(i,k) = tx1*dumns0*rho(i,k) ! m2/m3 -> cm2/cm3 - end if + endif - end do ! vertical k loop + enddo ! vertical k loop enddo do k=1,nlev do i=1,mgncol @@ -4307,9 +4307,9 @@ subroutine micro_mg_tend ( & if (qg(i,k)+qgtend(i,k)*deltat < qsmall) ngtend(i,k) = -ng(i,k) * oneodt !--ag - end do + enddo - end do + enddo ! DO STUFF FOR OUTPUT: !================================================== @@ -4459,10 +4459,10 @@ subroutine micro_mg_tend ( & else acsrfl(i,k) = zero fcsrfl(i,k) = zero - end if + endif - end do - end do + enddo + enddo do k=1,nlev do i = 1,mgncol @@ -4520,7 +4520,7 @@ subroutine calc_rercld(lamr, n0r, lamc, pgam, qric, qcic, ncic, rercld, mgncol,n if (Atmp > zero) then rercld(i,k) = rercld(i,k) + three *(qric(i,k) + qcic(i,k)) / (four * rhow * Atmp) - end if + endif enddo enddo end subroutine calc_rercld From c41612c39f57e570e9480ced8619c8a7c518edf8 Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Tue, 30 Nov 2021 01:10:52 +0000 Subject: [PATCH 07/15] adding/removing some blanks in sfc_sice --- physics/sfc_sice.f | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/physics/sfc_sice.f b/physics/sfc_sice.f index 312c35dfa..b88178702 100644 --- a/physics/sfc_sice.f +++ b/physics/sfc_sice.f @@ -579,10 +579,10 @@ subroutine ice3lay ! !===> ... begin here ! - dt2 = delt + delt - dt4 = dt2 + dt2 - dt6 = dt2 + dt4 - dt2i = one / dt2 + dt2 = delt + delt + dt4 = dt2 + dt2 + dt6 = dt2 + dt4 + dt2i = one / dt2 do i = 1, im if (flag(i)) then From 4c62ee2baef6c168cf4b0bda360e5521b7c1f88c Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Sun, 5 Dec 2021 03:07:05 +0000 Subject: [PATCH 08/15] cleanup/fix GFS_surface_composites --- physics/GFS_surface_composites.F90 | 51 +++-------------------------- physics/GFS_surface_composites.meta | 48 --------------------------- 2 files changed, 5 insertions(+), 94 deletions(-) diff --git a/physics/GFS_surface_composites.F90 b/physics/GFS_surface_composites.F90 index 14bc48cd7..20e7a57b1 100644 --- a/physics/GFS_surface_composites.F90 +++ b/physics/GFS_surface_composites.F90 @@ -30,10 +30,8 @@ end subroutine GFS_surface_composites_pre_finalize subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, frac_grid, & flag_cice, cplflx, cplice, cplwav2atm, landfrac, lakefrac, lakedepth, oceanfrac, frland, & dry, icy, lake, use_flake, wet, hice, cice, zorlo, zorll, zorli, & - snowd, snowd_lnd, snowd_ice, tprcp, tprcp_wat, & - tprcp_lnd, tprcp_ice, uustar, uustar_wat, uustar_lnd, uustar_ice, & - weasd, weasd_lnd, weasd_ice, ep1d_ice, tsfc, tsfco, tsfcl, tsfc_wat, & - tisfc, tsurf_wat, tsurf_lnd, tsurf_ice, & + tprcp, tprcp_wat, tprcp_lnd, tprcp_ice, uustar, uustar_wat, uustar_lnd, uustar_ice, & + ep1d_ice, tsfc, tsfco, tsfcl, tsfc_wat, tisfc, tsurf_wat, tsurf_lnd, tsurf_ice, & gflx_ice, tgice, islmsk, islmsk_cice, slmsk, qss, qss_wat, qss_lnd, qss_ice, & min_lakeice, min_seaice, kdt, huge, errmsg, errflg) @@ -47,12 +45,12 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra real(kind=kind_phys), dimension(:), intent(in ) :: landfrac, lakefrac, lakedepth, oceanfrac real(kind=kind_phys), dimension(:), intent(inout) :: cice, hice real(kind=kind_phys), dimension(:), intent( out) :: frland - real(kind=kind_phys), dimension(:), intent(in ) :: snowd, tprcp, uustar, weasd, qss + real(kind=kind_phys), dimension(:), intent(in ) :: tprcp, uustar, qss real(kind=kind_phys), dimension(:), intent(inout) :: tsfc, tsfco, tsfcl, tisfc - real(kind=kind_phys), dimension(:), intent(inout) :: snowd_lnd, snowd_ice, tprcp_wat, & + real(kind=kind_phys), dimension(:), intent(inout) :: tprcp_wat, & tprcp_lnd, tprcp_ice, tsfc_wat, tsurf_wat,tsurf_lnd, tsurf_ice, & - uustar_wat, uustar_lnd, uustar_ice, weasd_lnd, weasd_ice, & + uustar_wat, uustar_lnd, uustar_ice, & qss_wat, qss_lnd, qss_ice, ep1d_ice, gflx_ice real(kind=kind_phys), intent(in ) :: tgice integer, dimension(:), intent(inout) :: islmsk, islmsk_cice @@ -218,7 +216,6 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra endif if (dry(i)) then ! Land uustar_lnd(i) = uustar(i) - weasd_lnd(i) = weasd(i) tsurf_lnd(i) = tsfcl(i) ! DH* else @@ -229,7 +226,6 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra endif if (icy(i)) then ! Ice uustar_ice(i) = uustar(i) - weasd_ice(i) = weasd(i) tsurf_ice(i) = tisfc(i) ep1d_ice(i) = zero gflx_ice(i) = zero @@ -256,43 +252,6 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra endif enddo ! - if (.not. cplflx .or. kdt == 1) then - if (frac_grid) then - do i=1,im - if (dry(i)) then - if (icy(i)) then - tem = one / (cice(i)*(one-frland(i))) - snowd_ice(i) = max(zero, (snowd(i) - snowd_lnd(i)*frland(i)) * tem) - weasd_ice(i) = max(zero, (weasd(i) - weasd_lnd(i)*frland(i)) * tem) - endif - elseif (icy(i)) then - tem = one / cice(i) - snowd_lnd(i) = zero - snowd_ice(i) = snowd(i) * tem - weasd_lnd(i) = zero - weasd_ice(i) = weasd(i) * tem - endif - enddo - else - do i=1,im - if (dry(i)) then - snowd_lnd(i) = snowd(i) - weasd_lnd(i) = weasd(i) - snowd_ice(i) = zero - weasd_ice(i) = zero - elseif (icy(i)) then - snowd_lnd(i) = zero - weasd_lnd(i) = zero - tem = one / cice(i) - snowd_ice(i) = snowd(i) * tem - weasd_ice(i) = weasd(i) * tem - endif - enddo - endif - endif - -! write(0,*)' minmax of ice snow=',minval(snowd_ice),maxval(snowd_ice) - end subroutine GFS_surface_composites_pre_run end module GFS_surface_composites_pre diff --git a/physics/GFS_surface_composites.meta b/physics/GFS_surface_composites.meta index fde52ed23..35f5552e9 100644 --- a/physics/GFS_surface_composites.meta +++ b/physics/GFS_surface_composites.meta @@ -185,30 +185,6 @@ type = real kind = kind_phys intent = inout -[snowd] - standard_name = lwe_surface_snow - long_name = water equivalent snow depth - units = mm - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = in -[snowd_lnd] - standard_name = surface_snow_thickness_water_equivalent_over_land - long_name = water equivalent snow depth over land - units = mm - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = inout -[snowd_ice] - standard_name = surface_snow_thickness_water_equivalent_over_ice - long_name = water equivalent snow depth over ice - units = mm - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = inout [tprcp] standard_name = nonnegative_lwe_thickness_of_precipitation_amount_on_dynamics_timestep long_name = total precipitation amount in each time step @@ -273,30 +249,6 @@ type = real kind = kind_phys intent = inout -[weasd] - standard_name = lwe_thickness_of_surface_snow_amount - long_name = water equiv of acc snow depth over land and sea ice - units = mm - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = in -[weasd_lnd] - standard_name = water_equivalent_accumulated_snow_depth_over_land - long_name = water equiv of acc snow depth over land - units = mm - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = inout -[weasd_ice] - standard_name = water_equivalent_accumulated_snow_depth_over_ice - long_name = water equiv of acc snow depth over ice - units = mm - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = inout [ep1d_ice] standard_name = surface_upward_potential_latent_heat_flux_over_ice long_name = surface upward potential latent heat flux over ice From 7aa8caa4060e50df5d9807eae3d75bac3bfa6c86 Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Sat, 11 Dec 2021 01:19:40 +0000 Subject: [PATCH 09/15] some cosmetic change in a routine --- physics/cires_ugwp.F90 | 111 +++++++++++++++++++++-------------------- 1 file changed, 56 insertions(+), 55 deletions(-) diff --git a/physics/cires_ugwp.F90 b/physics/cires_ugwp.F90 index 6efce96f5..c4f0a255d 100644 --- a/physics/cires_ugwp.F90 +++ b/physics/cires_ugwp.F90 @@ -241,54 +241,54 @@ subroutine cires_ugwp_run(do_ugwp, me, master, im, levs, ntrac, dtp, kdt, lonr ! wrap everything in a do_ugwp 'if test' in order not to break the namelist functionality if (do_ugwp) then ! calling revised old GFS gravity wave drag - ! topo paras - ! w/ orographic effects - if(nmtvr == 14)then - ! calculate sgh30 for TOFD - sgh30 = abs(oro - oro_uf) - ! w/o orographic effects - else - sgh30 = 0. - endif - - zlwb(:) = 0. - - call GWDPS_V0(im, levs, lonr, do_tofd, Pdvdt, Pdudt, Pdtdt, Pkdis, & - ugrs, vgrs, tgrs, qgrs(:,:,1), kpbl, prsi,del,prsl, prslk, phii, phil, & - dtp, kdt, sgh30, hprime, oc, oa4, clx, theta, sigma, gamma, elvmax, & - dusfcg, dvsfcg, xlat_d, sinlat, coslat, area, cdmbgwd(1:2), & - me, master, rdxzb, con_g, con_omega, zmtb, zogw, tau_mtb, tau_ogw, & - tau_tofd, dudt_mtb, dudt_ogw, dudt_tms) - - else ! calling old GFS gravity wave drag as is - - do k=1,levs - do i=1,im - Pdvdt(i,k) = 0.0 - Pdudt(i,k) = 0.0 - Pdtdt(i,k) = 0.0 - Pkdis(i,k) = 0.0 - enddo - enddo - - if (cdmbgwd(1) > 0.0 .or. cdmbgwd(2) > 0.0) then - call gwdps_run(im, levs, Pdvdt, Pdudt, Pdtdt, & - ugrs, vgrs, tgrs, qgrs(:,:,1), & - kpbl, prsi, del, prsl, prslk, phii, phil, dtp, kdt, & - hprime, oc, oa4, clx, theta, sigma, gamma, & - elvmax, dusfcg, dvsfcg, & - con_g, con_cp, con_rd, con_rv, lonr, & - nmtvr, cdmbgwd, me, lprnt, ipr, rdxzb, & - errmsg, errflg) - if (errflg/=0) return + ! topo paras + ! w/ orographic effects + if(nmtvr == 14)then + ! calculate sgh30 for TOFD + sgh30 = abs(oro - oro_uf) + ! w/o orographic effects + else + sgh30 = 0. endif - tau_mtb = 0.0 ; tau_ogw = 0.0 ; tau_tofd = 0.0 - if (ldiag_ugwp) then - du3dt_mtb = 0.0 ; du3dt_ogw = 0.0 ; du3dt_tms= 0.0 - end if - - endif ! do_ugwp + zlwb(:) = 0. + + call GWDPS_V0(im, levs, lonr, do_tofd, Pdvdt, Pdudt, Pdtdt, Pkdis, & + ugrs, vgrs, tgrs, qgrs(:,:,1), kpbl, prsi,del,prsl, prslk, phii, phil, & + dtp, kdt, sgh30, hprime, oc, oa4, clx, theta, sigma, gamma, elvmax, & + dusfcg, dvsfcg, xlat_d, sinlat, coslat, area, cdmbgwd(1:2), & + me, master, rdxzb, con_g, con_omega, zmtb, zogw, tau_mtb, tau_ogw, & + tau_tofd, dudt_mtb, dudt_ogw, dudt_tms) + + else ! calling old GFS gravity wave drag as is + + do k=1,levs + do i=1,im + Pdvdt(i,k) = 0.0 + Pdudt(i,k) = 0.0 + Pdtdt(i,k) = 0.0 + Pkdis(i,k) = 0.0 + enddo + enddo + + if (cdmbgwd(1) > 0.0 .or. cdmbgwd(2) > 0.0) then + call gwdps_run(im, levs, Pdvdt, Pdudt, Pdtdt, & + ugrs, vgrs, tgrs, qgrs(:,:,1), & + kpbl, prsi, del, prsl, prslk, phii, phil, dtp, kdt, & + hprime, oc, oa4, clx, theta, sigma, gamma, & + elvmax, dusfcg, dvsfcg, & + con_g, con_cp, con_rd, con_rv, lonr, & + nmtvr, cdmbgwd, me, lprnt, ipr, rdxzb, & + errmsg, errflg) + if (errflg/=0) return + endif + + tau_mtb = 0.0 ; tau_ogw = 0.0 ; tau_tofd = 0.0 + if (ldiag_ugwp) then + du3dt_mtb = 0.0 ; du3dt_ogw = 0.0 ; du3dt_tms= 0.0 + endif + + endif ! do_ugwp if(ldiag3d .and. lssav .and. .not. flag_for_gwd_generic_tend) then @@ -348,19 +348,20 @@ subroutine cires_ugwp_run(do_ugwp, me, master, im, levs, ntrac, dtp, kdt, lonr endif call fv3_ugwp_solv2_v0(im, levs, dtp, tgrs, ugrs, vgrs,qgrs(:,:,1), & - prsl, prsi, phil, xlat_d, sinlat, coslat, gw_dudt, gw_dvdt, gw_dtdt, gw_kdis, & - tau_ngw, me, master, kdt) + prsl, prsi, phil, xlat_d, sinlat, coslat, & + gw_dudt, gw_dvdt, gw_dtdt, gw_kdis, tau_ngw, & + me, master, kdt) do k=1,levs do i=1,im - gw_dtdt(i,k) = pngw*gw_dtdt(i,k)+ pogw*Pdtdt(i,k) - gw_dudt(i,k) = pngw*gw_dudt(i,k)+ pogw*Pdudt(i,k) - gw_dvdt(i,k) = pngw*gw_dvdt(i,k)+ pogw*Pdvdt(i,k) - gw_kdis(i,k) = pngw*gw_kdis(i,k)+ pogw*Pkdis(i,k) + gw_dtdt(i,k) = pngw*gw_dtdt(i,k) + pogw*Pdtdt(i,k) + gw_dudt(i,k) = pngw*gw_dudt(i,k) + pogw*Pdudt(i,k) + gw_dvdt(i,k) = pngw*gw_dvdt(i,k) + pogw*Pdvdt(i,k) + gw_kdis(i,k) = pngw*gw_kdis(i,k) + pogw*Pkdis(i,k) ! accumulation of tendencies for CCPP to replicate EMC-physics updates (!! removed in latest code commit to VLAB) - !dudt(i,k) = dudt(i,k) +gw_dudt(i,k) - !dvdt(i,k) = dvdt(i,k) +gw_dvdt(i,k) - !dtdt(i,k) = dtdt(i,k) +gw_dtdt(i,k) + !dudt(i,k) = dudt(i,k) + gw_dudt(i,k) + !dvdt(i,k) = dvdt(i,k) + gw_dvdt(i,k) + !dtdt(i,k) = dtdt(i,k) + gw_dtdt(i,k) enddo enddo From 84343332fa236e91c86da1021a2a0ed68eeb43ec Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Sat, 25 Dec 2021 01:47:33 +0000 Subject: [PATCH 10/15] testing with zorli related change --- physics/GFS_surface_composites.F90 | 2 ++ 1 file changed, 2 insertions(+) diff --git a/physics/GFS_surface_composites.F90 b/physics/GFS_surface_composites.F90 index 20e7a57b1..b2c59d32c 100644 --- a/physics/GFS_surface_composites.F90 +++ b/physics/GFS_surface_composites.F90 @@ -209,6 +209,7 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra uustar_wat(i) = uustar(i) tsfc_wat(i) = tsfco(i) tsurf_wat(i) = tsfco(i) + zorlo(i) = max(1.0e-5, min(one, zorlo(i))) ! DH* else zorlo(i) = huge @@ -229,6 +230,7 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra tsurf_ice(i) = tisfc(i) ep1d_ice(i) = zero gflx_ice(i) = zero + zorli(i) = max(1.0e-5, min(one, zorli(i))) ! DH* else zorli(i) = huge From 88fcf4cbf8b41cfd5424782acf5edf5effc5ae45 Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Thu, 6 Jan 2022 14:48:46 +0000 Subject: [PATCH 11/15] updating GFS_surface_composites to address crash in cpld_debug_p7 --- physics/GFS_surface_composites.F90 | 56 ++++++++++++++++++++++++----- physics/GFS_surface_composites.meta | 55 ++++++++++++++++++++++++---- 2 files changed, 95 insertions(+), 16 deletions(-) diff --git a/physics/GFS_surface_composites.F90 b/physics/GFS_surface_composites.F90 index b2c59d32c..1175dac41 100644 --- a/physics/GFS_surface_composites.F90 +++ b/physics/GFS_surface_composites.F90 @@ -30,27 +30,29 @@ end subroutine GFS_surface_composites_pre_finalize subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, frac_grid, & flag_cice, cplflx, cplice, cplwav2atm, landfrac, lakefrac, lakedepth, oceanfrac, frland, & dry, icy, lake, use_flake, wet, hice, cice, zorlo, zorll, zorli, & - tprcp, tprcp_wat, tprcp_lnd, tprcp_ice, uustar, uustar_wat, uustar_lnd, uustar_ice, & - ep1d_ice, tsfc, tsfco, tsfcl, tsfc_wat, tisfc, tsurf_wat, tsurf_lnd, tsurf_ice, & + snowd, snowd_lnd, snowd_ice, tprcp, tprcp_wat, & + tprcp_lnd, tprcp_ice, uustar, uustar_wat, uustar_lnd, uustar_ice, & + weasd, weasd_lnd, weasd_ice, ep1d_ice, tsfc, tsfco, tsfcl, tsfc_wat, & + tisfc, tsurf_wat, tsurf_lnd, tsurf_ice, & gflx_ice, tgice, islmsk, islmsk_cice, slmsk, qss, qss_wat, qss_lnd, qss_ice, & - min_lakeice, min_seaice, kdt, huge, errmsg, errflg) + min_lakeice, min_seaice, huge, errmsg, errflg) implicit none ! Interface variables - integer, intent(in ) :: im, lkm, kdt + integer, intent(in ) :: im, lkm logical, intent(in ) :: flag_init, flag_restart, frac_grid, cplflx, cplice, cplwav2atm logical, dimension(:), intent(inout) :: flag_cice logical, dimension(:), intent(inout) :: dry, icy, lake, use_flake, wet real(kind=kind_phys), dimension(:), intent(in ) :: landfrac, lakefrac, lakedepth, oceanfrac real(kind=kind_phys), dimension(:), intent(inout) :: cice, hice real(kind=kind_phys), dimension(:), intent( out) :: frland - real(kind=kind_phys), dimension(:), intent(in ) :: tprcp, uustar, qss + real(kind=kind_phys), dimension(:), intent(in ) :: snowd, tprcp, uustar, weasd, qss real(kind=kind_phys), dimension(:), intent(inout) :: tsfc, tsfco, tsfcl, tisfc - real(kind=kind_phys), dimension(:), intent(inout) :: tprcp_wat, & + real(kind=kind_phys), dimension(:), intent(inout) :: snowd_lnd, snowd_ice, tprcp_wat, & tprcp_lnd, tprcp_ice, tsfc_wat, tsurf_wat,tsurf_lnd, tsurf_ice, & - uustar_wat, uustar_lnd, uustar_ice, & + uustar_wat, uustar_lnd, uustar_ice, weasd_lnd, weasd_ice, & qss_wat, qss_lnd, qss_ice, ep1d_ice, gflx_ice real(kind=kind_phys), intent(in ) :: tgice integer, dimension(:), intent(inout) :: islmsk, islmsk_cice @@ -62,6 +64,7 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra real(kind=kind_phys), parameter :: timin = 173.0_kind_phys ! minimum temperature allowed for snow/ice real(kind=kind_phys) :: tem + logical :: icy_old(im) ! CCPP error handling character(len=*), intent(out) :: errmsg @@ -74,6 +77,9 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra errmsg = '' errflg = 0 + do i=1,im + icy_old(i) = icy(i) + enddo if (frac_grid) then ! cice is ice fraction wrt water area do i=1,im frland(i) = landfrac(i) @@ -209,7 +215,7 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra uustar_wat(i) = uustar(i) tsfc_wat(i) = tsfco(i) tsurf_wat(i) = tsfco(i) - zorlo(i) = max(1.0e-5, min(one, zorlo(i))) + zorlo(i) = max(1.0e-5, min(one, zorlo(i))) ! DH* else zorlo(i) = huge @@ -217,6 +223,7 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra endif if (dry(i)) then ! Land uustar_lnd(i) = uustar(i) + weasd_lnd(i) = weasd(i) tsurf_lnd(i) = tsfcl(i) ! DH* else @@ -227,10 +234,11 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra endif if (icy(i)) then ! Ice uustar_ice(i) = uustar(i) + weasd_ice(i) = weasd(i) tsurf_ice(i) = tisfc(i) ep1d_ice(i) = zero gflx_ice(i) = zero - zorli(i) = max(1.0e-5, min(one, zorli(i))) + zorli(i) = max(1.0e-5, min(one, zorli(i))) ! DH* else zorli(i) = huge @@ -254,6 +262,36 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra endif enddo ! + if (frac_grid) then + do i=1,im + if (dry(i)) then + if (icy(i) .and. .not. icy_old(i)) then + tem = one / (cice(i)*(one-frland(i))) + snowd_ice(i) = max(zero, (snowd(i) - snowd_lnd(i)*frland(i)) * tem) + weasd_ice(i) = max(zero, (weasd(i) - weasd_lnd(i)*frland(i)) * tem) + endif + elseif (icy(i) .and. .not. icy_old(i)) then + tem = one / cice(i) + snowd_lnd(i) = zero + snowd_ice(i) = snowd(i) * tem + weasd_lnd(i) = zero + weasd_ice(i) = weasd(i) * tem + endif + enddo + else + do i=1,im + if (icy(i) .and. .not. icy_old(i)) then + snowd_lnd(i) = zero + weasd_lnd(i) = zero + tem = one / cice(i) + snowd_ice(i) = snowd(i) * tem + weasd_ice(i) = weasd(i) * tem + endif + enddo + endif + +! write(0,*)' minmax of ice snow=',minval(snowd_ice),maxval(snowd_ice) + end subroutine GFS_surface_composites_pre_run end module GFS_surface_composites_pre diff --git a/physics/GFS_surface_composites.meta b/physics/GFS_surface_composites.meta index 35f5552e9..2f7b451b1 100644 --- a/physics/GFS_surface_composites.meta +++ b/physics/GFS_surface_composites.meta @@ -185,6 +185,30 @@ type = real kind = kind_phys intent = inout +[snowd] + standard_name = lwe_surface_snow + long_name = water equivalent snow depth + units = mm + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in +[snowd_lnd] + standard_name = surface_snow_thickness_water_equivalent_over_land + long_name = water equivalent snow depth over land + units = mm + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = inout +[snowd_ice] + standard_name = surface_snow_thickness_water_equivalent_over_ice + long_name = water equivalent snow depth over ice + units = mm + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = inout [tprcp] standard_name = nonnegative_lwe_thickness_of_precipitation_amount_on_dynamics_timestep long_name = total precipitation amount in each time step @@ -249,6 +273,30 @@ type = real kind = kind_phys intent = inout +[weasd] + standard_name = lwe_thickness_of_surface_snow_amount + long_name = water equiv of acc snow depth over land and sea ice + units = mm + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in +[weasd_lnd] + standard_name = water_equivalent_accumulated_snow_depth_over_land + long_name = water equiv of acc snow depth over land + units = mm + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = inout +[weasd_ice] + standard_name = water_equivalent_accumulated_snow_depth_over_ice + long_name = water equiv of acc snow depth over ice + units = mm + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = inout [ep1d_ice] standard_name = surface_upward_potential_latent_heat_flux_over_ice long_name = surface upward potential latent heat flux over ice @@ -407,13 +455,6 @@ type = real kind = kind_phys intent = in -[kdt] - standard_name = index_of_timestep - long_name = current forecast iteration - units = index - dimensions = () - type = integer - intent = in [huge] standard_name = netcdf_float_fillvalue long_name = definition of NetCDF float FillValue From 01a91fa9f03b4e956a18423ccde6d6861ad8f682 Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Thu, 6 Jan 2022 16:47:34 +0000 Subject: [PATCH 12/15] reverting some change --- physics/GFS_surface_composites.F90 | 60 ++++++++++++++--------------- physics/GFS_surface_composites.meta | 7 ++++ 2 files changed, 36 insertions(+), 31 deletions(-) diff --git a/physics/GFS_surface_composites.F90 b/physics/GFS_surface_composites.F90 index 1175dac41..057ec4946 100644 --- a/physics/GFS_surface_composites.F90 +++ b/physics/GFS_surface_composites.F90 @@ -35,12 +35,12 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra weasd, weasd_lnd, weasd_ice, ep1d_ice, tsfc, tsfco, tsfcl, tsfc_wat, & tisfc, tsurf_wat, tsurf_lnd, tsurf_ice, & gflx_ice, tgice, islmsk, islmsk_cice, slmsk, qss, qss_wat, qss_lnd, qss_ice, & - min_lakeice, min_seaice, huge, errmsg, errflg) + min_lakeice, min_seaice, kdt, huge, errmsg, errflg) implicit none ! Interface variables - integer, intent(in ) :: im, lkm + integer, intent(in ) :: im, lkm, kdt logical, intent(in ) :: flag_init, flag_restart, frac_grid, cplflx, cplice, cplwav2atm logical, dimension(:), intent(inout) :: flag_cice logical, dimension(:), intent(inout) :: dry, icy, lake, use_flake, wet @@ -64,7 +64,6 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra real(kind=kind_phys), parameter :: timin = 173.0_kind_phys ! minimum temperature allowed for snow/ice real(kind=kind_phys) :: tem - logical :: icy_old(im) ! CCPP error handling character(len=*), intent(out) :: errmsg @@ -77,9 +76,6 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra errmsg = '' errflg = 0 - do i=1,im - icy_old(i) = icy(i) - enddo if (frac_grid) then ! cice is ice fraction wrt water area do i=1,im frland(i) = landfrac(i) @@ -262,32 +258,34 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra endif enddo ! - if (frac_grid) then - do i=1,im - if (dry(i)) then - if (icy(i) .and. .not. icy_old(i)) then - tem = one / (cice(i)*(one-frland(i))) - snowd_ice(i) = max(zero, (snowd(i) - snowd_lnd(i)*frland(i)) * tem) - weasd_ice(i) = max(zero, (weasd(i) - weasd_lnd(i)*frland(i)) * tem) + if (.not. cplflx .or. kdt == 1) then + if (frac_grid) then + do i=1,im + if (dry(i)) then + if (icy(i)) then + tem = one / (cice(i)*(one-frland(i))) + snowd_ice(i) = max(zero, (snowd(i) - snowd_lnd(i)*frland(i)) * tem) + weasd_ice(i) = max(zero, (weasd(i) - weasd_lnd(i)*frland(i)) * tem) + endif + elseif (icy(i)) then + tem = one / cice(i) + snowd_lnd(i) = zero + snowd_ice(i) = snowd(i) * tem + weasd_lnd(i) = zero + weasd_ice(i) = weasd(i) * tem endif - elseif (icy(i) .and. .not. icy_old(i)) then - tem = one / cice(i) - snowd_lnd(i) = zero - snowd_ice(i) = snowd(i) * tem - weasd_lnd(i) = zero - weasd_ice(i) = weasd(i) * tem - endif - enddo - else - do i=1,im - if (icy(i) .and. .not. icy_old(i)) then - snowd_lnd(i) = zero - weasd_lnd(i) = zero - tem = one / cice(i) - snowd_ice(i) = snowd(i) * tem - weasd_ice(i) = weasd(i) * tem - endif - enddo + enddo + else + do i=1,im + if (icy(i)) then + snowd_lnd(i) = zero + weasd_lnd(i) = zero + tem = one / cice(i) + snowd_ice(i) = snowd(i) * tem + weasd_ice(i) = weasd(i) * tem + endif + enddo + endif endif ! write(0,*)' minmax of ice snow=',minval(snowd_ice),maxval(snowd_ice) diff --git a/physics/GFS_surface_composites.meta b/physics/GFS_surface_composites.meta index 2f7b451b1..fde52ed23 100644 --- a/physics/GFS_surface_composites.meta +++ b/physics/GFS_surface_composites.meta @@ -455,6 +455,13 @@ type = real kind = kind_phys intent = in +[kdt] + standard_name = index_of_timestep + long_name = current forecast iteration + units = index + dimensions = () + type = integer + intent = in [huge] standard_name = netcdf_float_fillvalue long_name = definition of NetCDF float FillValue From 30a7afd6973547a78b5f1f8524e72b16c7674f71 Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Thu, 6 Jan 2022 17:32:30 +0000 Subject: [PATCH 13/15] still more mods --- physics/GFS_surface_composites.F90 | 27 +++++++++++++++------------ 1 file changed, 15 insertions(+), 12 deletions(-) diff --git a/physics/GFS_surface_composites.F90 b/physics/GFS_surface_composites.F90 index 057ec4946..ae30aa21f 100644 --- a/physics/GFS_surface_composites.F90 +++ b/physics/GFS_surface_composites.F90 @@ -258,34 +258,37 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra endif enddo ! - if (.not. cplflx .or. kdt == 1) then - if (frac_grid) then - do i=1,im - if (dry(i)) then - if (icy(i)) then + if (frac_grid) then + do i=1,im + if (dry(i)) then + if (icy(i)) then + if (kdt == 1 .or. (.not. cplflx .and. lakefrac(i) > zero)) then tem = one / (cice(i)*(one-frland(i))) snowd_ice(i) = max(zero, (snowd(i) - snowd_lnd(i)*frland(i)) * tem) weasd_ice(i) = max(zero, (weasd(i) - weasd_lnd(i)*frland(i)) * tem) endif - elseif (icy(i)) then + endif + elseif (icy(i)) then + if (kdt == 1 .or. (.not. cplflx .and. lakefrac(i) > zero)) then tem = one / cice(i) snowd_lnd(i) = zero snowd_ice(i) = snowd(i) * tem weasd_lnd(i) = zero weasd_ice(i) = weasd(i) * tem endif - enddo - else - do i=1,im - if (icy(i)) then + endif + enddo + do i=1,im + if (icy(i)) then + if (kdt == 1 .or. (.not. cplflx .and. lakefrac(i) > zero)) then snowd_lnd(i) = zero weasd_lnd(i) = zero tem = one / cice(i) snowd_ice(i) = snowd(i) * tem weasd_ice(i) = weasd(i) * tem endif - enddo - endif + endif + enddo endif ! write(0,*)' minmax of ice snow=',minval(snowd_ice),maxval(snowd_ice) From 7a708756e8f4cde4166191c72f559e1c847db33a Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Thu, 6 Jan 2022 17:44:25 +0000 Subject: [PATCH 14/15] yet more mods --- physics/GFS_surface_composites.F90 | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/physics/GFS_surface_composites.F90 b/physics/GFS_surface_composites.F90 index ae30aa21f..bd61a0c4d 100644 --- a/physics/GFS_surface_composites.F90 +++ b/physics/GFS_surface_composites.F90 @@ -262,14 +262,14 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra do i=1,im if (dry(i)) then if (icy(i)) then - if (kdt == 1 .or. (.not. cplflx .and. lakefrac(i) > zero)) then + if (kdt == 1 .or. (.not. cplflx .or. lakefrac(i) > zero)) then tem = one / (cice(i)*(one-frland(i))) snowd_ice(i) = max(zero, (snowd(i) - snowd_lnd(i)*frland(i)) * tem) weasd_ice(i) = max(zero, (weasd(i) - weasd_lnd(i)*frland(i)) * tem) endif endif elseif (icy(i)) then - if (kdt == 1 .or. (.not. cplflx .and. lakefrac(i) > zero)) then + if (kdt == 1 .or. (.not. cplflx .or. lakefrac(i) > zero)) then tem = one / cice(i) snowd_lnd(i) = zero snowd_ice(i) = snowd(i) * tem @@ -280,7 +280,7 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra enddo do i=1,im if (icy(i)) then - if (kdt == 1 .or. (.not. cplflx .and. lakefrac(i) > zero)) then + if (kdt == 1 .or. (.not. cplflx .or. lakefrac(i) > zero)) then snowd_lnd(i) = zero weasd_lnd(i) = zero tem = one / cice(i) From 8d97f46a76c2f185cc365b6a380024bcc4e0efd7 Mon Sep 17 00:00:00 2001 From: "Shrinivas.Moorthi" Date: Fri, 7 Jan 2022 00:34:10 +0000 Subject: [PATCH 15/15] fix a typo --- physics/GFS_surface_composites.F90 | 1 + 1 file changed, 1 insertion(+) diff --git a/physics/GFS_surface_composites.F90 b/physics/GFS_surface_composites.F90 index 14fc3b520..510b3f427 100644 --- a/physics/GFS_surface_composites.F90 +++ b/physics/GFS_surface_composites.F90 @@ -278,6 +278,7 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, fra endif endif enddo + else do i=1,im if (icy(i)) then if (kdt == 1 .or. (.not. cplflx .or. lakefrac(i) > zero)) then