SUBROUTINE wrt_rst (ng) ! !svn $Id: wrt_rst.F 331 2009-03-12 00:34:51Z arango $ !================================================== Hernan G. Arango === ! Copyright (c) 2002-2009 The ROMS/TOMS Group ! ! Licensed under a MIT/X style license ! ! See License_ROMS.txt ! !======================================================================= ! ! ! This routine writes fields into restart NetCDF file. ! ! ! !======================================================================= ! USE mod_param USE mod_parallel USE mod_grid USE mod_iounits USE mod_mixing USE mod_ncparam USE mod_netcdf USE mod_ocean USE mod_scalars USE mod_sediment USE mod_stepping ! USE nf_fwrite2d_mod, ONLY : nf_fwrite2d USE nf_fwrite3d_mod, ONLY : nf_fwrite3d ! implicit none ! ! Imported variable declarations. ! integer, intent(in) :: ng ! ! Local variable declarations. ! integer :: LBi, UBi, LBj, UBj integer :: gfactor, gtype, i, itrc, status, varid integer :: ntmp(1) real(r8) :: scale ! LBi=LBOUND(GRID(ng)%h,DIM=1) UBi=UBOUND(GRID(ng)%h,DIM=1) LBj=LBOUND(GRID(ng)%h,DIM=2) UBj=UBOUND(GRID(ng)%h,DIM=2) ! SourceFile='wrt_rst.F' ! !----------------------------------------------------------------------- ! Write out restart fields. !----------------------------------------------------------------------- ! IF (exit_flag.ne.NoError) RETURN ! ! Set grid type factor to write full (gfactor=1) fields or water ! points (gfactor=-1) fields only. ! gfactor=1 ! ! Set time record index. ! tRSTindx(ng)=tRSTindx(ng)+1 NrecRST(ng)=NrecRST(ng)+1 ! ! If requested, set time index to recycle time records in restart ! file. ! IF (LcycleRST(ng)) THEN tRSTindx(ng)=MOD(tRSTindx(ng)-1,2)+1 END IF ! ! Write out model time (s). ! CALL netcdf_put_fvar (ng, iNLM, RSTname(ng), & & TRIM(Vname(idtime,ng)), time(ng:), & & (/tRSTindx(ng)/), (/1/), & & ncid = ncRSTid(ng), & & varid = rstVid(idtime,ng)) IF (exit_flag.ne.NoError) RETURN ! ! Write out wet/dry mask at RHO-points. ! scale=1.0_r8 gtype=gfactor*r2dvar status=nf_fwrite2d(ng, iNLM, ncRSTid(ng), rstVid(idRwet,ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % rmask, & & GRID(ng) % rmask_wet, & & SetFillVal = .FALSE.) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idRwet)), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF ! ! Write out wet/dry mask at U-points. ! scale=1.0_r8 gtype=gfactor*u2dvar status=nf_fwrite2d(ng, iNLM, ncRSTid(ng), rstVid(idUwet,ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % umask, & & GRID(ng) % umask_wet, & & SetFillVal = .FALSE.) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idUwet)), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF ! ! Write out wet/dry mask at V-points. ! scale=1.0_r8 gtype=gfactor*v2dvar status=nf_fwrite2d(ng, iNLM, ncRSTid(ng), rstVid(idVwet,ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % vmask, & & GRID(ng) % vmask_wet, & & SetFillVal = .FALSE.) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idVwet)), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF ! ! Write out free-surface (m). ! scale=1.0_r8 gtype=gfactor*r2dvar status=nf_fwrite2d(ng, iNLM, ncRSTid(ng), rstVid(idFsur,ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % rmask, & & OCEAN(ng) % zeta(:,:,kstp(ng)), & & SetFillVal = .FALSE.) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idFsur)), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF ! ! Write out 2D momentum component (m/s) in the XI-direction. ! scale=1.0_r8 gtype=gfactor*u2dvar status=nf_fwrite2d(ng, iNLM, ncRSTid(ng), rstVid(idUbar,ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % umask, & & OCEAN(ng) % ubar(:,:,kstp(ng))) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idUbar)), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF ! ! Write out 2D momentum component (m/s) in the ETA-direction. ! scale=1.0_r8 gtype=gfactor*v2dvar status=nf_fwrite2d(ng, iNLM, ncRSTid(ng), rstVid(idVbar,ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % vmask, & & OCEAN(ng) % vbar(:,:,kstp(ng))) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idVbar)), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF ! ! Write out 3D momentum component (m/s) in the XI-direction. ! scale=1.0_r8 gtype=gfactor*u3dvar status=nf_fwrite3d(ng, iNLM, ncRSTid(ng), rstVid(idUvel,ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, 1, N(ng), scale, & & GRID(ng) % umask, & & OCEAN(ng) % u(:,:,:,nrhs(ng))) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idUvel)), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF ! ! Write out momentum component (m/s) in the ETA-direction. ! scale=1.0_r8 gtype=gfactor*v3dvar status=nf_fwrite3d(ng, iNLM, ncRSTid(ng), rstVid(idVvel,ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, 1, N(ng), scale, & & GRID(ng) % vmask, & & OCEAN(ng) % v(:,:,:,nrhs(ng))) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idVvel)), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF ! ! Write out tracer type variables. ! DO itrc=1,NT(ng) scale=1.0_r8 gtype=gfactor*r3dvar status=nf_fwrite3d(ng, iNLM, ncRSTid(ng), rstTid(itrc,ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, 1, N(ng), scale, & & GRID(ng) % rmask, & & OCEAN(ng) % t(:,:,:,nrhs(ng),itrc)) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idTvar(itrc))), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF END DO ! ! Write out density anomaly. ! scale=1.0_r8 gtype=gfactor*r3dvar status=nf_fwrite3d(ng, iNLM, ncRSTid(ng), rstVid(idDano,ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, 1, N(ng), scale, & & GRID(ng) % rmask, & & OCEAN(ng) % rho) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idDano)), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF ! ! Write out vertical viscosity coefficient. ! scale=1.0_r8 gtype=gfactor*w3dvar status=nf_fwrite3d(ng, iNLM, ncRSTid(ng), rstVid(idVvis,ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, 0, N(ng), scale, & & GRID(ng) % rmask, & & MIXING(ng) % Akv) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idVvis)), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF ! ! Write out vertical diffusion coefficient for potential temperature. ! scale=1.0_r8 gtype=gfactor*w3dvar status=nf_fwrite3d(ng, iNLM, ncRSTid(ng), rstVid(idTdif,ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, 0, N(ng), scale, & & GRID(ng) % rmask, & & MIXING(ng) % Akt(:,:,:,itemp)) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idTdif)), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF ! ! Write out vertical diffusion coefficient for salinity. ! scale=1.0_r8 gtype=gfactor*w3dvar status=nf_fwrite3d(ng, iNLM, ncRSTid(ng), rstVid(idSdif,ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, 0, N(ng), scale, & & GRID(ng) % rmask, & & MIXING(ng) % Akt(:,:,:,isalt)) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idSdif)), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF ! ! Write out sediment fraction of each size class in each bed layer. ! DO i=1,NST scale=1.0_r8 gtype=gfactor*b3dvar status=nf_fwrite3d(ng, iNLM, ncRSTid(ng), rstVid(idfrac(i),ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, 1, Nbed, scale, & & GRID(ng) % rmask, & & OCEAN(ng) % bed_frac(:,:,:,i)) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idfrac(i))), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF END DO ! ! Write out sediment mass of each size class in each bed layer. ! DO i=1,NST scale=1.0_r8 gtype=gfactor*b3dvar status=nf_fwrite3d(ng, iNLM, ncRSTid(ng), rstVid(idBmas(i),ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, 1, Nbed, scale, & & GRID(ng) % rmask, & & OCEAN(ng) % bed_mass(:,:,:,nrhs(ng),i)) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idBmas(i))), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF END DO ! ! Write out sediment properties in each bed layer. ! DO i=1,MBEDP IF (i.eq.itauc) THEN scale=rho0 ELSE scale=1.0_r8 END IF gtype=gfactor*b3dvar status=nf_fwrite3d(ng, iNLM, ncRSTid(ng), rstVid(idSbed(i),ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, 1, Nbed, scale, & & GRID(ng) % rmask, & & OCEAN(ng) % bed(:,:,:,i)) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idSbed(i))), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF END DO ! ! Write out exposed sediment layer properties. Notice that only the ! first four properties (mean grain diameter, mean grain density, ! mean settling velocity, mean critical erosion stress, ! ripple length and ripple height) are written. ! DO i=1,6 scale=1.0_r8 gtype=gfactor*r2dvar status=nf_fwrite2d(ng, iNLM, ncRSTid(ng), rstVid(idBott(i),ng), & & tRSTindx(ng), gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % rmask, & & OCEAN(ng) % bottom(:,:,i)) IF (status.ne.nf90_noerr) THEN IF (Master) THEN WRITE (stdout,10) TRIM(Vname(1,idBott(i))), tRSTindx(ng) END IF exit_flag=3 ioerror=status RETURN END IF END DO ! !----------------------------------------------------------------------- ! Synchronize restart NetCDF file to disk. !----------------------------------------------------------------------- ! CALL netcdf_sync (ng, iNLM, RSTname(ng), ncRSTid(ng)) IF (exit_flag.ne.NoError) RETURN IF (Master) WRITE (stdout,20) kstp(ng), nrhs(ng), tRSTindx(ng) ! 10 FORMAT (/,' WRT_RST - error while writing variable: ',a,/,11x, & & 'into restart NetCDF file for time record: ',i4) 20 FORMAT (6x,'WRT_RST - wrote re-start fields (Index=', i1, & & ',',i1,') into time record = ',i7.7) RETURN END SUBROUTINE wrt_rst