Skip to content

Commit

Permalink
Merge pull request ufs-community#8 from junwang-noaa/bugfixes
Browse files Browse the repository at this point in the history 
Bugfixes in inline post, netcdf time stamp and gfs physics for coupled system. GFS physics update to reduce the cold bias in lower atmosphere layer. VLAB tickets: #69814 and #69735. fv3atm issue: ufs-community#4, ufs-community#5, ufs-community#6 and ufs-community#7.
  • Loading branch information
@junwang-noaa
junwang-noaa authored Nov 15, 2019
2 parents 45dbc34 + 97d12bc commit b18739d
Show file tree
Hide file tree
Showing 10 changed files with 136 additions and 98 deletions.
2 changes: 1 addition & 1 deletion atmos_cubed_sphere
Submodule atmos_cubed_sphere updated from 786447 to f967c7
2 changes: 1 addition & 1 deletion ccpp/physics
Submodule physics updated from d4b1cd to c0bec7
2 changes: 1 addition & 1 deletion gfsphysics/GFS_layer/GFS_physics_driver.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1149,7 +1149,7 @@ subroutine GFS_physics_driver &
if (fice(i) < one) then
wet(i) = .true.
! Sfcprop%tsfco(i) = tgice
Sfcprop%tsfco(i) = max(Sfcprop%tisfc(i), tgice)
if (.not. Model%cplflx) Sfcprop%tsfco(i) = max(Sfcprop%tisfc(i), tgice)
! Sfcprop%tsfco(i) = max((Sfcprop%tsfc(i) - fice(i)*sfcprop%tisfc(i)) &
! / (one - fice(i)), tgice)
endif
Expand Down
11 changes: 10 additions & 1 deletion gfsphysics/GFS_layer/GFS_restart.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -117,6 +117,9 @@ subroutine GFS_restart_populate (Restart, Model, Statein, Stateout, Sfcprop, &
#endif

Restart%num3d = Model%ntot3d
if(Model%lrefres) then
Restart%num3d = Model%ntot3d+1
endif
#ifdef CCPP
! GF
if (Model%imfdeepcnv == 3) then
Expand Down Expand Up @@ -252,7 +255,13 @@ subroutine GFS_restart_populate (Restart, Model, Statein, Stateout, Sfcprop, &
Restart%data(nb,num)%var3p => Tbd(nb)%phy_f3d(:,:,num)
enddo
enddo

if (Model%lrefres) then
num = Model%ntot3d+1
restart%name3d(num) = 'ref_f3d'
do nb = 1,nblks
Restart%data(nb,num)%var3p => IntDiag(nb)%refl_10cm(:,:)
enddo
endif
#ifdef CCPP
!--- RAP/HRRR-specific variables, 3D
num = Model%ntot3d
Expand Down
12 changes: 11 additions & 1 deletion gfsphysics/GFS_layer/GFS_typedefs.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -701,6 +701,9 @@ module GFS_typedefs
!--- GFDL microphysical paramters
logical :: lgfdlmprad !< flag for GFDL mp scheme and radiation consistency

!--- Thompson,GFDL mp parameter
logical :: lrefres !< flag for radar reflectivity in restart file

!--- land/surface model parameters
integer :: lsm !< flag for land surface model lsm=1 for noah lsm
integer :: lsm_noah=1 !< flag for NOAH land surface model
Expand Down Expand Up @@ -2740,6 +2743,9 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, &
!--- GFDL microphysical parameters
logical :: lgfdlmprad = .false. !< flag for GFDLMP radiation interaction

!--- Thompson,GFDL microphysical parameter
logical :: lrefres = .false. !< flag for radar reflectivity in restart file

!--- land/surface model parameters
integer :: lsm = 1 !< flag for land surface model to use =0 for osu lsm; =1 for noah lsm; =2 for noah mp lsm; =3 for RUC lsm
integer :: lsoil = 4 !< number of soil layers
Expand Down Expand Up @@ -3023,7 +3029,7 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, &
mg_do_graupel, mg_do_hail, mg_nccons, mg_nicons, mg_ngcons, &
mg_ncnst, mg_ninst, mg_ngnst, sed_supersat, do_sb_physics, &
mg_alf, mg_qcmin, mg_do_ice_gmao, mg_do_liq_liu, &
ltaerosol, lradar, ttendlim, lgfdlmprad, &
ltaerosol, lradar, lrefres, ttendlim, lgfdlmprad, &
!--- max hourly
avg_max_length, &
!--- land/surface model control
Expand Down Expand Up @@ -3312,6 +3318,8 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, &
Model%ttendlim = ttendlim
!--- gfdl MP parameters
Model%lgfdlmprad = lgfdlmprad
!--- Thompson,GFDL MP parameter
Model%lrefres = lrefres

!--- land/surface model parameters
Model%lsm = lsm
Expand Down Expand Up @@ -4310,6 +4318,7 @@ subroutine control_print(Model)
print *, ' Thompson microphysical parameters'
print *, ' ltaerosol : ', Model%ltaerosol
print *, ' lradar : ', Model%lradar
print *, ' lrefres : ', Model%lrefres
print *, ' ttendlim : ', Model%ttendlim
print *, ' '
endif
Expand All @@ -4327,6 +4336,7 @@ subroutine control_print(Model)
if (Model%imp_physics == Model%imp_physics_gfdl) then
print *, ' GFDL microphysical parameters'
print *, ' GFDL MP radiation inter: ', Model%lgfdlmprad
print *, ' lrefres : ', Model%lrefres
print *, ' '
endif

Expand Down
8 changes: 4 additions & 4 deletions gfsphysics/physics/gfdl_cloud_microphys.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4688,7 +4688,7 @@ subroutine cloud_diagnosis (is, ie, ks, ke, den, delp, lsm, qmw, qmi, qmr, qms,
real :: n0r = 8.0e6, n0s = 3.0e6, n0g = 4.0e6
real :: alphar = 0.8, alphas = 0.25, alphag = 0.5
real :: gammar = 17.837789, gammas = 8.2850630, gammag = 11.631769
real :: qmin = 1.0e-12, beta = 1.22
real :: qmin = 1.0e-12, beta = 1.22, qmin1 = 9.e-6

do k = ks, ke
do i = is, ie
Expand Down Expand Up @@ -4718,7 +4718,7 @@ subroutine cloud_diagnosis (is, ie, ks, ke, den, delp, lsm, qmw, qmi, qmr, qms,
! cloud ice (Heymsfield and Mcfarquhar, 1996)
! -----------------------------------------------------------------------

if (qmi (i, k) .gt. qmin) then
if (qmi (i, k) .gt. qmin1) then
qci (i, k) = dpg * qmi (i, k) * 1.0e3
rei_fac = log (1.0e3 * qmi (i, k) * den (i, k))
if (t (i, k) - tice .lt. - 50) then
Expand All @@ -4744,7 +4744,7 @@ subroutine cloud_diagnosis (is, ie, ks, ke, den, delp, lsm, qmw, qmi, qmr, qms,
! cloud ice (Wyser, 1998)
! -----------------------------------------------------------------------

if (qmi (i, k) .gt. qmin) then
if (qmi (i, k) .gt. qmin1) then
qci (i, k) = dpg * qmi (i, k) * 1.0e3
bw = - 2. + 1.e-3 * log10 (den (i, k) * qmi (i, k) / rho_0) * max (0.0, tice - t (i, k)) ** 1.5
rei (i, k) = 377.4 + bw * (203.3 + bw * (37.91 + 2.3696 * bw))
Expand Down Expand Up @@ -4774,7 +4774,7 @@ subroutine cloud_diagnosis (is, ie, ks, ke, den, delp, lsm, qmw, qmi, qmr, qms,
! snow (Lin et al., 1983)
! -----------------------------------------------------------------------

if (qms (i, k) .gt. qmin) then
if (qms (i, k) .gt. qmin1) then
qcs (i, k) = dpg * qms (i, k) * 1.0e3
lambdas = exp (0.25 * log (pi * rhos * n0s / qms (i, k) / den (i, k)))
res (i, k) = 0.5 * exp (log (gammas / 6) / alphas) / lambdas * 1.0e6
Expand Down
53 changes: 27 additions & 26 deletions gfsphysics/physics/satmedmfvdifq.f
View file
Original file line number Diff line number Diff line change
Expand Up @@ -151,7 +151,7 @@ subroutine satmedmfvdifq(ix,im,km,ntrac,ntcw,ntiw,ntke,
& rlmn, rlmn1, rlmx, elmx,
& ttend, utend, vtend, qtend,
& zfac, zfmin, vk, spdk2,
& tkmin, xkzinv, xkgdx,
& tkmin, tkminx, xkzinv, xkgdx,
& zlup, zldn, bsum,
& tem, tem1, tem2,
& ptem, ptem0, ptem1, ptem2
Expand All @@ -176,11 +176,11 @@ subroutine satmedmfvdifq(ix,im,km,ntrac,ntcw,ntiw,ntke,
parameter(prmin=0.25,prmax=4.0)
parameter(pr0=1.0,prtke=1.0,prscu=0.67)
parameter(f0=1.e-4,crbmin=0.15,crbmax=0.35)
parameter(tkmin=1.e-9,dspmax=10.0)
parameter(tkmin=1.e-9,tkminx=0.2,dspmax=10.0)
parameter(qmin=1.e-8,qlmin=1.e-12,zfmin=1.e-8)
parameter(aphi5=5.,aphi16=16.)
parameter(elmfac=1.0,elefac=1.0,cql=100.)
parameter(dw2min=1.e-4,dkmax=1000.,xkgdx=25000.)
parameter(dw2min=1.e-4,dkmax=1000.,xkgdx=5000.)
parameter(qlcr=3.5e-5,zstblmax=2500.,xkzinv=0.1)
parameter(h1=0.33333333)
parameter(ck0=0.4,ck1=0.15,ch0=0.4,ch1=0.15)
Expand Down Expand Up @@ -273,20 +273,20 @@ subroutine satmedmfvdifq(ix,im,km,ntrac,ntcw,ntiw,ntke,
xkzo(i,k) = 0.0
xkzmo(i,k) = 0.0
if (k < kinver(i)) then
! vertical background diffusivity
ptem = prsi(i,k+1) * tx1(i)
tem1 = 1.0 - ptem
tem2 = tem1 * tem1 * 10.0
tem2 = min(1.0, exp(-tem2))
xkzo(i,k) = xkzm_hx(i) * tem2
!
! minimum turbulent mixing length
ptem = prsl(i,k) * tx1(i)
tem1 = 1.0 - ptem
tem2 = tem1 * tem1 * 2.5
tem2 = min(1.0, exp(-tem2))
rlmnz(i,k)= rlmn * tem2
rlmnz(i,k)= max(rlmnz(i,k), rlmn1)
! vertical background diffusivity for momentum
! vertical background diffusivity
ptem = prsi(i,k+1) * tx1(i)
tem1 = 1.0 - ptem
tem2 = tem1 * tem1 * 10.0
tem2 = min(1.0, exp(-tem2))
xkzo(i,k) = xkzm_hx(i) * tem2
! vertical background diffusivity for momentum
if (ptem >= xkzm_s) then
xkzmo(i,k) = xkzm_mx(i)
kx1(i) = k + 1
Expand Down Expand Up @@ -674,20 +674,20 @@ subroutine satmedmfvdifq(ix,im,km,ntrac,ntcw,ntiw,ntke,
!
! background diffusivity decreasing with increasing surface layer stability
!
do i = 1, im
if(.not.sfcflg(i)) then
tem = (1. + 5. * rbsoil(i))**2.
! tem = (1. + 5. * zol(i))**2.
frik(i) = 0.1 + 0.9 / tem
endif
enddo
!
do k = 1,km1
do i=1,im
xkzo(i,k) = frik(i) * xkzo(i,k)
xkzmo(i,k)= frik(i) * xkzmo(i,k)
enddo
enddo
! do i = 1, im
! if(.not.sfcflg(i)) then
! tem = (1. + 5. * rbsoil(i))**2.
!! tem = (1. + 5. * zol(i))**2.
! frik(i) = 0.1 + 0.9 / tem
! endif
! enddo
!
! do k = 1,km1
! do i=1,im
! xkzo(i,k) = frik(i) * xkzo(i,k)
! xkzmo(i,k)= frik(i) * xkzmo(i,k)
! enddo
! enddo
!
! The background vertical diffusivities in the inversion layers are limited
! to be less than or equal to xkzminv
Expand Down Expand Up @@ -867,13 +867,14 @@ subroutine satmedmfvdifq(ix,im,km,ntrac,ntcw,ntiw,ntke,
do i = 1, im
if(k == 1) then
tem = ckz(i,1)
tem1 = xkzmo(i,1)
tem1 = 0.5 * xkzmo(i,1)
else
tem = 0.5 * (ckz(i,k-1) + ckz(i,k))
tem1 = 0.5 * (xkzmo(i,k-1) + xkzmo(i,k))
endif
ptem = tem1 / (tem * elm(i,k))
tkmnz(i,k) = ptem * ptem
tkmnz(i,k) = min(tkmnz(i,k), tkminx)
tkmnz(i,k) = max(tkmnz(i,k), tkmin)
enddo
enddo
Expand Down
9 changes: 2 additions & 7 deletions io/FV3GFS_io.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1119,15 +1119,10 @@ subroutine sfc_prop_restart_read (Sfcprop, Atm_block, Model, fv_domain)
Sfcprop(nb)%zorll(ix) = Sfcprop(nb)%zorlo(ix)
Sfcprop(nb)%zorl(ix) = Sfcprop(nb)%zorlo(ix)
Sfcprop(nb)%tsfc(ix) = Sfcprop(nb)%tsfco(ix)
if (Sfcprop(nb)%slmsk(ix) < 0.1 .or. Sfcprop(nb)%slmsk(ix) > 1.9) then
if (Sfcprop(nb)%slmsk(ix) > 1.9) then
Sfcprop(nb)%landfrac(ix) = 0.0
if (Sfcprop(nb)%oro_uf(ix) > 0.01) then
Sfcprop(nb)%lakefrac(ix) = 1.0 ! lake
else
Sfcprop(nb)%lakefrac(ix) = 0.0 ! ocean
endif
else
Sfcprop(nb)%landfrac(ix) = 1.0 ! land
Sfcprop(nb)%landfrac(ix) = Sfcprop(nb)%slmsk(ix)
endif
enddo
enddo
Expand Down
67 changes: 40 additions & 27 deletions io/module_wrt_grid_comp.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -71,6 +71,7 @@ module module_wrt_grid_comp
logical,save :: first_init=.false.
logical,save :: first_run=.false.
logical,save :: first_getlatlon=.true.
logical,save :: change_wrtidate=.false.
!
!-----------------------------------------------------------------------
!
Expand Down Expand Up @@ -173,7 +174,8 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc)
real(ESMF_KIND_R4) :: valueR4
real(ESMF_KIND_R8) :: valueR8

integer :: attCount, axeslen, jidx, noutfile
integer :: attCount, axeslen, jidx, idx, noutfile
character(19) :: newdate
character(128) :: FBlist_outfilename(100), outfile_name
character(128),dimension(:,:), allocatable :: outfilename
real(8), dimension(:), allocatable :: slat
Expand All @@ -183,7 +185,6 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc)
real(ESMF_KIND_R8) :: geo_lon, geo_lat
real(ESMF_KIND_R8) :: lon1_r8, lat1_r8
real(ESMF_KIND_R8) :: x1, y1, x, y
type(ESMF_Time) :: IO_BASETIME_IAU
type(ESMF_TimeInterval) :: IAU_offsetTI
type(ESMF_DataCopy_Flag) :: copyflag=ESMF_DATACOPY_REFERENCE
! real(8),parameter :: PI=3.14159265358979d0
Expand Down Expand Up @@ -469,8 +470,32 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc)
call ESMF_Finalize(endflag=ESMF_END_ABORT)

endif
!
!-----------------------------------------------------------------------
!*** get write grid component initial time from clock
!-----------------------------------------------------------------------
!
call ESMF_ClockGet(clock =CLOCK & !<-- The ESMF Clock
,startTime=wrt_int_state%IO_BASETIME & !<-- The Clock's starting time
,rc =RC)


call ESMF_TimeGet(time=wrt_int_state%IO_BASETIME,yy=idate(1),mm=idate(2),dd=idate(3),h=idate(4), &
m=idate(5),s=idate(6),rc=rc)
! if (lprnt) write(0,*) 'in wrt initial, io_baseline time=',idate,'rc=',rc
idate(7) = 1
wrt_int_state%idate = idate
wrt_int_state%fdate = idate
! update IO-BASETIME and idate on write grid comp when IAU is enabled
if(iau_offset > 0 ) then
call ESMF_TimeIntervalSet(IAU_offsetTI, h=iau_offset, rc=rc)
wrt_int_state%IO_BASETIME = wrt_int_state%IO_BASETIME + IAU_offsetTI
call ESMF_TimeGet(time=wrt_int_state%IO_BASETIME,yy=idate(1),mm=idate(2),dd=idate(3),h=idate(4), &
m=idate(5),s=idate(6),rc=rc)
wrt_int_state%idate = idate
change_wrtidate = .true.
if (lprnt) print *,'in wrt initial, with iau, io_baseline time=',idate,'rc=',rc
endif
!
! Create field bundle
!-------------------------------------------------------------------
!
Expand Down Expand Up @@ -867,6 +892,17 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc)

if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return

! update the time:units when idate on write grid component is changed
if ( index(trim(attNameList(i)),'time:units')>0) then
if ( change_wrtidate ) then
idx = index(trim(valueS),' since ')
if(lprnt) print *,'in write grid comp, time:unit=',trim(valueS)
write(newdate,'(I4.4,a,I2.2,a,I2.2,a,I2.2,a,I2.2,a,I2.2)') idate(1),'-', &
idate(2),'-',idate(3),' ',idate(4),':',idate(5),':',idate(6)
valueS = valueS(1:idx+6)//newdate
if(lprnt) print *,'in write grid comp, new time:unit=',trim(valueS)
endif
endif
call ESMF_AttributeSet(wrtgrid, convention="NetCDF", purpose="FV3", &
name=trim(attNameList(i)), value=valueS, rc=rc)

Expand Down Expand Up @@ -1036,28 +1072,6 @@ subroutine wrt_initialize(wrt_comp, imp_state_write, exp_state_write, clock, rc)
if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
!
!-----------------------------------------------------------------------
!*** SET THE IO_BaseTime TO THE INITIAL CLOCK TIME.
!-----------------------------------------------------------------------
!
call ESMF_ClockGet(clock =CLOCK & !<-- The ESMF Clock
,startTime=wrt_int_state%IO_BASETIME & !<-- The Clock's starting time
,rc =RC)

call ESMF_TimeGet(time=wrt_int_state%IO_BASETIME,yy=idate(1),mm=idate(2),dd=idate(3),h=idate(4), &
m=idate(5),s=idate(6),rc=rc)
! if (lprnt) write(0,*) 'in wrt initial, io_baseline time=',idate,'rc=',rc
idate(7) = 1
wrt_int_state%idate = idate
wrt_int_state%fdate = idate
if(iau_offset > 0 ) then
call ESMF_TimeIntervalSet(IAU_offsetTI, h=iau_offset, rc=rc)
IO_BASETIME_IAU = wrt_int_state%IO_BASETIME + IAU_offsetTI
call ESMF_TimeGet(time=IO_BASETIME_IAU,yy=idate(1),mm=idate(2),dd=idate(3),h=idate(4), &
m=idate(5),s=idate(6),rc=rc)
! if (lprnt) write(0,*) 'in wrt initial, with iau, io_baseline time=',idate,'rc=',rc
endif
!
!-----------------------------------------------------------------------
!*** SET THE FIRST HISTORY FILE'S TIME INDEX.
!-----------------------------------------------------------------------
!
Expand Down Expand Up @@ -1250,10 +1264,9 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc)
! 'nseconds_num=',nseconds_num,nseconds_den,'mype=',mype
!
nf_seconds = nf_hours*3600+nf_minuteS*60+nseconds+real(nseconds_num)/real(nseconds_den)
! shift forecast hour by iau_offset if iau is on.
nf_seconds = nf_seconds - iau_offset*3600
wrt_int_state%nfhour = nf_seconds/3600.
nf_hours = int(nf_seconds/3600.)
if(mype == lead_write_task) print *,'in write grid comp, nf_hours=',nf_hours
! if iau_offset > nf_hours, don't write out anything
if (nf_hours < 0) return

Expand Down
Loading

0 comments on commit b18739d

Please sign in to comment.