Prate_avg fix (NOAA-EMC#979), UPP call freq and single restart var (NOAA-EMC#973) combo (NOAA-EMC#996)

dpsarmie · web-flow · commit 40187b2172e9 · 2025-08-15T12:21:25.000-04:00
* Add fixes for Prate and UPP Call

* Add changes to unit test due to fundamental change to output_fh

* Reduce restart control variables to one
diff --git a/ccpp/CCPP_driver.F90 b/ccpp/CCPP_driver.F90
@@ -59,6 +59,8 @@ subroutine CCPP_step (step, nblks, ierr)
     ! Local variables
     integer :: nb, nt, ntX
     integer :: ierr2
+    integer :: kdt_iau
+    logical :: iauwindow_center
     ! DH* 20210104 - remove kdt_rad when code to clear diagnostic buckets is removed
     integer :: kdt_rad
 
@@ -167,8 +169,17 @@ subroutine CCPP_step (step, nblks, ierr)
       endif
 
       !--- determine if physics diagnostics buckets need to be cleared
+      iauwindow_center = .false.
+      if (GFS_control%iau_offset > 0) then
+        kdt_iau = nint(GFS_control%iau_offset*3600./GFS_control%dtp)
+        if (GFS_control%kdt-1 == kdt_iau) then
+          iauwindow_center = .true.
+          if( GFS_control%me == 0)print *,'in ccpp step vary, iauwindow_center=',iauwindow_center,&
+            'kdt=',GFS_control%kdt,'dtp=',GFS_control%dtp,'iau_offset=',GFS_control%iau_offset
+        endif
+      endif
       if ((mod(GFS_control%kdt-1,GFS_control%nszero)) == 0) then
-        call GFS_Intdiag%phys_zero(GFS_control)
+        call GFS_Intdiag%phys_zero(GFS_control, iauwindow_center=iauwindow_center)
       endif
 
       !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
diff --git a/fv3/atmos_model.F90 b/fv3/atmos_model.F90
@@ -1026,8 +1026,10 @@ subroutine update_atmos_model_state (Atmos, rc)
     if (ANY(nint(output_fh(:)*3600.0) == seconds) .or. (GFS_control%kdt == first_kdt)) then
       if (mpp_pe() == mpp_root_pe()) write(6,*) "---isec,seconds",isec,seconds
       time_int = real(isec)
-      call InitTimeFromIAUOffset(Atmos, time_int, time_intfull, seconds)
-      if (mpp_pe() == mpp_root_pe()) write(6,*) ' gfs diags time since last bucket empty: ',time_int/3600.,'hrs'
+      time_intfull = real(seconds)
+      call InitTimeFromIAUOffset(Atmos, time_int, time_intfull)
+      if (mpp_pe() == mpp_root_pe()) write(6,*) 'gfs diags time since last bucket empty: ',time_int,' time_intfull=', &
+         time_intfull,' kdt=',GFS_control%kdt
       call atmosphere_nggps_diag(Atmos%Time)
       call fv3atm_diag_output(Atmos%Time, GFS_Diag, Atm_block, GFS_control%nx, GFS_control%ny, &
                             GFS_control%levs, 1, 1, 1.0_GFS_kind_phys, time_int, time_intfull, &
@@ -1069,27 +1071,17 @@ end subroutine update_atmos_model_state
   !> @param[inout] atmos the main atmos model configurations 
   !> @param[inout] time_init model initialization time
   !> @param[inout] time_intfull model time remaining
-  !> @param seconds time since model initialization
   !>
   !> @author Daniel Sarmiento @date May 16, 2025
- subroutine InitTimeFromIAUOffset(Atmos, time_int, time_intfull, seconds)
+ subroutine InitTimeFromIAUOffset(Atmos, time_int, time_intfull)
 
    type (atmos_data_type),   intent(inout)  :: Atmos
    real(kind=GFS_kind_phys), intent(inout)  :: time_int, time_intfull
-   integer,                  intent(inout)  :: seconds
-   integer                                  :: isec_fhzero
 
    if(Atmos%iau_offset > zero) then
      if( time_int - Atmos%iau_offset*3600. > zero ) then
        time_int = time_int - Atmos%iau_offset*3600.
-     else if(seconds == Atmos%iau_offset*3600) then
-       call get_time (Atmos%Time - diag_time_fhzero, isec_fhzero)
-       time_int = real(isec_fhzero)
-       if (mpp_pe() == mpp_root_pe()) write(6,*) "---iseczero",isec_fhzero
      endif
-   endif
-   time_intfull = real(seconds)
-   if(Atmos%iau_offset > zero) then
      if( time_intfull - Atmos%iau_offset*3600. > zero) then
        time_intfull = time_intfull - Atmos%iau_offset*3600.
      endif
diff --git a/fv3/fv3_cap.F90 b/fv3/fv3_cap.F90
@@ -1005,11 +1005,7 @@ subroutine OutputHours_FrequencyInput(nfhmax, output_startfh, outputfh2)
     if( nfhmax>output_startfh) nfh = nint((nfhmax-output_startfh)/outputfh2(1)) + 1
     if( nfh > 0) then
       allocate(output_fh(nfh))
-      if( output_startfh == 0) then
-        output_fh(1) = dt_atmos/3600.
-      else
-        output_fh(1) = output_startfh
-      endif
+      output_fh(1) = output_startfh + dt_atmos/3600.
       do i=2,nfh
         output_fh(i) = (i-1)*outputfh2(1) + output_startfh
         ! Except fh000, which is the first time output, if any other of the
diff --git a/fv3/io/module_wrt_grid_comp.F90 b/fv3/io/module_wrt_grid_comp.F90
@@ -1723,6 +1723,7 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc)
       logical                               :: use_parallel_netcdf
       real, allocatable                     :: output_fh(:)
       logical                               :: is_restart_bundle, restart_written
+      logical                               :: lupp_history, lrestart
       integer                               :: tileCount
       type(ESMF_Info)                       :: fcstInfo, wrtInfo
       character(len=ESMF_MAXSTR)            :: output_grid_name
@@ -1783,7 +1784,30 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc)
       if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
 
       if (nf_hours < 0) return
-
+      !
+      !*** set up output time on write grid comp:
+      !
+      call ESMF_TimeIntervalGet(timeinterval=io_currtimediff, s=fcst_seconds, rc=rc)
+      if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
+      
+      ! fcst_seconds is number of seconds in io_currtimediff, which is time interval between currenttime and io_basetime.
+      ! io_basetime has been adjusted by iau_offset in initialize phase.
+      ! Since output_fh and frestart and NOT adjusted by iau_offset, in order to compare
+      ! them with fcst_seconds, we must also adjust fcst_seconds by iau_offset
+      if (iau_offset > 0) then
+        fcst_seconds = fcst_seconds + iau_offset*3600
+      endif
+      !
+      !--- set up logic variable to run upp/history and restart
+      !
+      lupp_history = .false.
+      lrestart = .false.
+      if ( ANY(nint(output_fh(:)*3600) == fcst_seconds) ) lupp_history = .true.
+      if ( ANY(frestart(:) == fcst_seconds) ) lrestart = .true.
+      if ( .not. (lupp_history .or. lrestart ) ) return
+      !
+      !--- set up current forecast hours
+      !
       if (lflname_fulltime) then
         ndig = max(log10(nf_hours+0.5)+1., 3.)
         write(cform, '("(I",I1,".",I1,",A1,I2.2,A1,I2.2)")') ndig, ndig
@@ -1794,7 +1818,9 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc)
         write(cfhour, cform) nf_hours
       endif
 !
-       if(lprnt) print *,'in wrt run, nfhour=',nfhour,' cfhour=',trim(cfhour)
+      if(lprnt) print *,'in wrt run, cdate=',cdate(1:4),'fcst_seconds=',fcst_seconds/3600.,'nfhour=',nfhour,&
+      'lupp_history=', lupp_history, 'lrestart=',lrestart,'write grid comp not return,cfhour=',trim(cfhour)
+!
 !
 !-----------------------------------------------------------------------
 !*** loop on the "output_" FieldBundles, i.e. files that need to write out
@@ -2020,7 +2046,7 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc)
 !*** do post
 !-----------------------------------------------------------------------
       lmask_fields = .false.
-      if( wrt_int_state%write_dopost ) then
+      if( wrt_int_state%write_dopost .and. lupp_history ) then
 #ifdef INLINE_POST
         wbeg = MPI_Wtime()
         do n=1,ngrids
@@ -2081,20 +2107,7 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc)
 ! ** now loop through output field bundle
 !-----------------------------------------------------------------------
 
-      call ESMF_TimeIntervalGet(timeinterval=io_currtimediff, s=fcst_seconds, rc=rc)
-      if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
-
-      ! fcst_seconds is number of seconds in io_currtimediff, which is time interval between currenttime and io_basetime.
-      ! io_basetime has been adjusted by iau_offset in initialize phase.
-      ! Since output_fh and frestart and NOT adjusted by iau_offset, in order to compare
-      ! them with fcst_seconds, we must also adjust fcst_seconds by iau_offset
-      if (iau_offset > 0) then
-        fcst_seconds = fcst_seconds + iau_offset*3600
-      endif
-
-      if ( (wrt_int_state%output_history .and. ANY(nint(output_fh(:)*3600.0) == fcst_seconds)) .or. ANY(frestart(:) == fcst_seconds) ) then
-
-        ! if (lprnt) write(0,*)'wrt_run: loop over wrt_int_state%FBCount ',wrt_int_state%FBCount, ' nfhour ',  nfhour, ' cdate ', cdate(1:6)
+      if ( (wrt_int_state%output_history .and. lupp_history) .or. lrestart ) then
         two_phase_loop: do out_phase = 1, 2
 
           restart_written = .false.
@@ -2106,9 +2119,9 @@ subroutine wrt_run(wrt_comp, imp_state_write, exp_state_write,clock,rc)
             is_restart_bundle = .false.
             if (wrtFBName(1:8) == 'restart_') then
               is_restart_bundle = .true.
-              if (.not.(ANY(frestart(:) == fcst_seconds))) cycle
+              if (.not.lrestart) cycle
             else
-              if (.not.(wrt_int_state%output_history .and. ANY(nint(output_fh(:)*3600.0) == fcst_seconds))) cycle
+              if (.not.(wrt_int_state%output_history .and. lupp_history)) cycle
             endif
 
             if (out_phase == 1 .and. is_restart_bundle) cycle
diff --git a/fv3/module_fcst_grid_comp.F90 b/fv3/module_fcst_grid_comp.F90
@@ -595,8 +595,8 @@ subroutine fcst_initialize(fcst_comp, importState, exportState, clock, rc)
     logical               :: top_parent_is_global
     logical               :: history_file_on_native_grid
 
-    integer                       :: num_restart_interval, restart_starttime
-    real,dimension(:),allocatable :: restart_interval
+    integer                       :: num_restart_fh, restart_starttime
+    real,dimension(:),allocatable :: restart_fh
 
     integer           :: urc
     type(ESMF_State)  :: tempState
@@ -634,16 +634,16 @@ subroutine fcst_initialize(fcst_comp, 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
 
-    num_restart_interval = ESMF_ConfigGetLen(config=CF, label ='restart_interval:',rc=rc)
+    num_restart_fh = 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 ufs 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 (mype == 0) print *,'af ufs config,num_restart_fh=',num_restart_fh
+    if (num_restart_fh<=0) num_restart_fh = 1
+    allocate(restart_fh(num_restart_fh))
+    restart_fh = 0
+    call ESMF_ConfigGetAttribute(CF,valueList=restart_fh,label='restart_interval:', &
+                                 count=num_restart_fh, rc=rc)
     if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return
-    if (mype == 0) print *,'af ufs config,restart_interval=',restart_interval
+    if (mype == 0) print *,'af ufs config,restart_fh=',restart_fh
 !
     call fms_init(fcst_mpi_comm%mpi_val)
     call mpp_init()
@@ -762,8 +762,8 @@ subroutine fcst_initialize(fcst_comp, importState, exportState, clock, rc)
     if (mype == 0) write(*,*)'time_init=', date_init,'time=',date,'time_end=',date_end,'dt_atmos=',dt_atmos
     
     call fcst_time_array_setup(Time_init, Time_end, Time_step_restart, &
-                                   Time_restart, num_restart_interval, &
-                                   restart_interval)
+                                   Time_restart, num_restart_fh, &
+                                   restart_fh)
 
 !------ initialize component models ------
 
@@ -1158,51 +1158,41 @@ end subroutine fcst_initialize
   !>
   !> @param[inout] Time_init model initialization time
   !> @param[inout] Time_end model end time
-  !> @param[inout] Time_step_restart restart time based on restart_interval
+  !> @param[inout] Time_step_restart restart time based on restart_fh
   !> @param[inout] Time_restart calculated restart time
-  !> @param[inout] num_restart_interval user defined restart interval
-  !> @param[inout] restart_interval restart interval, allocatable
+  !> @param[inout] num_restart_fh user defined restart interval
+  !> @param[inout] restart_fh restart interval, allocatable
   !>
   !> @author Daniel Sarmiento @date May 16, 2025
   subroutine fcst_time_array_setup(Time_init, Time_end, Time_step_restart, &
-                                   Time_restart, num_restart_interval, &
-                                   restart_interval)
+                                   Time_restart, num_restart_fh, &
+                                   restart_fh)
 
     type(time_type), intent(inout)                 :: Time_init, Time_end, &
                                                       Time_step_restart, &
                                                       Time_restart
     type(time_type)                                :: iautime
-    integer,         intent(inout)                 :: num_restart_interval
+    integer,         intent(inout)                 :: num_restart_fh
     integer                                        :: total_inttime, tmpvar, &
                                                       i, restart_starttime
     logical                                        :: freq_restart
-    real, dimension(:), allocatable, intent(inout) :: restart_interval
+    real, dimension(:), allocatable, intent(inout) :: restart_fh
                                 
     ! set up forecast time array that controls when to write out restart files
     frestart = 0
     call get_time(Time_end - Time_init, total_inttime)
-    ! set iau offset time
-    Atmos%iau_offset    = iau_offset
-    if(iau_offset > 0 ) then
-      iautime =  set_time(iau_offset * 3600, 0)
-    endif
     ! if the second item is -1, the first number is frequency
     freq_restart = .false.
-    if(num_restart_interval == 2) then
-      if(restart_interval(2)== -1) freq_restart = .true.
+    if(num_restart_fh == 2) then
+      if(restart_fh(2)== -1) freq_restart = .true.
     endif
     if(freq_restart) then
-      if(restart_interval(1) >= 0) then
-        tmpvar = restart_interval(1) * 3600
+      if(restart_fh(1) >= 0) then
+        tmpvar = restart_fh(1) * 3600
         Time_step_restart = set_time (tmpvar, 0)
-        if(iau_offset > 0 ) then
-          Time_restart = Time_init + iautime + Time_step_restart
-          frestart(1) = tmpvar + iau_offset *3600
-        else
-          Time_restart = Time_init + Time_step_restart
-          frestart(1) = tmpvar
-        endif
-        if(restart_interval(1) > 0) then
+        Time_restart = Time_init + Time_step_restart
+        frestart(1) = tmpvar
+        if(restart_fh(1) > 0) then
           i = 2
           do while ( Time_restart < Time_end )
             frestart(i) = frestart(i-1) + tmpvar
@@ -1212,17 +1202,13 @@ subroutine fcst_time_array_setup(Time_init, Time_end, Time_step_restart, &
         endif
       endif
     ! otherwise it is an array with forecast time at which the restart files will be written out
-    else if(num_restart_interval >= 1) then
-      if(num_restart_interval == 1 .and. restart_interval(1) == 0 ) then
+    else if(num_restart_fh >= 1) then
+      if(num_restart_fh == 1 .and. restart_fh(1) == 0 ) then
         frestart(1) = total_inttime
       else
-        if(iau_offset > 0 ) then
-          restart_starttime = iau_offset *3600
-        else
-          restart_starttime = 0
-        endif
-        do i=1,num_restart_interval
-          frestart(i) = restart_interval(i) * 3600. + restart_starttime
+        restart_starttime = 0
+        do i=1,num_restart_fh
+          frestart(i) = restart_fh(i) * 3600. + restart_starttime
         enddo
       endif
     endif
diff --git a/tests/test_fv3_cap.F90 b/tests/test_fv3_cap.F90
@@ -3,6 +3,11 @@ program test_output_hours
   use module_fv3_config, only: dt_atmos, output_fh
   use module_fv3_io_def, only: lflname_fulltime
 
+  ! Test changelog:
+  ! 2025/07/29 - D. Sarmiento: output_fh changed. When first hour is not 0 
+  !              (i.e. IAU), then output will now be created at the first
+  !              timestep. This is now aligns with the coldstart functionality.
+
   implicit none
   
   ! Test variables
@@ -73,8 +78,8 @@ subroutine test_frequency_input()
       stop 4
     end if
     
-    ! First value (should be 6.0)
-    if (abs(output_fh(1) - 6.0) > 1e-6) then
+    ! First value (should be 6.5)
+    if (abs(output_fh(1) - (output_startfh + dt_atmos/3600.)) > 1e-6) then
       print *, "First output time is incorrect"
       stop 5
     end if
