Merge pull request NOAA-EMC#24 from climbfuji/update_dtc_develop_from…

…_develop

dtc/develop: update from develop 2020/01/27

atmos_model.F90

@@ -113,7 +113,8 @@ module atmos_model_mod
                               FV3GFS_diag_register, FV3GFS_diag_output,  &
                               DIAG_SIZE
 use fv_iau_mod,         only: iau_external_data_type,getiauforcing,iau_initialize
-use module_fv3_config,  only: output_1st_tstep_rst, first_kdt, nsout
+use module_fv3_config,  only: output_1st_tstep_rst, first_kdt, nsout,    &
+                              frestart, restart_endfcst
 
 !-----------------------------------------------------------------------
 
@@ -221,7 +222,8 @@ module atmos_model_mod
   logical,parameter :: flip_vc = .true.
 #endif
 
-  real(kind=IPD_kind_phys), parameter :: zero=0.0, one=1.0
+  real(kind=IPD_kind_phys), parameter :: zero = 0.0_IPD_kind_phys, &
+                                         one  = 1.0_IPD_kind_phys
 
 contains
 
@@ -944,17 +946,19 @@ end subroutine update_atmos_model_state
 subroutine atmos_model_end (Atmos)
   type (atmos_data_type), intent(inout) :: Atmos
 !---local variables
-  integer :: idx
+  integer :: idx, seconds
 #ifdef CCPP
   integer :: ierr
 #endif
 
 !-----------------------------------------------------------------------
 !---- termination routine for atmospheric model ----
 
-    call atmosphere_end (Atmos % Time, Atmos%grid)
-    call FV3GFS_restart_write (IPD_Data, IPD_Restart, Atm_block, &
-                               IPD_Control, Atmos%domain)
+    call atmosphere_end (Atmos % Time, Atmos%grid, restart_endfcst)
+    if(restart_endfcst) then
+      call FV3GFS_restart_write (IPD_Data, IPD_Restart, Atm_block, &
+                                 IPD_Control, Atmos%domain)
+    endif
 
 #ifdef CCPP
 !   Fast physics (from dynamics) are finalized in atmosphere_end above;
@@ -1457,6 +1461,24 @@ subroutine update_atmos_chemistry(state, rc)
         enddo
       enddo
 
+      ! -- zero out accumulated fields
+!$OMP parallel do default (none) &
+!$OMP             shared  (nj, ni, Atm_block, IPD_Control, IPD_Data) &
+!$OMP             private (j, jb, i, ib, nb, ix)
+      do j = 1, nj
+        jb = j + Atm_block%jsc - 1
+        do i = 1, ni
+          ib = i + Atm_block%isc - 1
+          nb = Atm_block%blkno(ib,jb)
+          ix = Atm_block%ixp(ib,jb)
+          IPD_Data(nb)%coupling%rainc_cpl(ix)  = zero
+          if (.not.IPD_Control%cplflx) then
+            IPD_Data(nb)%coupling%rain_cpl(ix) = zero
+            IPD_Data(nb)%coupling%snow_cpl(ix) = zero
+          end if
+        enddo
+      enddo
+
       if (IPD_Control%debug) then
         ! -- diagnostics
         write(6,'("update_atmos: prsi   - min/max/avg",3g16.6)') minval(prsi),   maxval(prsi),   sum(prsi)/size(prsi)
@@ -1698,7 +1720,7 @@ subroutine assign_importdata(rc)
                   IPD_Data(nb)%Coupling%ficein_cpl(ix)   = zero
                   if (IPD_Data(nb)%Sfcprop%oceanfrac(ix) > zero) then
                     if (datar8(i,j) >= IPD_control%min_seaice*IPD_Data(nb)%Sfcprop%oceanfrac(ix)) then
-                      IPD_Data(nb)%Coupling%ficein_cpl(ix) = datar8(i,j)
+                      IPD_Data(nb)%Coupling%ficein_cpl(ix) = max(zero, min(datar8(i,j),one))
 !                     if (IPD_Data(nb)%Sfcprop%oceanfrac(ix) == one) IPD_Data(nb)%Sfcprop%slmsk(ix) = 2. !slmsk=2 crashes in gcycle on partial land points
                       IPD_Data(nb)%Sfcprop%slmsk(ix)         = 2.                                        !slmsk=2 crashes in gcycle on partial land points
                       IPD_Data(nb)%Coupling%slimskin_cpl(ix) = 4.

fv3_cap.F90

@@ -30,7 +30,7 @@ module fv3gfs_cap_mod
                                     calendar, calendar_type, cpl,            &
                                     force_date_from_configure,               &
                                     cplprint_flag,output_1st_tstep_rst,      &
-                                    first_kdt                            
+                                    first_kdt,num_restart_interval       
 
   use module_fv3_io_def,      only: num_pes_fcst,write_groups,app_domain,    &
                                     num_files, filename_base,                &
@@ -278,9 +278,16 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)
     CALL ESMF_ConfigLoadFile(config=CF ,filename='model_configure' ,rc=RC)
     if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
 !
-    CALL ESMF_ConfigGetAttribute(config=CF,value=restart_interval, &
-                                 label ='restart_interval:',rc=rc)
+    num_restart_interval = ESMF_ConfigGetLen(config=CF, label ='restart_interval:',rc=rc)
     if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
+    if(mype == 0) print *,'af nems config,num_restart_interval=',num_restart_interval
+    if (num_restart_interval<=0) num_restart_interval = 1
+    allocate(restart_interval(num_restart_interval))
+    restart_interval = 0
+    CALL  ESMF_ConfigGetAttribute(CF,valueList=restart_interval,label='restart_interval:', &
+      count=num_restart_interval, rc=RC)
+    if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
+    if(mype == 0) print *,'af nems config,restart_interval=',restart_interval
 !
     CALL ESMF_ConfigGetAttribute(config=CF,value=calendar, &
                                  label ='calendar:',rc=rc)
@@ -326,9 +333,8 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc)
                                    label ='app_domain:',rc=rc)
       if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
 
-      if(mype == 0) print *,'af nems config,restart_interval=',restart_interval, &
-      'quilting=',quilting,'write_groups=',write_groups,wrttasks_per_group,      &
-      'calendar=',trim(calendar),'calendar_type=',calendar_type
+      if(mype == 0) print *,'af nems config,quilting=',quilting,'write_groups=', &
+        write_groups,wrttasks_per_group,'calendar=',trim(calendar),'calendar_type=',calendar_type
 !
       CALL ESMF_ConfigGetAttribute(config=CF,value=num_files, &
                                    label ='num_files:',rc=rc)

gfsphysics/GFS_layer/GFS_physics_driver.F90

@@ -1103,7 +1103,7 @@ subroutine GFS_physics_driver                                 &
 !*## CCPP ##
       enddo
 !
-!## CCPP ##* note: this block is not yet in CCPP
+!## CCPP ##* GFS_surface_generic.F90/GFS_surface_generic_pre_run
       if (Model%cplflx) then
         do i=1,im
           islmsk_cice(i) = nint(Coupling%slimskin_cpl(i))
@@ -1273,7 +1273,7 @@ subroutine GFS_physics_driver                                 &
           dtdt(i,k)  = zero
           dtdtc(i,k) = zero
 
-!## CCPP ##* note: this block is not yet in CCPP
+!## CCPP ##* GFS_typedefs.F90/interstitial_phys_reset
 !vay-2018
 ! Pure tendency arrays w/o accumulation of Phys-tendencies from each
 !      chain of GFS-physics (later add container for species)
@@ -1911,14 +1911,16 @@ subroutine GFS_physics_driver                                 &
 !    &,' stsoil=',stsoil(ipr,:)
 
 !  --- ...  surface energy balance over seaice
-!## CCPP ##* This section is not in the CCPP yet
+!## CCPP ##* sfc_sice.f/sfc_sice_run (local adjustment to avoid resetting islmsk after call to sfc_sice_run)
         if (Model%cplflx) then
           do i=1,im
             if (flag_cice(i)) then
                islmsk (i) = islmsk_cice(i)
             endif
           enddo
+!*## CCPP ##
 
+!## CCPP ##* sfc_cice.f/sfc_cice_run
 ! call sfc_cice for sea ice points in the coupled model (i.e. islmsk=4)
 !
           call sfc_cice                                                  &
@@ -1954,7 +1956,7 @@ subroutine GFS_physics_driver                                 &
             snowd3(:,2), qss3(:,2), snowmt, gflx3(:,2), cmm3(:,2), chh3(:,2),    &
             evap3(:,2),  hflx3(:,2))
 !*## CCPP ##
-!## CCPP ##* This section is not in the CCPP yet.
+!## CCPP ##* This section is not needed for CCPP.
         if (Model%cplflx) then
           do i = 1, im
             if (flag_cice(i)) then
@@ -2805,7 +2807,7 @@ subroutine GFS_physics_driver                                 &
       endif
 !*## CCPP ##
 
-!## CCPP ##* This block is not yet in CCPP
+!## CCPP ##* GFS_PBL_generic.F90/GFS_PBL_generic_post_run
       if (Model%cplchm) then
         do i = 1, im
           tem1 = max(Diag%q1(i), 1.e-8)
@@ -2814,7 +2816,6 @@ subroutine GFS_physics_driver                                 &
         enddo
         Coupling%dkt     (:,:) = dkt (:,:)
       endif
-!*## CCPP ##
 
 !     if (lprnt) then
 !       write(0,*) ' dusfc1=',dusfc1(ipr),' kdt=',kdt,' lat=',lat
@@ -2827,8 +2828,6 @@ subroutine GFS_physics_driver                                 &
 
 !  --- ...  coupling insertion
 
-!## CCPP ## This block is not in the CCPP yet. It should probably be put in 
-! GFS_PBL_generic.F90/GFS_PBL_generic_post_run.
       if (Model%cplflx) then
         do i=1,im
           if (Sfcprop%oceanfrac(i) > zero) then               ! Ocean only, NO LAKES
@@ -3182,7 +3181,7 @@ subroutine GFS_physics_driver                                 &
       Stateout%gq0(1:im,:,:) = Statein%qgrs(1:im,:,:) + dqdt(1:im,:,:) * dtp
 !*## CCPP ##
 
-! DH* TODO - WHERE IS THIS IN CCPP?
+!## CCPP ##* This is not in the CCPP yet.
 !================================================================================
 !     above: updates of the state by UGWP oro-GWS and RF-damp
 !  Diag%tav_ugwp & Diag%uav_ugwp(i,k)-Updated U-T state before moist/micro !  physics
@@ -3197,7 +3196,7 @@ subroutine GFS_physics_driver                                 &
           enddo
         enddo
       endif
-! *DH
+!*## CCPP ##
 
 !================================================================================
 ! It is not clear Do we need it, "ideaca_up", having stability check inside UGWP-module
@@ -3308,9 +3307,13 @@ subroutine GFS_physics_driver                                 &
         dtdt(1:im,:) = Stateout%gt0(1:im,:)
       endif   ! end if_ldiag3d/cnvgwd
 
-      if (Model%ldiag3d) then
+      if (Model%ldiag3d .or. Model%cplchm) then
         dqdt(1:im,:,1) = Stateout%gq0(1:im,:,1)
-      endif   ! end if_ldiag3d
+      endif   ! end if_ldiag3d/cplchm
+
+      if (Model%cplchm) then
+        Coupling%dqdti(1:im,:) = zero
+      endif   ! end if_cplchm
 !*## CCPP ##
 
 !## CCPP ## Only get_prs_fv3.F90/get_phi_fv3_run is a scheme (GFS_HYDRO is assumed to be undefined)

gfsphysics/GFS_layer/GFS_radiation_driver.F90

@@ -2135,18 +2135,40 @@ subroutine GFS_radiation_driver                             &
               Diag%fluxr(i,11-j) = Diag%fluxr(i,11-j) + tem0d * Statein%prsi(i,itop+kt)
               Diag%fluxr(i,14-j) = Diag%fluxr(i,14-j) + tem0d * Statein%prsi(i,ibtc+kb)
               Diag%fluxr(i,17-j) = Diag%fluxr(i,17-j) + tem0d * Statein%tgrs(i,itop)
+            enddo
+          enddo
 
 !       Anning adds optical depth and emissivity output
-              tem1 = 0.
-              tem2 = 0.
-              do k=ibtc,itop
-                tem1 = tem1 + cldtausw(i,k)      ! approx .55 mu channel
-                tem2 = tem2 + cldtaulw(i,k)      ! approx 10. mu channel
+          if (Model%lsswr .and. (nday > 0)) then
+            do j = 1, 3
+              do i = 1, IM
+                tem0d = raddt * cldsa(i,j)
+                itop  = mtopa(i,j) - kd
+                ibtc  = mbota(i,j) - kd
+                tem1 = 0.
+                do k=ibtc,itop
+                  tem1 = tem1 + cldtausw(i,k)      ! approx .55 um channel
+                enddo
+                Diag%fluxr(i,43-j) = Diag%fluxr(i,43-j) + tem0d * tem1
               enddo
-              Diag%fluxr(i,43-j) = Diag%fluxr(i,43-j) + tem0d * tem1
-              Diag%fluxr(i,46-j) = Diag%fluxr(i,46-j) + tem0d * (1.0-exp(-tem2))
             enddo
-          enddo
+          endif
+
+          if (Model%lslwr) then
+            do j = 1, 3
+              do i = 1, IM
+                tem0d = raddt * cldsa(i,j)
+                itop  = mtopa(i,j) - kd
+                ibtc  = mbota(i,j) - kd
+                tem2 = 0.
+                do k=ibtc,itop
+                  tem2 = tem2 + cldtaulw(i,k)      ! approx 10. um channel
+                enddo
+                Diag%fluxr(i,46-j) = Diag%fluxr(i,46-j) + tem0d * (1.0-exp(-tem2))
+              enddo
+            enddo
+          endif
+
         endif
 
       endif                                ! end_if_lssav

gfsphysics/physics/gfdl_cloud_microphys.F90

@@ -3266,7 +3266,7 @@ subroutine fall_speed (ktop, kbot, den, qs, qi, qg, ql, tk, vts, vti, vtg)
             else
                 tc (k) = tk (k) - tice
                 vti (k) = (3. + log10 (qi (k) * den (k))) * (tc (k) * (aa * tc (k) + bb) + cc) + dd * tc (k) + ee
-                vti (k) = vi0 * exp (log_10 * vti (k)) * 0.8 
+                vti (k) = vi0 * exp (log_10 * vti (k)) * 0.9
                 vti (k) = min (vi_max, max (vf_min, vti (k)))
             endif
         enddo

gfsphysics/physics/moninedmf_hafs.f

@@ -1360,7 +1360,11 @@ subroutine moninedmf_hafs(ix,im,km,ntrac,ntcw,dv,du,tau,rtg,      &
           tem = 0.5 * (diss(i,k-1)+diss(i,k))
           tem  = max(tem, 0.)
           ttend = tem / cp
-          tau(i,k) = tau(i,k) + 0.5*ttend
+          if (alpha .gt. 0.0) then
+            tau(i,k) = tau(i,k) + 0.5*ttend
+          else
+            tau(i,k) = tau(i,k) + 0.7*ttend  ! in HWRF/HMON, use 0.7
+          endif
         enddo
       enddo
 !

gfsphysics/physics/samfdeepcnv.f

@@ -1547,22 +1547,22 @@ subroutine samfdeepcnv(im,ix,km,delt,itc,ntc,ntk,ntr,delp,
       enddo
       enddo
       do i = 1, im
-        betamn = betas
-        if(islimsk(i) == 1) betamn = betal
-        if(ntk > 0) then
-          betamx = betamn + dbeta
-          if(tkemean(i) > tkemx) then
-            beta = betamn
-          else if(tkemean(i) < tkemn) then
-            beta = betamx
+        if(cnvflg(i)) then
+          betamn = betas
+          if(islimsk(i) == 1) betamn = betal
+          if(ntk > 0) then
+            betamx = betamn + dbeta
+            if(tkemean(i) > tkemx) then
+              beta = betamn
+            else if(tkemean(i) < tkemn) then
+              beta = betamx
+            else
+              tem = (betamx - betamn) * (tkemean(i) - tkemn)
+              beta = betamx - tem  / dtke
+            endif
           else
-            tem = (betamx - betamn) * (tkemean(i) - tkemn)
-            beta = betamx - tem  / dtke
+            beta = betamn
           endif
-        else
-          beta = betamn
-        endif
-        if(cnvflg(i)) then
           dz  = (sumx(i)+zi(i,1))/float(kbcon(i))
           tem = 1./float(kbcon(i))
           xlamd(i) = (1.-beta**tem)/dz

gfsphysics/physics/satmedmfvdifq.f

@@ -148,7 +148,8 @@ subroutine satmedmfvdifq(ix,im,km,ntrac,ntcw,ntiw,ntke,
      &                     epsi,    beta,   chx,    cqx,
      &                     rdt,     rdz,    qmin,   qlmin,
      &                     rimin,   rbcr,   rbint,  tdzmin,
-     &                     rlmn,    rlmn1,  rlmx,   elmx,
+     &                     rlmn,    rlmn1,  rlmn2,  
+     &                     rlmx,    elmx,
      &                     ttend,   utend,  vtend,  qtend,
      &                     zfac,    zfmin,  vk,     spdk2,
      &                     tkmin,   tkminx, xkzinv, xkgdx,
@@ -172,7 +173,8 @@ subroutine satmedmfvdifq(ix,im,km,ntrac,ntcw,ntiw,ntke,
       parameter(gamcrt=3.,gamcrq=0.,sfcfrac=0.1)
       parameter(vk=0.4,rimin=-100.)
       parameter(rbcr=0.25,zolcru=-0.02,tdzmin=1.e-3)
-      parameter(rlmn=30.,rlmn1=5.,rlmx=300.,elmx=300.)
+      parameter(rlmn=30.,rlmn1=5.,rlmn2=10.)
+      parameter(rlmx=300.,elmx=300.)
       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)
@@ -698,8 +700,9 @@ subroutine satmedmfvdifq(ix,im,km,ntrac,ntcw,ntiw,ntke,
 !         if(tem1 > 1.e-5) then
           tem1 = tvx(i,k+1)-tvx(i,k)
           if(tem1 > 0.) then
-             xkzo(i,k)  = min(xkzo(i,k),xkzinv)
-             xkzmo(i,k) = min(xkzmo(i,k),xkzinv)
+             xkzo(i,k)  = min(xkzo(i,k), xkzinv)
+             xkzmo(i,k) = min(xkzmo(i,k), xkzinv)
+             rlmnz(i,k) = min(rlmnz(i,k), rlmn2)
           endif
         enddo
       enddo

gfsphysics/physics/sflx.f

@@ -251,6 +251,7 @@ subroutine sflx                                                   &
       runoff2 = 0.0
       runoff3 = 0.0
       snomlt  = 0.0
+      rc      = 0.0
 
 !  --- ...  define local variable ice to achieve:
 !             sea-ice case,          ice =  1