initial.F

Go to the documentation of this file.

#include "cppdefs.h"
      SUBROUTINE initial
!
!git $Id$
!================================================== Hernan G. Arango ===
!  Copyright (c) 2002-2025 The ROMS Group                              !
!    Licensed under a MIT/X style license                              !
!    See License_ROMS.md                                               !
!=======================================================================
!                                                                      !
!  This routine initializes all model variables.                       !
!                                                                      !
!=======================================================================
!
      USE mod_param
      USE mod_parallel
#ifdef BBL_MODEL
      USE mod_bbl
#endif
#ifdef FOUR_DVAR
      USE mod_fourdvar
#endif
      USE mod_grid
      USE mod_iounits
      USE mod_ncparam
#ifdef NESTING
      USE mod_nesting
#endif
      USE mod_ocean
      USE mod_scalars
      USE mod_stepping
!
      USE analytical_mod
      USE close_io_mod,      ONLY : close_inp
      USE dateclock_mod,     ONLY : time_string
#if defined FOUR_DVAR && !defined CORRELATION
      USE def_ini_mod,       ONLY : def_ini
#endif
#ifdef DISTRIBUTE
      USE distribute_mod,    ONLY : mp_bcasti
#endif
      USE get_state_mod,     ONLY : get_state
#ifdef WET_DRY
      USE get_wetdry_mod,    ONLY : get_wetdry
#endif
#ifdef TLM_CHECK
      USE ini_adjust_mod,    ONLY : ini_perturb
#endif
#if defined SOLVE3D && defined NONLINEAR
      USE ini_fields_mod,    ONLY : set_zeta_timeavg
#endif
      USE ini_hmixcoef_mod,  ONLY : ini_hmixcoef
#if defined WAV_COUPLING && defined MCT_LIB
      USE mct_coupler_mod,   ONLY : ocn2wav_coupling
#endif
#ifdef NESTING
      USE nesting_mod,       ONLY : nesting
#endif
#ifdef SOLVE3D
      USE set_depth_mod,     ONLY : set_depth0, set_depth
      USE omega_mod,         ONLY : omega
      USE rho_eos_mod,       ONLY : rho_eos
      USE set_massflux_mod,  ONLY : set_massflux
#endif
#if defined OBSERVATIONS && !defined RBL4DVAR_FCT_SENSITIVITY
      USE obs_initial_mod,   ONLY : obs_initial
#endif
#ifdef MASKING
      USE set_masks_mod,     ONLY : set_masks
#endif
      USE stiffness_mod,     ONLY : stiffness
      USE strings_mod,       ONLY : founderror
#ifdef WET_DRY
      USE wetdry_mod,        ONLY : wetdry
#endif
#if defined PROPAGATOR || \
    (defined masking    && (defined read_water || defined write_water))
      USE wpoints_mod,       ONLY : wpoints
#endif
!
      implicit none
!
!  Local variable declarations.
!
      logical :: update = .false.
!
      integer :: Fcount
      integer :: ng, thread, tile
#ifdef NESTING
      integer :: ig, nl
      integer :: cr, i, m
#endif
 
      integer, dimension(Ngrids) :: IniRec, Tindex
 
#if defined ADJUST_BOUNDARY || \
    defined adjust_stflux   || defined adjust_wstress
      integer :: irec
#endif
!
#ifdef SP4DVAR
      real(dp) :: my_dstart
!
#endif
!
      character (len=*), parameter :: MyFile =                          &
     &  __FILE__
!
!=======================================================================
!   Initialize model variables.
!=======================================================================
!
!$OMP MASTER
      IF (master) THEN
#if defined PERTURBATION
        WRITE (stdout,10) nrun
 10     FORMAT (/,' <<<< Ensemble/Perturbation Run: ',i5.5,' >>>>',/)
#elif defined I4DVAR      || defined RBL4DVAR          || \
      defined r4dvar      || defined sensitivity_4dvar || \
      defined tl_rbl4dvar || defined tl_r4dvar
        WRITE (stdout,10) outer, inner
 10     FORMAT (/,' <<<< 4D Variational Data Assimilation, ',           &
     &          'Outer = ',i3.3, ', Inner = ',i3.3,' >>>>',/)
#endif
        WRITE (stdout,20) 'INITIAL: Configuring and initializing ',     &
     &                    'forward nonlinear model ...'
 20     FORMAT (/,1x,a,a,/,1x,'*******')
      END IF
!$OMP END MASTER
!
!-----------------------------------------------------------------------
!  Initialize time stepping indices and counters.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        iif(ng)=1
        indx1(ng)=1
        next_kstp(ng)=1
        kstp(ng)=1
        krhs(ng)=1
        knew(ng)=1
        predictor_2d_step(ng)=.false.
!
        iic(ng)=0
        nstp(ng)=1
        nrhs(ng)=1
        nnew(ng)=1
#ifdef FLOATS
        nf(ng)=0
        nfp1(ng)=1
        nfm1(ng)=4
        nfm2(ng)=3
        nfm3(ng)=2
#endif
!
        inirec(ng)=nrrec(ng)
        tindex(ng)=1
!
        synchro_flag(ng)=.true.
        first_time(ng)=0
#ifdef SP4DVAR
        IF (initime(ng).lt.0.0_dp) THEN
          my_dstart=dstart
        ELSE
          my_dstart=initime(ng)/86400.0_dp
        END IF
#else
        tdays(ng)=dstart
#endif
        time(ng)=tdays(ng)*day2sec
!$OMP MASTER
        ntstart(ng)=int((time(ng)-dstart*day2sec)/dt(ng))+1
        ntend(ng)=ntstart(ng)+ntimes(ng)-1
        ntfirst(ng)=ntstart(ng)
!$OMP END MASTER
!$OMP BARRIER
        step_counter(ng)=0
      END DO
!
!  Initialize global diagnostics variables.
!
      avgke=0.0_dp
      avgpe=0.0_dp
      avgkp=0.0_dp
      volume=0.0_dp
 
#ifdef PROFILE
!
!-----------------------------------------------------------------------
!  Start time wall clocks.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        DO thread=thread_range
          CALL wclock_on (ng, inlm, 2, __line__, myfile)
        END DO
      END DO
!$OMP BARRIER
#endif
 
#ifdef FOUR_DVAR
!
!-----------------------------------------------------------------------
!  If variational data assimilation, reset several IO switches and
!  variables.
!-----------------------------------------------------------------------
!
!  Set initial conditions record to process. If applicable open existing
!  nonlinear model initial conditions NetCDF file and, if needed, define
!  new variables. Then, inquire about available variables.
!
      DO ng=1,ngrids
        IF (ldefini(ng)) THEN
          ldefini(ng)=.false.          ! needed to inquire variables IDs
# ifndef CORRELATION
          CALL def_ini (ng)
# endif
# ifdef DISTRIBUTE
          CALL mp_bcasti (ng, inlm, exit_flag)
# endif
          IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
          inirec(ng)=nrrec(ng)
          ini(ng)%Rindex=inirec(ng)
        ELSE
          inirec(ng)=ini(ng)%Rindex
        END IF
      END DO
 
# ifdef ADJUST_BOUNDARY
!
!  Initialize open boundary counter for storage arrays.
!
      DO ng=1,ngrids
        obccount(ng)=0
      END DO
# endif
# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
!
!  Initialize surface forcing counter for storage arrays.
!
      DO ng=1,ngrids
        sfcount(ng)=0
      END DO
# endif
!
!  Reset nonlinear history time record counters. These counters are
!  reset on every iteration pass. This file is created on the first
!  iteration pass.
!
      DO ng=1,ngrids
        his(ng)%Rindex=0
        fcount=his(ng)%Fcount
        his(ng)%Nrec(fcount)=0
      END DO
 
# ifdef I4DVAR
!
!  Activate switches to writting data into average, history and
!  restart files.
!
      DO ng=1,ngrids
        lwrtavg(ng)=.true.
        lwrthis(ng)=.true.
        lwrtrst(ng)=.true.
      END DO
# endif
!$OMP BARRIER
#endif
!
!-----------------------------------------------------------------------
!  Set application grid, metrics, and associated variables and
!  parameters.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        IF (setgridconfig(ng)) THEN
          CALL set_grid (ng, inlm)
          setgridconfig(ng)=.false.
        END IF
      END DO
!
!-----------------------------------------------------------------------
!  Initialize horizontal mixing coefficients. If applicable, scale
!  mixing coefficients according to the grid size (smallest area).
#ifndef ANA_SPONGE
!  Also increase their values in sponge areas using the "visc_factor"
!  and/or "diff_factor" read from input Grid NetCDF file.
#endif
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        DO tile=first_tile(ng),last_tile(ng),+1
          CALL ini_hmixcoef (ng, tile, inlm)
        END DO
!$OMP BARRIER
      END DO
 
#ifdef ANA_SPONGE
!
!-----------------------------------------------------------------------
!  Increase horizontal mixing coefficients in sponge areas using
!  analytical functions.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        IF (lsponge(ng)) THEN
          DO tile=first_tile(ng),last_tile(ng),+1
            CALL ana_sponge (ng, tile, inlm)
          END DO
!$OMP BARRIER
        END IF
      END DO
#endif
!
!=======================================================================
!  Initialize model state variables and forcing.  This part is
!  executed for each ensemble/perturbation/iteration run.
!=======================================================================
 
#ifdef TLM_CHECK
!
!  Clear state variables.
!
      DO ng=1,ngrids
        DO tile=first_tile(ng),last_tile(ng),+1
          CALL initialize_ocean (ng, tile, inlm)
        END DO
!$OMP BARRIER
      END DO
#endif
 
#if defined SOLVE3D && !defined INI_FILE
!
!-----------------------------------------------------------------------
!  If analytical initial conditions, compute initial time-evolving
!  depths with zero free-surface.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        DO tile=first_tile(ng),last_tile(ng),+1
          CALL set_depth (ng, tile, inlm)
        END DO
!$OMP BARRIER
      END DO
#endif
!
!-----------------------------------------------------------------------
!  Set primitive variables initial conditions.
!-----------------------------------------------------------------------
 
#ifdef ANA_INITIAL
!
!  Analytical initial conditions for momentum and active tracers.
!
      DO ng=1,ngrids
        IF (nrrec(ng).eq.0) THEN
          DO tile=first_tile(ng),last_tile(ng),+1
            CALL ana_initial (ng, tile, inlm)
          END DO
!$OMP BARRIER
        END IF
      END DO
#endif
 
#if defined ANA_PASSIVE && defined SOLVE3D
!
!  Analytical initial conditions for inert passive tracers.
!
      DO ng=1,ngrids
        IF (nrrec(ng).eq.0) THEN
          DO tile=first_tile(ng),last_tile(ng),+1
            CALL ana_passive (ng, tile, inlm)
          END DO
!$OMP BARRIER
        END IF
      END DO
#endif
 
#if defined ANA_BIOLOGY && defined SOLVE3D
!
!  Analytical initial conditions for biology tracers.
!
      DO ng=1,ngrids
        IF (nrrec(ng).eq.0) THEN
          DO tile=first_tile(ng),last_tile(ng),+1
            CALL ana_biology (ng, tile, inlm)
          END DO
!$OMP BARRIER
        END IF
      END DO
#endif
 
#if defined ANA_SEDIMENT && defined SOLVE3D
!
!  Analytical initial conditions for sediment tracers.
!
 
      DO ng=1,ngrids
        IF (nrrec(ng).eq.0) THEN
          DO tile=first_tile(ng),last_tile(ng),+1
            CALL ana_sediment (ng, tile, inlm)
          END DO
!$OMP BARRIER
        END IF
      END DO
#endif
 
#if defined INI_FILE && !defined RBL4DVAR_FCT_SENSITIVITY
!
!  Read in initial conditions from initial NetCDF file.
!
      DO ng=1,ngrids
!$OMP MASTER
        CALL get_state (ng, inlm, 1, ini(ng), inirec(ng), tindex(ng))
!$OMP END MASTER
# ifdef DISTRIBUTE
        CALL mp_bcasti (ng, inlm, exit_flag)
# endif
!$OMP BARRIER
        IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
        time(ng)=io_time                     ! needed for shared-memory
      END DO
#else
!
!  If restart, read in initial conditions restart NetCDF file.
!
      DO ng=1,ngrids
        IF (nrrec(ng).ne.0) THEN
!$OMP MASTER
          CALL get_state (ng, 0, 1, ini(ng), inirec(ng), tindex(ng))
!$OMP END MASTER
# ifdef DISTRIBUTE
          CALL mp_bcasti (ng, inlm, exit_flag)
# endif
!$OMP BARRIER
          IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
          time(ng)=io_time                   ! needed for shared-memory
        END IF
      END DO
#endif
 
#ifdef WET_DRY
!
!-----------------------------------------------------------------------
!  Process initial wet/dry masks.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
!
!  If restart, read in wet/dry masks.
!
        IF (nrrec(ng).ne.0) THEN
!$OMP MASTER
# ifdef DISTRIBUTE
          CALL get_wetdry (ng, myrank, inlm, inirec(ng))
# else
          CALL get_wetdry (ng, -1, inlm, inirec(ng))
# endif
!$OMP END MASTER
# ifdef DISTRIBUTE
          CALL mp_bcasti (ng, inlm, exit_flag)
# endif
!$OMP BARRIER
          IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
        ELSE
          DO tile=first_tile(ng),last_tile(ng),+1
            CALL wetdry (ng, tile, tindex(ng), .true.)
          END DO
!$OMP BARRIER
        END IF
      END DO
#endif
 
#if defined OBSERVATIONS && !defined RBL4DVAR_FCT_SENSITIVITY
!
!-----------------------------------------------------------------------
!  Open observations NetCDF file and initialize various variables
!  needed for processing the nonlinear state solution at observation
!  locations. Need to be done after processing initial conditions since
!  the correct initial time is needed to determine the first "ObsTime"
!  to process.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
!$OMP MASTER
        CALL obs_initial (ng, inlm, .false.)
!$OMP END MASTER
        IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
      END DO
!$OMP BARRIER
#endif
 
#if (defined ADJUST_BOUNDARY          || \
     defined adjust_stflux            || \
     defined adjust_wstress)          && \
    (defined rbl4dvar                 || \
     defined rbl4dvar_ana_sensitivity || \
     defined rbl4dvar_fct_sensitivity || \
     defined tl_rbl4dvar)
!
!-----------------------------------------------------------------------
!  Read in the surface forcing and or open boundary conditions
!  increments for RBL4D-Var from record IniRec of the NLM initial
!  NetCDF file.
!-----------------------------------------------------------------------
!
      IF (nrun.gt.1) THEN
        DO ng=1,ngrids
!$OMP MASTER
          CALL get_state (ng, 5, 5, ini(ng), inirec(ng), tindex(ng))
!$OMP END MASTER
# ifdef DISTRIBUTE
          CALL mp_bcasti (ng, inlm, exit_flag)
# endif
!$OMP BARRIER
          IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
        END DO
      END IF
#endif
 
#ifdef TLM_CHECK
!
!-----------------------------------------------------------------------
!  Add a perturbation to nonlinear state variable according to the outer
!  loop iteration with the steepest descent direction of the gradient
!  (adjoint state).
!-----------------------------------------------------------------------
!
      IF (outer.ge.1) THEN
        DO ng=1,ngrids
          DO tile=first_tile(ng),last_tile(ng),+1
            CALL ini_perturb (ng, tile, lnew(ng), tindex(ng))
          END DO
!$OMP BARRIER
        END DO
      END IF
#endif
 
#ifdef SOLVE3D
!
!-----------------------------------------------------------------------
!  Compute time independent (Zt_avg1=0). Set Zt_avg1 to the initial
!  free-surface value (no lateral boundary conditions applied). Then,
!  set initial time dependent depths and level thicknesses.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        DO tile=first_tile(ng),last_tile(ng),+1
          CALL set_depth0 (ng, tile, inlm)
          CALL set_zeta_timeavg (ng, tile, inlm)
          CALL set_depth  (ng, tile, inlm)
        END DO
!$OMP BARRIER
      END DO
!
!-----------------------------------------------------------------------
!  Compute initial horizontal mass fluxes, Hz*u/n and Hz*v/m.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        DO tile=first_tile(ng),last_tile(ng),+1
          CALL set_massflux (ng, tile, inlm)
        END DO
!$OMP BARRIER
      END DO
!
!-----------------------------------------------------------------------
!  Compute initial S-coordinates vertical velocity. Compute initial
!  density anomaly from potential temperature and salinity via equation
!  of state for seawater.  Also compute other equation of state related
!  quatities.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        DO tile=first_tile(ng),last_tile(ng),+1
          CALL omega (ng, tile, inlm)
          CALL rho_eos (ng, tile, inlm)
        END DO
!$OMP BARRIER
      END DO
#endif
 
#ifdef ANA_PSOURCE
!
!-----------------------------------------------------------------------
!  Set point Sources/Sinks position, direction, special flag, and mass
!  transport nondimensional shape profile with analytcal expressions.
!  Point sources are at U- and V-points. We need to get their positions
!  to process internal Land/Sea masking arrays during initialization.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        IF (luvsrc(ng).or.lwsrc(ng).or.any(ltracersrc(:,ng))) THEN
          DO tile=first_tile(ng),last_tile(ng),+1
            CALL ana_psource (ng, tile, inlm)
          END DO
        END IF
!$OMP BARRIER
      END DO
#endif
 
#if defined FOUR_DVAR || !defined TANGENT || !defined ADJOINT
!
!-----------------------------------------------------------------------
!  Read in initial forcing, climatology and assimilation data from
!  input NetCDF files.  It loads the first relevant data record for
!  the time-interpolation between snapshots.
!-----------------------------------------------------------------------
 
# ifdef ADJUST_BOUNDARY
!
!  If first pass of iteration loop, set time of open boundary
!  adjustment.
!
!$OMP MASTER
      IF (nrun.eq.erstr) THEN
        DO ng=1,ngrids
          obc_time(1,ng)=time(ng)
          DO irec=2,nbrec(ng)
            obc_time(irec,ng)=obc_time(irec-1,ng)+nobc(ng)*dt(ng)
          END DO
        END DO
      END IF
!$OMP END MASTER
!$OMP BARRIER
# endif
# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
!
!  If first pass of iteration loop, set time of surface forcing
!  adjustment.
!
!$OMP MASTER
      IF (nrun.eq.erstr) THEN
        DO ng=1,ngrids
          sf_time(1,ng)=time(ng)
          DO irec=2,nfrec(ng)
            sf_time(irec,ng)=sf_time(irec-1,ng)+nsff(ng)*dt(ng)
          END DO
        END DO
      END IF
!$OMP END MASTER
!$OMP BARRIER
# endif
# if !defined CORRELATION
!
!  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.
!
      DO ng=1,ngrids
!$OMP MASTER
        CALL close_inp (ng, inlm)
        CALL check_multifile (ng, inlm)
!$OMP END MASTER
#  ifdef DISTRIBUTE
        CALL mp_bcasti (ng, inlm, exit_flag)
#  endif
!$OMP BARRIER
        IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
      END DO
!
!  If applicable, read in input data.
!
      DO ng=1,ngrids
!$OMP MASTER
        CALL get_idata (ng)
        CALL get_data (ng)
!$OMP END MASTER
#  ifdef DISTRIBUTE
        CALL mp_bcasti (ng, inlm, exit_flag)
#  endif
!$OMP BARRIER
        IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
      END DO
# endif
#endif
 
#ifdef MASKING
!
!-----------------------------------------------------------------------
!  Set internal I/O mask arrays.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        DO tile=first_tile(ng),last_tile(ng),+1
          CALL set_masks (ng, tile, inlm)
        END DO
!$OMP BARRIER
      END DO
#endif
 
#if !defined CORRELATION
# ifdef NESTING
#  if defined MASKING || defined WET_DRY
!
!-----------------------------------------------------------------------
!  If nesting and Land/Sea masking, scale horizontal interpolation
!  weights to account for land contact points.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        CALL nesting (ng, inlm, nmask)
      END DO
#  endif
!
!-----------------------------------------------------------------------
!  If nesting, process state fields initial conditions in the contact
!  regions.
!-----------------------------------------------------------------------
!
!  Free-surface and 2D-momentum.
!
      DO nl=1,nestlayers
        DO ig=1,gridsinlayer(nl)
          ng=gridnumber(ig,nl)
          IF (any(compositegrid(:,ng))) THEN
            CALL nesting (ng, inlm, nfsic)        ! free-surface
#  ifndef SOLVE3D
            CALL nesting (ng, inlm, n2dic)        ! 2d momentum
#  endif
          END IF
        END DO
      END DO
 
#  ifdef SOLVE3D
!
!  Determine vertical indices and vertical interpolation weights in
!  the contact zone using initial unperturbed depth arrays.
!
      DO ng=1,ngrids
        CALL nesting (ng, inlm, nzwgt)
      END DO
!
!  3D-momentum and tracers.
!
      DO nl=1,nestlayers
        DO ig=1,gridsinlayer(nl)
          ng=gridnumber(ig,nl)
          IF (any(compositegrid(:,ng))) THEN
            CALL nesting (ng, inlm, n3dic)        ! 3D momentum
            CALL nesting (ng, inlm, ntvic)        ! Tracer variables
          END IF
        END DO
      END DO
#  endif
# endif
#endif
 
#if defined PROPAGATOR || \
   (defined masking    && (defined read_water || defined write_water ))
!
!-----------------------------------------------------------------------
!  Set variables associated with the processing water points and/or
!  size of packed state arrays.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        DO tile=first_tile(ng),last_tile(ng),+1
          CALL wpoints (ng, tile, inlm)
        END DO
!$OMP BARRIER
      END DO
#endif
 
#if defined NLM_OUTER                || \
    defined rbl4dvar                 || \
    defined rbl4dvar_ana_sensitivity || \
    defined rbl4dvar_fct_sensitivity || \
    defined tl_rbl4dvar
!
!-----------------------------------------------------------------------
!  Read in convolved adjoint impulse forcing (first record) and its
!  application time.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        IF (sporadicimpulse(ng)) THEN
          frcrec(ng)=1
!$OMP MASTER
          CALL get_state (ng, 7, 7, tlf(ng), frcrec(ng), 1)
!$OMP END MASTER
!$OMP BARRIER
# ifdef DISTRIBUTE
          CALL mp_bcasti (ng, itlm, exit_flag)
# endif
          IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
        END IF
      END DO
#endif
 
#if defined ANA_DRAG && defined UV_DRAG_GRID
!
!-----------------------------------------------------------------------
!  Set analytical spatially varying bottom friction parameter.
!-----------------------------------------------------------------------
!
      IF (nrun.eq.erstr) THEN
        DO ng=1,ngrids
          DO tile=first_tile(ng),last_tile(ng),+1
            CALL ana_drag (ng, tile, inlm)
          END DO
!$OMP BARRIER
        END DO
      END IF
#endif
!
!-----------------------------------------------------------------------
!  Compute grid stiffness.
!-----------------------------------------------------------------------
!
      IF (lstiffness) THEN
        lstiffness=.false.
        DO ng=1,ngrids
          DO tile=first_tile(ng),last_tile(ng),+1
            CALL stiffness (ng, tile, inlm)
          END DO
!$OMP BARRIER
        END DO
      END IF
 
#if defined FLOATS || defined STATIONS
!
!-----------------------------------------------------------------------
!  If applicable, convert initial locations to fractional grid
!  coordinates.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
!$OMP MASTER
        CALL grid_coords (ng, inlm)
!$OMP END MASTER
!$OMP BARRIER
      END DO
#endif
 
# if defined WAV_COUPLING && defined MCT_LIB
!
!-----------------------------------------------------------------------
!  Read in initial forcing from coupled wave model.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        DO tile=first_tile(ng),last_tile(ng),+1
          CALL ocn2wav_coupling (ng, tile)
        END DO
!$OMP BARRIER
        IF (master) WRITE (stdout,'(/)')
      END DO
# endif
!
!-----------------------------------------------------------------------
!  Initialize time-stepping counter and date/time string. Save NLM
!  initial conditions time.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        initime(ng)=time(ng)
        iic(ng)=ntstart(ng)
        CALL time_string (time(ng), time_code(ng))
      END DO
 
#ifdef PROFILE
!
!-----------------------------------------------------------------------
!  Turn off initialization time wall clock.
!-----------------------------------------------------------------------
!
      DO ng=1,ngrids
        DO thread=thread_range
          CALL wclock_off (ng, inlm, 2, __line__, myfile)
        END DO
!$OMP BARRIER
      END DO
#endif
!
      RETURN
      END SUBROUTINE initial