Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Bug fix in GFS_surface_composites.F90 from Moorthi #520

Closed
Closed
Show file tree
Hide file tree
Changes from 1 commit
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
171 changes: 101 additions & 70 deletions physics/GFS_surface_composites.F90
Original file line number Diff line number Diff line change
Expand Up @@ -24,15 +24,15 @@ end subroutine GFS_surface_composites_pre_finalize
!> \section arg_table_GFS_surface_composites_pre_run Argument Table
!! \htmlinclude GFS_surface_composites_pre_run.html
!!
subroutine GFS_surface_composites_pre_run (im, lkm, frac_grid, flag_cice, cplflx, cplwav2atm, &
landfrac, lakefrac, lakedepth, oceanfrac, frland, &
dry, icy, lake, ocean, wet, cice, cimin, zorl, zorlo, zorll, zorli, zorl_wat, &
zorl_lnd, zorl_ice, snowd, snowd_wat, snowd_lnd, snowd_ice, tprcp, tprcp_wat, &
tprcp_lnd, tprcp_ice, uustar, uustar_wat, uustar_lnd, uustar_ice, &
weasd, weasd_wat, weasd_lnd, weasd_ice, ep1d_ice, tsfc, tsfco, tsfcl, tsfc_wat,&
tsfc_lnd, tsfc_ice, tisfc, tice, tsurf, tsurf_wat, tsurf_lnd, tsurf_ice, &
gflx_ice, tgice, islmsk, semis_rad, semis_wat, semis_lnd, semis_ice, &
qss, qss_wat, qss_lnd, qss_ice, hflx, hflx_wat, hflx_lnd, hflx_ice, &
subroutine GFS_surface_composites_pre_run (im, lkm, frac_grid, flag_cice, cplflx, cplwav2atm, &
landfrac, lakefrac, lakedepth, oceanfrac, frland, &
dry, icy, lake, ocean, wet, hice, cice, zorl, zorlo, zorll, zorli, zorl_wat, &
zorl_lnd, zorl_ice, snowd, snowd_wat, snowd_lnd, snowd_ice, tprcp, tprcp_wat, &
tprcp_lnd, tprcp_ice, uustar, uustar_wat, uustar_lnd, uustar_ice, &
weasd, weasd_wat, weasd_lnd, weasd_ice, ep1d_ice, tsfc, tsfco, tsfcl, tsfc_wat, &
tsfc_lnd, tsfc_ice, tisfc, tice, tsurf, tsurf_wat, tsurf_lnd, tsurf_ice, &
gflx_ice, tgice, islmsk, islmsk_cice, slmsk, semis_rad, semis_wat, semis_lnd, semis_ice, &
qss, qss_wat, qss_lnd, qss_ice, hflx, hflx_wat, hflx_lnd, hflx_ice, &
min_lakeice, min_seaice, errmsg, errflg)

implicit none
Expand All @@ -42,9 +42,8 @@ subroutine GFS_surface_composites_pre_run (im, lkm, frac_grid, flag_cice, cplflx
logical, intent(in ) :: frac_grid, cplflx, cplwav2atm
logical, dimension(im), intent(inout) :: flag_cice
logical, dimension(im), intent(inout) :: dry, icy, lake, ocean, wet
real(kind=kind_phys), intent(in ) :: cimin
real(kind=kind_phys), dimension(im), intent(in ) :: landfrac, lakefrac, lakedepth, oceanfrac
real(kind=kind_phys), dimension(im), intent(inout) :: cice
real(kind=kind_phys), dimension(im), intent(inout) :: cice, hice
real(kind=kind_phys), dimension(im), intent( out) :: frland
real(kind=kind_phys), dimension(im), intent(in ) :: zorl, snowd, tprcp, uustar, weasd, qss, hflx

Expand All @@ -55,11 +54,13 @@ subroutine GFS_surface_composites_pre_run (im, lkm, frac_grid, flag_cice, cplflx
qss_wat, qss_lnd, qss_ice, hflx_wat, hflx_lnd, hflx_ice, ep1d_ice, gflx_ice
real(kind=kind_phys), dimension(im), intent( out) :: tice
real(kind=kind_phys), intent(in ) :: tgice
integer, dimension(im), intent(inout) :: islmsk
integer, dimension(im), intent(inout) :: islmsk, islmsk_cice
real(kind=kind_phys), dimension(im), intent(in ) :: semis_rad
real(kind=kind_phys), dimension(im), intent(inout) :: semis_wat, semis_lnd, semis_ice
real(kind=kind_phys), dimension(im), intent(inout) :: semis_wat, semis_lnd, semis_ice, slmsk
real(kind=kind_phys), intent(in ) :: min_lakeice, min_seaice

real(kind=kind_phys), parameter :: timin = 173.0_kind_phys ! minimum temperature allowed for snow/ice

! CCPP error handling
character(len=*), intent(out) :: errmsg
integer, intent(out) :: errflg
Expand All @@ -76,37 +77,49 @@ subroutine GFS_surface_composites_pre_run (im, lkm, frac_grid, flag_cice, cplflx
frland(i) = landfrac(i)
if (frland(i) > zero) dry(i) = .true.
if (frland(i) < one) then
if (flag_cice(i)) then
if (oceanfrac(i) > zero) then
if (cice(i) >= min_seaice) then
icy(i) = .true.
if (cice(i) < one) wet(i) = .true. ! some open ocean/lake water exists
tisfc(i) = max(timin, min(tisfc(i), tgice))
if (cplflx) then
islmsk_cice(i) = 4
flag_cice(i) = .true.
else
islmsk_cice(i) = 2
endif
islmsk(i) = 2
else
cice(i) = zero
hice(i) = zero
flag_cice(i) = .false.
! islmsk_cice(i) = 0
! islmsk(i) = 0
wet(i) = .true. ! some open ocean/lake water exists
islmsk_cice(i) = 0
islmsk(i) = 0
endif
if (cice(i) < one) then
wet(i) = .true. ! some open ocean
if (.not. cplflx .and. icy(i)) tsfco(i) = max(tisfc(i), tgice)
endif
else
if (cice(i) >= min_lakeice) then
icy(i) = .true.
if (cice(i) < one) wet(i) = .true. ! some open ocean/lake water exists
islmsk(i) = 2
tisfc(i) = max(timin, min(tisfc(i), tgice))
else
cice(i) = zero
! islmsk(i) = 0
wet(i) = .true. ! some open ocean/lake water exists
hice(i) = zero
islmsk(i) = 0
endif
endif
if (wet(i) .and. .not. cplflx) then
if (oceanfrac(i) > zero) then
tsfco(i) = max(tsfco(i), tisfc(i), tgice)
elseif (icy(i)) then
tsfco(i) = max(tisfc(i), tgice)
islmsk_cice(i) = islmsk(i)
if (cice(i) < one) then
wet(i) = .true. ! some open lake
if (icy(i)) tsfco(i) = max(tisfc(i), tgice)
endif
endif
else
else ! all land
cice(i) = zero
hice(i) = zero
islmsk_cice(i) = 1
islmsk(i) = 1
endif
enddo

Expand All @@ -118,27 +131,39 @@ subroutine GFS_surface_composites_pre_run (im, lkm, frac_grid, flag_cice, cplflx
dry(i) = .true.
frland(i) = one
cice(i) = zero
hice(i) = zero
else
frland(i) = zero
if (flag_cice(i)) then
if (cice(i) > min_seaice) then
icy(i) = .true.
if (oceanfrac(i) > zero) then
if (cice(i) >= min_seaice) then
icy(i) = .true.
tisfc(i) = max(timin, min(tisfc(i), tgice))
else
cice(i) = zero
hice(i) = zero
flag_cice(i) = .false.
islmsk(i) = 0
islmsk_cice(i) = 0
endif
if (cice(i) < one) then
wet(i) = .true. ! some open ocean
if (.not. cplflx .and. icy(i)) tsfco(i) = max(tisfc(i), tgice)
endif
else
if (cice(i) > min_lakeice) then
if (cice(i) >= min_lakeice) then
icy(i) = .true.
tisfc(i) = max(timin, min(tisfc(i), tgice))
else
cice(i) = zero
hice(i) = zero
flag_cice(i) = .false.
islmsk(i) = 0
endif
endif
if (cice(i) < one) then
wet(i) = .true. ! some open ocean/lake water exists
if (.not. cplflx .and. icy(i)) tsfco(i) = max(tisfc(i), tgice)
islmsk_cice(i) = islmsk(i)
if (cice(i) < one) then
wet(i) = .true. ! some open lake
if (icy(i)) tsfco(i) = max(tisfc(i), tgice)
endif
endif
endif
enddo
Expand All @@ -162,17 +187,17 @@ subroutine GFS_surface_composites_pre_run (im, lkm, frac_grid, flag_cice, cplflx
tprcp_lnd(i) = tprcp(i)
tprcp_ice(i) = tprcp(i)
if (wet(i)) then ! Water
uustar_wat(i) = uustar(i)
! uustar_wat(i) = uustar(i)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why is this line commented out? Does it break anything keeping it? Some of this was added specifically for other physics options.

zorl_wat(i) = zorlo(i)
tsfc_wat(i) = tsfco(i)
tsurf_wat(i) = tsfco(i)
! weasd_wat(i) = weasd(i)
! snowd_wat(i) = snowd(i)
weasd_wat(i) = zero
snowd_wat(i) = zero
semis_wat(i) = 0.984d0
qss_wat(i) = qss(i)
hflx_wat(i) = hflx(i)
semis_wat(i) = 0.984_kind_phys
! qss_wat(i) = qss(i)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why are these lines commented out? They were added specifically for MYNN, unless they break your results then please do not change this. All qss_* and hflx_* are required by MYNN and possibly RUC LSM, as far as I remember. uustar_wat possibly as well.

! hflx_wat(i) = hflx(i)
endif
if (dry(i)) then ! Land
uustar_lnd(i) = uustar(i)
Expand All @@ -182,8 +207,8 @@ subroutine GFS_surface_composites_pre_run (im, lkm, frac_grid, flag_cice, cplflx
tsurf_lnd(i) = tsfcl(i)
snowd_lnd(i) = snowd(i)
semis_lnd(i) = semis_rad(i)
qss_lnd(i) = qss(i)
hflx_lnd(i) = hflx(i)
! qss_lnd(i) = qss(i)
! hflx_lnd(i) = hflx(i)
end if
if (icy(i)) then ! Ice
uustar_ice(i) = uustar(i)
Expand All @@ -195,9 +220,10 @@ subroutine GFS_surface_composites_pre_run (im, lkm, frac_grid, flag_cice, cplflx
ep1d_ice(i) = zero
gflx_ice(i) = zero
semis_ice(i) = 0.95_kind_phys
qss_ice(i) = qss(i)
hflx_ice(i) = hflx(i)
! qss_ice(i) = qss(i)
! hflx_ice(i) = hflx(i)
endif
if (nint(slmsk(i)) /= 1) slmsk(i) = islmsk(i)
enddo

! to prepare to separate lake from ocean under water category
Expand Down Expand Up @@ -364,7 +390,7 @@ subroutine GFS_surface_composites_post_run (

! Local variables
integer :: i, k
real(kind=kind_phys) :: txl, txi, txo, tem
real(kind=kind_phys) :: txl, txi, txo, wfrac

! Initialize CCPP error handling variables
errmsg = ''
Expand All @@ -377,9 +403,10 @@ subroutine GFS_surface_composites_post_run (
do i=1, im

! Three-way composites (fields from sfc_diff)
txl = landfrac(i)
txi = cice(i)*(one - txl) ! txi = ice fraction wrt whole cell
txo = max(zero, one - txl - txi)
txl = landfrac(i) ! land fraction
wfrac = one - txl ! ocean fraction
txi = cice(i) * wfrac ! txi = ice fraction wrt whole cell
txo = max(zero, wfrac-txi) ! txo = open water fraction

zorl(i) = txl*zorl_lnd(i) + txi*zorl_ice(i) + txo*zorl_wat(i)
cd(i) = txl*cd_lnd(i) + txi*cd_ice(i) + txo*cd_wat(i)
Expand All @@ -404,11 +431,10 @@ subroutine GFS_surface_composites_post_run (
!tprcp(i) = txl*tprcp_lnd(i) + txi*tprcp_ice(i) + txo*tprcp_wat(i)

if (.not. flag_cice(i) .and. islmsk(i) == 2) then
tem = one - txl
evap(i) = txl*evap_lnd(i) + tem*evap_ice(i)
hflx(i) = txl*hflx_lnd(i) + tem*hflx_ice(i)
qss(i) = txl*qss_lnd(i) + tem*qss_ice(i)
gflx(i) = txl*gflx_lnd(i) + tem*gflx_ice(i)
evap(i) = txl*evap_lnd(i) + wfrac*evap_ice(i)
hflx(i) = txl*hflx_lnd(i) + wfrac*hflx_ice(i)
qss(i) = txl*qss_lnd(i) + wfrac*qss_ice(i)
gflx(i) = txl*gflx_lnd(i) + wfrac*gflx_ice(i)
else
evap(i) = txl*evap_lnd(i) + txi*evap_ice(i) + txo*evap_wat(i)
hflx(i) = txl*hflx_lnd(i) + txi*hflx_ice(i) + txo*hflx_wat(i)
Expand Down Expand Up @@ -451,14 +477,18 @@ subroutine GFS_surface_composites_post_run (
! tisfc(i) = tsfc_ice(i) ! over ice when uncoupled
! endif

if (.not. flag_cice(i)) then
if (islmsk(i) == 2) then ! return updated lake ice thickness & concentration to global array
tisfc(i) = tice(i)
else ! this would be over open ocean or land (no ice fraction)
hice(i) = zero
cice(i) = zero
tisfc(i) = tsfc(i)
endif
! if (.not. flag_cice(i)) then
! if (islmsk(i) == 2) then ! return updated lake ice thickness & concentration to global array
! tisfc(i) = tice(i)
! else ! this would be over open ocean or land (no ice fraction)
! hice(i) = zero
! cice(i) = zero
! tisfc(i) = tsfc(i)
! endif
! endif
if (.not. icy(i)) then
hice(i) = zero
cice(i) = zero
endif
enddo

Expand All @@ -478,6 +508,9 @@ subroutine GFS_surface_composites_post_run (
fh2(i) = fh2_lnd(i)
!tsurf(i) = tsurf_lnd(i)
tsfcl(i) = tsfc_lnd(i) ! over land
tsfc(i) = tsfcl(i)
tsfco(i) = tsfc(i)
tisfc(i) = tsfc(i)
cmm(i) = cmm_lnd(i)
chh(i) = chh_lnd(i)
gflx(i) = gflx_lnd(i)
Expand All @@ -488,11 +521,8 @@ subroutine GFS_surface_composites_post_run (
evap(i) = evap_lnd(i)
hflx(i) = hflx_lnd(i)
qss(i) = qss_lnd(i)
tsfc(i) = tsfc_lnd(i)
hice(i) = zero
cice(i) = zero
tisfc(i) = tsfc(i)
tsfco(i) = tsfc(i)
elseif (islmsk(i) == 0) then
zorl(i) = zorl_wat(i)
cd(i) = cd_wat(i)
Expand All @@ -506,7 +536,9 @@ subroutine GFS_surface_composites_post_run (
fh2(i) = fh2_wat(i)
!tsurf(i) = tsurf_wat(i)
tsfco(i) = tsfc_wat(i) ! over lake (and ocean when uncoupled)
tsfc(i) = tsfco(i)
tsfcl(i) = tsfc(i)
tisfc(i) = tsfc(i)
cmm(i) = cmm_wat(i)
chh(i) = chh_wat(i)
gflx(i) = gflx_wat(i)
Expand All @@ -517,10 +549,8 @@ subroutine GFS_surface_composites_post_run (
evap(i) = evap_wat(i)
hflx(i) = hflx_wat(i)
qss(i) = qss_wat(i)
tsfc(i) = tsfc_wat(i)
hice(i) = zero
cice(i) = zero
tisfc(i) = tsfc(i)
else ! islmsk(i) == 2
zorl(i) = zorl_ice(i)
cd(i) = cd_ice(i)
Expand All @@ -544,12 +574,13 @@ subroutine GFS_surface_composites_post_run (
evap(i) = evap_ice(i)
hflx(i) = hflx_ice(i)
qss(i) = qss_ice(i)
tisfc(i) = tice(i)
if (.not. flag_cice(i)) then
tisfc(i) = tice(i) ! over lake ice (and sea ice when uncoupled)
! tisfc(i) = tice(i) ! over lake ice (and sea ice when uncoupled)
zorl(i) = cice(i) * zorl_ice(i) + (one - cice(i)) * zorl_wat(i)
tsfc(i) = tsfc_ice(i)
tsfc(i) = tsfc_ice(i) ! over lake (and ocean when uncoupled)
elseif (wet(i)) then
if (cice(i) > min_seaice) then ! this was already done for lake ice in sfc_sice
if (cice(i) >= min_seaice) then ! this was already done for lake ice in sfc_sice
txi = cice(i)
txo = one - txi
evap(i) = txi * evap_ice(i) + txo * evap_wat(i)
Expand All @@ -576,7 +607,7 @@ subroutine GFS_surface_composites_post_run (
endif
tsfcl(i) = tsfc(i)
do k=1,kice ! store tiice in stc to reduce output in the nonfrac grid case
stc(i,k)=tiice(i,k)
stc(i,k) = tiice(i,k)
end do
endif

Expand Down
Loading