MODULE ad_wrt_his_mod ! !git $Id$ !svn $Id: ad_wrt_his.F 1189 2023-08-15 21:26:58Z arango $ !================================================== Hernan G. Arango === ! Copyright (c) 2002-2023 The ROMS/TOMS Group ! ! Licensed under a MIT/X style license ! ! See License_ROMS.md ! !======================================================================= ! ! ! This routine writes requested adjoint model fields into adjoint ! ! history NetCDF file. ! ! ! ! Notice that only momentum is affected by the full time-averaged ! ! masks. If applicable, these mask contains information about ! ! river runoff and time-dependent wetting and drying variations. ! ! ! !======================================================================= ! USE mod_param USE mod_parallel USE mod_forces USE mod_fourdvar USE mod_grid USE mod_iounits USE mod_mixing USE mod_ncparam USE mod_ocean USE mod_scalars USE mod_stepping ! USE nf_fwrite2d_mod, ONLY : nf_fwrite2d USE nf_fwrite3d_mod, ONLY : nf_fwrite3d USE omega_mod, ONLY : scale_omega USE strings_mod, ONLY : FoundError ! implicit none ! PUBLIC :: ad_wrt_his PRIVATE :: ad_wrt_his_nf90 ! CONTAINS ! !*********************************************************************** SUBROUTINE ad_wrt_his (ng, tile) !*********************************************************************** ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile ! ! Local variable declarations. ! integer :: LBi, UBi, LBj, UBj ! character (len=*), parameter :: MyFile = & & "ROMS/Adjoint/ad_wrt_his.F" ! !----------------------------------------------------------------------- ! Write out history fields according to IO type. !----------------------------------------------------------------------- ! LBi=BOUNDS(ng)%LBi(tile) UBi=BOUNDS(ng)%UBi(tile) LBj=BOUNDS(ng)%LBj(tile) UBj=BOUNDS(ng)%UBj(tile) ! SELECT CASE (ADM(ng)%IOtype) CASE (io_nf90) CALL ad_wrt_his_nf90 (ng, tile, & & LBi, UBi, LBj, UBj) CASE DEFAULT IF (Master) WRITE (stdout,10) ADM(ng)%IOtype exit_flag=3 END SELECT IF (FoundError(exit_flag, NoError, 116, MyFile)) RETURN ! 10 FORMAT (' AD_WRT_HIS - Illegal output type, io_type = ',i0, & & /,14x,'Check KeyWord ''OUT_LIB'' in ''roms.in''.') ! RETURN END SUBROUTINE ad_wrt_his ! !*********************************************************************** SUBROUTINE ad_wrt_his_nf90 (ng, tile, & & LBi, UBi, LBj, UBj) !*********************************************************************** ! USE mod_netcdf ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile integer, intent(in) :: LBi, UBi, LBj, UBj ! ! Local variable declarations. ! integer :: Fcount, i, j, gfactor, gtype, status integer :: kout integer :: kfout integer :: itrc, k, nout ! real(dp) :: scale real(r8) :: Tval(1) ! real(r8), allocatable :: Wr3d(:,:,:) ! character (len=*), parameter :: MyFile = & & "ROMS/Adjoint/ad_wrt_his.F"//", ad_wrt_his_nf90" ! SourceFile=MyFile ! !----------------------------------------------------------------------- ! Write out adjoint fields. !----------------------------------------------------------------------- ! IF (FoundError(exit_flag, NoError, 169, MyFile)) RETURN ! ! Set grid type factor to write full (gfactor=1) fields or water ! points (gfactor=-1) fields only. ! gfactor=1 ! ! Determine time index to write. The "nout" index is updated to the ! version of "ad_main3d" that updates the "iic" counter at the bottom. ! Therefore, we need to change the conditional "iic(ng).ne.ntend(ng)" ! to "iic(ng).gt.ntend(ng)" to get identical solutions. ! kout=kstp(ng) kfout=2 IF (iic(ng).gt.ntend(ng)) THEN nout=nnew(ng) ELSE nout=nstp(ng) END IF ! ! Set time record index. ! ADM(ng)%Rindex=ADM(ng)%Rindex+1 Fcount=ADM(ng)%load ADM(ng)%Nrec(Fcount)=ADM(ng)%Nrec(Fcount)+1 ! ! Report. ! IF (Master) WRITE (stdout,10) kout, nout, ADM(ng)%Rindex ! ! If requested, set time index to recycle time records in the adjoint ! file. ! IF (LcycleADJ(ng)) THEN ADM(ng)%Rindex=MOD(ADM(ng)%Rindex-1,2)+1 END IF ! ! Write out model time (s). ! IF (LwrtTime(ng)) THEN IF (LwrtPER(ng)) THEN Tval(1)=REAL(ADM(ng)%Rindex,r8)*day2sec ELSE Tval(1)=ForceTime(ng) END IF CALL netcdf_put_fvar (ng, iADM, ADM(ng)%name, & & TRIM(Vname(1,idtime)), tval, & & (/ADM(ng)%Rindex/), (/1/), & & ncid = ADM(ng)%ncid, & & varid = ADM(ng)%Vid(idtime)) IF (FoundError(exit_flag, NoError, 253, MyFile)) RETURN END IF ! ! Write out bathymetry. ! IF (Hout(idbath,ng)) THEN scale=1.0_dp gtype=gfactor*r2dvar status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idbath, & & ADM(ng)%Vid(idbath), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % rmask, & & GRID(ng)% ad_h, & & SetFillVal = .FALSE.) IF (FoundError(status, nf90_noerr, 346, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idbath)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF ! ! Write out free-surface (m). ! IF (Hout(idFsur,ng)) THEN IF (WRTforce(ng)) THEN scale=1.0_dp gtype=gfactor*r2dvar status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idFsur, & & ADM(ng)%Vid(idFsur), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % rmask, & & OCEAN(ng)% f_zetaG(:,:,kfout)) IF (FoundError(status, nf90_noerr, 371, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idFsur)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF ELSE scale=1.0_dp gtype=gfactor*r2dvar IF (LwrtState2d(ng)) THEN status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idFsur, & & ADM(ng)%Vid(idFsur), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % rmask, & & OCEAN(ng)% ad_zeta(:,:,kout)) ELSE status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idFsur, & & ADM(ng)%Vid(idFsur), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % rmask, & & OCEAN(ng)% ad_zeta_sol) ENDIF IF (FoundError(status, nf90_noerr, 412, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idFsur)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF END IF ! ! Write out 2D U-momentum component (m/s). ! IF (Hout(idUbar,ng)) THEN IF (WRTforce(ng)) THEN scale=1.0_dp gtype=gfactor*u2dvar status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idUbar, & & ADM(ng)%Vid(idUbar), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % umask_full, & & OCEAN(ng) % f_ubarG(:,:,kfout)) IF (FoundError(status, nf90_noerr, 465, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idUbar)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF ELSE scale=1.0_dp gtype=gfactor*u2dvar IF (LwrtState2d(ng)) THEN status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idUbar, & & ADM(ng)%Vid(idUbar), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % umask_full, & & OCEAN(ng) % ad_ubar(:,:,kout)) ELSE status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idUbar, & & ADM(ng)%Vid(idUbar), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % umask_full, & & OCEAN(ng) % ad_ubar_sol) END IF IF (FoundError(status, nf90_noerr, 496, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idUbar)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF END IF ! ! Write out 2D V-momentum component (m/s). ! IF (Hout(idVbar,ng)) THEN IF (WRTforce(ng)) THEN scale=1.0_dp gtype=gfactor*v2dvar status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idVbar, & & ADM(ng)%Vid(idVbar), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % vmask_full, & & OCEAN(ng) % f_vbarG(:,:,kfout)) IF (FoundError(status, nf90_noerr, 549, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idVbar)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF ELSE scale=1.0_dp gtype=gfactor*v2dvar IF (LwrtState2d(ng)) THEN status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idVbar, & & ADM(ng)%Vid(idVbar), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % vmask_full, & & OCEAN(ng) % ad_vbar(:,:,kout)) ELSE status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idVbar, & & ADM(ng)%Vid(idVbar), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % vmask_full, & & OCEAN(ng) % ad_vbar_sol) END IF IF (FoundError(status, nf90_noerr, 580, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idVbar)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF END IF ! ! Write out 3D U-momentum component (m/s). ! IF (Hout(idUvel,ng)) THEN IF (WRTforce(ng)) THEN scale=1.0_dp gtype=gfactor*u3dvar status=nf_fwrite3d(ng, iADM, ADM(ng)%ncid, idUvel, & & ADM(ng)%Vid(idUvel), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, 1, N(ng), scale, & & GRID(ng) % umask_full, & & OCEAN(ng) % f_uG(:,:,:,kfout)) IF (FoundError(status, nf90_noerr, 635, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idUvel)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF ELSE scale=1.0_dp gtype=gfactor*u3dvar status=nf_fwrite3d(ng, iADM, ADM(ng)%ncid, idUvel, & & ADM(ng)%Vid(idUvel), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, 1, N(ng), scale, & & GRID(ng) % umask_full, & & OCEAN(ng) % ad_u(:,:,:,nout)) IF (FoundError(status, nf90_noerr, 670, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idUvel)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF END IF ! ! Write out 3D V-momentum component (m/s). ! IF (Hout(idVvel,ng)) THEN IF (WRTforce(ng)) THEN scale=1.0_dp gtype=gfactor*v3dvar status=nf_fwrite3d(ng, iADM, ADM(ng)%ncid, idVvel, & & ADM(ng)%Vid(idVvel), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, 1, N(ng), scale, & & GRID(ng) % vmask_full, & & OCEAN(ng) % f_vG(:,:,:,kfout)) IF (FoundError(status, nf90_noerr, 723, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idVvel)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF ELSE scale=1.0_dp gtype=gfactor*v3dvar status=nf_fwrite3d(ng, iADM, ADM(ng)%ncid, idVvel, & & ADM(ng)%Vid(idVvel), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, 1, N(ng), scale, & & GRID(ng) % vmask_full, & & OCEAN(ng) % ad_v(:,:,:,nout)) IF (FoundError(status, nf90_noerr, 758, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idVvel)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF END IF ! ! Write out S-coordinate omega vertical velocity (m/s). ! IF (Hout(idOvel,ng)) THEN IF (.not.allocated(Wr3d)) THEN allocate (Wr3d(LBi:UBi,LBj:UBj,0:N(ng))) Wr3d(LBi:UBi,LBj:UBj,0:N(ng))=0.0_r8 END IF scale=1.0_dp gtype=gfactor*w3dvar CALL scale_omega (ng, tile, LBi, UBi, LBj, UBj, 0, N(ng), & & GRID(ng) % pm, & & GRID(ng) % pn, & & OCEAN(ng) % ad_W_sol, & & Wr3d) status=nf_fwrite3d(ng, iADM, ADM(ng)%ncid, idOvel, & & ADM(ng)%Vid(idOvel), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, 0, N(ng), scale, & & GRID(ng) % rmask, & & Wr3d) IF (FoundError(status, nf90_noerr, 818, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idOvel)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF deallocate (Wr3d) END IF ! ! Write out tracer type variables. ! DO itrc=1,NT(ng) IF (Hout(idTvar(itrc),ng)) THEN IF (WRTforce(ng)) THEN scale=1.0_dp gtype=gfactor*r3dvar status=nf_fwrite3d(ng, iADM, ADM(ng)%ncid, idTvar(itrc), & & ADM(ng)%Tid(itrc), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, 1, N(ng), scale, & & GRID(ng) % rmask, & & OCEAN(ng) % f_tG(:,:,:,kfout,itrc)) IF (FoundError(status, nf90_noerr, 845, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idTvar(itrc))), & & ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF ELSE scale=1.0_dp gtype=gfactor*r3dvar status=nf_fwrite3d(ng, iADM, ADM(ng)%ncid, idTvar(itrc), & & ADM(ng)%Tid(itrc), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, 1, N(ng), scale, & & GRID(ng) % rmask, & & OCEAN(ng) % ad_t(:,:,:,nout,itrc)) IF (FoundError(status, nf90_noerr, 881, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idTvar(itrc))), & & ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF END IF END DO ! ! Write out density anomaly. ! IF (Hout(idDano,ng)) THEN scale=1.0_dp gtype=gfactor*r3dvar status=nf_fwrite3d(ng, iADM, ADM(ng)%ncid, idDano, & & ADM(ng)%Vid(idDano), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, 1, N(ng), scale, & & GRID(ng) % rmask, & & OCEAN(ng) % ad_rho) IF (FoundError(status, nf90_noerr, 937, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idDano)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF ! ! Write out vertical viscosity coefficient. ! IF (Hout(idVvis,ng)) THEN scale=1.0_dp gtype=gfactor*w3dvar status=nf_fwrite3d(ng, iADM, ADM(ng)%ncid, idVvis, & & ADM(ng)%Vid(idVvis), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, 0, N(ng), scale, & & GRID(ng) % rmask, & & MIXING(ng) % ad_Akv, & & SetFillVal = .FALSE.) IF (FoundError(status, nf90_noerr, 961, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idVvis)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF ! ! Write out vertical diffusion coefficient for potential temperature. ! IF (Hout(idTdif,ng)) THEN scale=1.0_dp gtype=gfactor*w3dvar status=nf_fwrite3d(ng, iADM, ADM(ng)%ncid, idTdif, & & ADM(ng)%Vid(idTdif), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, 0, N(ng), scale, & & GRID(ng) % rmask, & & MIXING(ng) % ad_Akt(:,:,:,itemp), & & SetFillVal = .FALSE.) IF (FoundError(status, nf90_noerr, 985, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idTdif)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF ! ! Write out vertical diffusion coefficient for salinity. ! IF (Hout(idSdif,ng)) THEN scale=1.0_dp gtype=gfactor*w3dvar status=nf_fwrite3d(ng, iADM, ADM(ng)%ncid, idSdif, & & ADM(ng)%Vid(idSdif), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, 0, N(ng), scale, & & GRID(ng) % rmask, & & MIXING(ng) % ad_Akt(:,:,:,isalt), & & SetFillVal = .FALSE.) IF (FoundError(status, nf90_noerr, 1010, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idSdif)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF ! ! Write out net surface active tracer fluxes. ! DO itrc=1,NT(ng) IF (Hout(idTsur(itrc),ng)) THEN scale=1.0_dp gtype=gfactor*r2dvar status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idTsur(itrc), & & ADM(ng)%Vid(idTsur(itrc)), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % rmask, & & FORCES(ng) % ad_stflx(:,:,itrc)) IF (FoundError(status, nf90_noerr, 1046, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idTsur(itrc))), & & ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF END DO ! ! Write out surface U-momentum stress. ! IF (Hout(idUsms,ng)) THEN scale=1.0_dp gtype=gfactor*u2dvar status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idUsms, & & ADM(ng)%Vid(idUsms), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % umask, & & FORCES(ng) % ad_sustr) IF (FoundError(status, nf90_noerr, 1080, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idUsms)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF ! ! Write out surface V-momentum stress. ! IF (Hout(idVsms,ng)) THEN scale=1.0_dp gtype=gfactor*v2dvar status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idVsms, & & ADM(ng)%Vid(idVsms), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % vmask, & & FORCES(ng) % ad_svstr) IF (FoundError(status, nf90_noerr, 1109, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idVsms)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF ! ! Write out bottom U-momentum stress. ! IF (Hout(idUbms,ng)) THEN scale=1.0_dp gtype=gfactor*u2dvar status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idUbms, & & ADM(ng)%Vid(idUbms), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % umask, & & FORCES(ng) % ad_bustr_sol) IF (FoundError(status, nf90_noerr, 1134, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idUbms)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF ! ! Write out bottom V-momentum stress. ! IF (Hout(idVbms,ng)) THEN scale=1.0_dp gtype=gfactor*v2dvar status=nf_fwrite2d(ng, iADM, ADM(ng)%ncid, idVbms, & & ADM(ng)%Vid(idVbms), & & ADM(ng)%Rindex, gtype, & & LBi, UBi, LBj, UBj, scale, & & GRID(ng) % vmask, & & FORCES(ng) % ad_bvstr_sol) IF (FoundError(status, nf90_noerr, 1158, MyFile)) THEN IF (Master) THEN WRITE (stdout,20) TRIM(Vname(1,idVbms)), ADM(ng)%Rindex END IF exit_flag=3 ioerror=status RETURN END IF END IF ! !----------------------------------------------------------------------- ! Synchronize adjoint history NetCDF file to disk to allow other ! processes to access data immediately after it is written. !----------------------------------------------------------------------- ! CALL netcdf_sync (ng, iADM, ADM(ng)%name, ADM(ng)%ncid) IF (FoundError(exit_flag, NoError, 1174, MyFile)) RETURN ! 10 FORMAT (2x,'AD_WRT_HIS_NF90 - writing adjoint', t42, & & 'fields (Index=',i1,',',i1,') in record = ',i0) 20 FORMAT (/,' AD_WRT_HIS_NF90 - error while writing variable: ',a, & & /,19x,'into adjoint NetCDF file for time record: ',i0) ! RETURN END SUBROUTINE ad_wrt_his_nf90 END MODULE ad_wrt_his_mod