Very IMPORTANT: Corrected initialization of input multi-file options

[arango]

I tracked a very nasty bug when using the multifile options in iterative algorithms (4D-Var and friends). The multifile NetCDF option is used for open boundaries, climatology, and forcing data. For example, we can have in ocean.in:

    BRYNAME == ../om/doppio_bry_MercatorV3_2014_atmpress_NAVD88.nc |
               ../om/doppio_bry_MercatorV3_2015_atmpress_NAVD88.nc |
               ../om/doppio_bry_MercatorV3_2016_atmpress_NAVD88.nc

    CLMNAME == ../om/doppio_clm_MercatorV3_2014_atmpress_NAVD88.nc |
               ../om/doppio_clm_MercatorV3_2015_atmpress_NAVD88.nc |
               ../om/doppio_clm_MercatorV3_2016_atmpress_NAVD88.nc

    ...

    NFFILES == 8                          ! number of unique forcing files

    FRCNAME == ../om/lwrad_down_nam_3hourly_MAB_and_GoM_2014.nc |
               ../om/lwrad_down_nam_3hourly_MAB_and_GoM_2015.nc |
               ../om/lwrad_down_nam_3hourly_MAB_and_GoM_2016.nc \

               ../om/Pair_nam_3hourly_MAB_and_GoM_2014.nc |
               ../om/Pair_nam_3hourly_MAB_and_GoM_2015.nc |
               ../om/Pair_nam_3hourly_MAB_and_GoM_2016.nc \

               ../om/Qair_nam_3hourly_MAB_and_GoM_2014.nc |
               ../om/Qair_nam_3hourly_MAB_and_GoM_2015.nc |
               ../om/Qair_nam_3hourly_MAB_and_GoM_2016.nc \

               ../om/rain_nam_3hourly_MAB_and_GoM_2014.nc |
               ../om/rain_nam_3hourly_MAB_and_GoM_2015.nc |
               ../om/rain_nam_3hourly_MAB_and_GoM_2016.nc \

               ../om/swrad_daily_nam_3hourly_MAB_and_GoM_2014.nc |
               ../om/swrad_daily_nam_3hourly_MAB_and_GoM_2015.nc |
               ../om/swrad_daily_nam_3hourly_MAB_and_GoM_2016.nc \

               ../om/Tair_nam_3hourly_MAB_and_GoM_2014.nc |
               ../om/Tair_nam_3hourly_MAB_and_GoM_2015.nc |
               ../om/Tair_nam_3hourly_MAB_and_GoM_2016.nc \

               ../om/Uwind_nam_3hourly_MAB_and_GoM_2014.nc |
               ../om/Uwind_nam_3hourly_MAB_and_GoM_2015.nc |
               ../om/Uwind_nam_3hourly_MAB_and_GoM_2016.nc \

               ../om/Vwind_nam_3hourly_MAB_and_GoM_2014.nc |
               ../om/Vwind_nam_3hourly_MAB_and_GoM_2015.nc |
               ../om/Vwind_nam_3hourly_MAB_and_GoM_2016.nc

Here, the input data is split into annual NetCDF files.

The problem was the location of close_inp in initial.F, tl_initial.F, rp_initial.F, and ad_initial.F. It needs to be before the call to check_multifile and not after because it erased all the initialization in BRY(:), CLM(:), and FRC(:,:) structure that is essential for processing that appropriate multifile. For example, in tl_initia.F we need to have instead:

!
!-----------------------------------------------------------------------
!  If applicable, close all input boundary, climatology, and forcing
!  NetCDF files and set associated parameters to the closed state. This
!  step is essential in iterative algorithms that run the full TLM
!  repetitively. Then, Initialize several parameters in their file
!  structure, so the appropriate input single or multi-file is selected
!  during initialization/restart.
!-----------------------------------------------------------------------
!
!$OMP MASTER
      CALL close_inp (ng, iTLM)
      CALL check_multifile (ng, iTLM)
!$OMP END MASTER
# ifdef DISTRIBUTE
      CALL mp_bcasti (ng, iTLM, exit_flag)
# endif
!$OMP BARRIER
      IF (FoundError(exit_flag, NoError, __LINE__,                      &
     &               __FILE__)) RETURN

The call to close_inp is necessary for iterative algorithms that call the nonlinear, tangent, representer, and adjoint models repetitively during a simulation because the input data needs to be processed from start (end) records during forward (backward) time-stepping.

Also, I corrected an out bounds error in routine close_inp of file close_io.F when closing the boundary, climatology, and forcing NetCDF files. For example, we needed to have:

         FRC(i,ng)%Fcount=Fcount
         FRC(i,ng)%name=TRIM(FRC(i,ng)%files(Fcount))
         lstr=LEN_TRIM(FRC(i,ng)%name)
         FRC(i,ng)%base=FRC(i,ng)%name(1:lstr-3)
         CALL netcdf_close (ng, model, FRC(i,ng)%ncid,                 &
    &                       FRC(i,ng)%name, .FALSE.)
         IF (FoundError(exit_flag, NoError, __LINE__,                  &
    &                   __FILE__)) RETURN

         ...

         BRY(ng)%Fcount=Fcount
         BRY(ng)%name=TRIM(BRY(ng)%files(Fcount))
         lstr=LEN_TRIM(BRY(ng)%name)
         BRY(ng)%base=BRY(ng)%name(1:lstr-3)
         CALL netcdf_close (ng, model, BRY(ng)%ncid,                   &
    &                       BRY(ng)%name, .FALSE.)
         IF (FoundError(exit_flag, NoError, __LINE__,                  &
    &                   __FILE__)) RETURN

         ...

         CLM(ng)%Fcount=Fcount
         CLM(ng)%name=TRIM(CLM(ng)%files(Fcount))
         lstr=LEN_TRIM(CLM(ng)%name)
         CLM(ng)%base=CLM(ng)%name(1:lstr-3)
         CALL netcdf_close (ng, model, CLM(ng)%ncid,                   &
    &                       CLM(ng)%name, .FALSE.)
         IF (FoundError(exit_flag, NoError, __LINE__,                  &
    &                   __FILE__)) RETURN

We need to pass FRC(i,ng)%name, BRY(ng)%name, and CLM(ng)%name arguments to routine netcdf_close.

Many thanks to Andy Moore for bringing this to my attention and helping me to set-up a case that I can follow in the debugger.