MODULE roms_kernel_mod ! !git $Id$ !svn $Id: so_semi_roms.h 1166 2023-05-17 20:11:58Z arango $ !================================================== Hernan G. Arango === ! Copyright (c) 2002-2023 The ROMS/TOMS Group Andrew M. Moore ! ! Licensed under a MIT/X style license ! ! See License_ROMS.md ! !======================================================================= ! ! ! ROMS/TOMS Stochastic Optimals, Seminorm Estimation Driver: ! ! ! ! This driver computes the eigenvectors of the stochastic optimals ! ! operator with respect the seminorm of the chosen functional. The ! ! stochastic optimals provide information about the influence of ! ! stochastic variations (biases) in ocean forcing. They can be ! ! used to build forecast ensembles. ! ! ! ! These routines control the initialization, time-stepping, and ! ! finalization of ROMS/TOMS model following ESMF conventions: ! ! ! ! ROMS_initialize ! ! ROMS_run ! ! ROMS_finalize ! ! ! ! Reference: ! ! ! ! Moore, A.M. et al., 2004: A comprehensive ocean prediction and ! ! analysis system based on the tangent linear and adjoint of a ! ! regional ocean model, Ocean Modelling, 7, 227-258. ! ! ! !======================================================================= ! USE mod_param USE mod_parallel USE mod_arrays USE mod_forces USE mod_iounits USE mod_ncparam USE mod_netcdf #if defined PIO_LIB && defined DISTRIBUTE USE mod_pio_netcdf #endif USE mod_scalars USE mod_stepping USE mod_storage ! USE propagator_mod ! USE ad_def_his_mod, ONLY : ad_def_his USE ad_wrt_his_mod, ONLY : ad_wrt_his USE close_io_mod, ONLY : close_file, close_inp, close_out #ifdef CHECKPOINTING USE def_gst_mod, ONLY : def_gst USE get_gst_mod, ONLY : get_gst #endif USE inp_par_mod, ONLY : inp_par #ifdef MCT_LIB # ifdef ATM_COUPLING USE mct_coupler_mod, ONLY : initialize_ocn2atm_coupling # endif # ifdef WAV_COUPLING USE mct_coupler_mod, ONLY : initialize_ocn2wav_coupling # endif #endif USE packing_mod, ONLY : ad_unpack USE packing_mod, ONLY : r_norm2 USE strings_mod, ONLY : FoundError #ifdef CHECKPOINTING USE wrt_gst_mod, ONLY : wrt_gst #endif USE wrt_rst_mod, ONLY : wrt_rst ! implicit none ! PRIVATE :: IRAM_error PUBLIC :: ROMS_initialize PUBLIC :: ROMS_run PUBLIC :: ROMS_finalize ! CONTAINS ! SUBROUTINE ROMS_initialize (first, mpiCOMM) ! !======================================================================= ! ! ! This routine allocates and initializes ROMS/TOMS state variables ! ! and internal and external parameters. ! ! ! !======================================================================= ! ! Imported variable declarations. ! logical, intent(inout) :: first ! integer, intent(in), optional :: mpiCOMM ! ! Local variable declarations. ! logical :: allocate_vars = .TRUE. ! #ifdef DISTRIBUTE integer :: MyError, MySize #endif integer :: chunk_size, ng, thread #ifdef _OPENMP integer :: my_threadnum #endif ! character (len=*), parameter :: MyFile = & & __FILE__//", ROMS_initialize" #ifdef DISTRIBUTE ! !----------------------------------------------------------------------- ! Set distribute-memory (mpi) world communictor. !----------------------------------------------------------------------- ! IF (PRESENT(mpiCOMM)) THEN OCN_COMM_WORLD=mpiCOMM ELSE OCN_COMM_WORLD=MPI_COMM_WORLD END IF CALL mpi_comm_rank (OCN_COMM_WORLD, MyRank, MyError) CALL mpi_comm_size (OCN_COMM_WORLD, MySize, MyError) #endif ! !----------------------------------------------------------------------- ! On first pass, initialize model parameters a variables for all ! nested/composed grids. Notice that the logical switch "first" ! is used to allow multiple calls to this routine during ensemble ! configurations. !----------------------------------------------------------------------- ! IF (first) THEN first=.FALSE. ! ! Initialize parallel control switches. These scalars switches are ! independent from standard input parameters. ! CALL initialize_parallel ! ! Read in model tunable parameters from standard input. Allocate and ! initialize variables in several modules after the number of nested ! grids and dimension parameters are known. ! CALL inp_par (iADM) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN ! ! Set domain decomposition tile partition range. This range is ! computed only once since the "first_tile" and "last_tile" values ! are private for each parallel thread/node. ! #if defined _OPENMP MyThread=my_threadnum() #elif defined DISTRIBUTE MyThread=MyRank #else MyThread=0 #endif DO ng=1,Ngrids chunk_size=(NtileX(ng)*NtileE(ng)+numthreads-1)/numthreads first_tile(ng)=MyThread*chunk_size last_tile (ng)=first_tile(ng)+chunk_size-1 END DO ! ! Initialize internal wall clocks. Notice that the timings does not ! includes processing standard input because several parameters are ! needed to allocate clock variables. ! IF (Master) THEN WRITE (stdout,10) 10 FORMAT (/,' Process Information:',/) END IF ! DO ng=1,Ngrids DO thread=THREAD_RANGE CALL wclock_on (ng, iADM, 0, __LINE__, MyFile) END DO END DO ! ! Allocate and initialize all model state arrays. ! CALL ROMS_allocate_arrays (allocate_vars) CALL ROMS_initialize_arrays IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF #if defined MCT_LIB && (defined ATM_COUPLING || defined WAV_COUPLING) ! !----------------------------------------------------------------------- ! Initialize coupling streams between model(s). !----------------------------------------------------------------------- ! DO ng=1,Ngrids # ifdef ATM_COUPLING CALL initialize_ocn2atm_coupling (ng, MyRank) # endif # ifdef WAV_COUPLING CALL initialize_ocn2wav_coupling (ng, MyRank) # endif END DO #endif ! !----------------------------------------------------------------------- ! Initialize adjoint model for all grids first in order to compute ! the size of the state vector, Nstate. This size is computed in ! routine "wpoints". !----------------------------------------------------------------------- #ifdef FORWARD_FLUXES ! ! Set the BLK structure to contain the nonlinear model surface fluxes ! needed by the tangent linear and adjoint models. Also, set switches ! to process that structure in routine "check_multifile". Notice that ! it is possible to split the solution into multiple NetCDF files to ! reduce their size. ! ! The switch LreadFRC is deactivated because all the atmospheric ! forcing, including shortwave radiation, is read from the NLM ! surface fluxes or is assigned during ESM coupling. Such fluxes ! are available from the QCK structure. There is no need for reading ! and processing from the FRC structure input forcing-files. ! CALL edit_multifile ('QCK2BLK') IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN DO ng=1,Ngrids LreadBLK(ng)=.TRUE. LreadFRC(ng)=.FALSE. END DO #endif ! ! Initialize perturbation tangent linear model. ! DO ng=1,Ngrids LreadFWD(ng)=.TRUE. CALL ad_initial (ng) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END DO ! ! Allocate arrays associated with Generalized Stability Theory (GST) ! analysis. ! CALL allocate_storage ! ! Initialize various IO flags. ! Nrun=0 DO ng=1,Ngrids LdefADJ(ng)=.TRUE. LwrtADJ(ng)=.TRUE. LdefTLM(ng)=.TRUE. LwrtTLM(ng)=.TRUE. LwrtPER(ng)=.FALSE. LcycleTLM(ng)=.FALSE. LcycleADJ(ng)=.FALSE. END DO ! ! Initialize ARPACK parameters. ! Lrvec=.TRUE. ! Compute Ritz vectors bmat='I' ! standard eigenvalue problem which='LM' ! compute NEV largest eigenvalues howmany='A' ! compute NEV Ritz vectors DO ng=1,Ngrids ido(ng)=0 ! reverse communication flag info(ng)=0 ! random initial residual vector iparam(1,ng)=1 ! exact shifts iparam(3,ng)=MaxIterGST ! maximum number of Arnoldi iterations iparam(4,ng)=1 ! block size in the recurrence iparam(7,ng)=1 ! type of eigenproblem being solved END DO ! ! ARPACK debugging parameters. ! logfil=stdout ! output logical unit ndigit=-3 ! number of decimal digits msaupd=1 ! iterations, timings, Ritz msaup2=1 ! norms, Ritz values msaitr=0 mseigt=0 msapps=0 msgets=0 mseupd=0 ! ! Determine size of the eigenproblem (Nsize) and size of work space ! array SworkL (LworkL). ! DO ng=1,Ngrids Nconv(ng)=0 Nsize(ng)=Nend(ng)-Nstr(ng)+1 END DO #ifdef CHECKPOINTING ! ! If restart, read in checkpointing data GST restart NetCDF file. ! Otherwise, create checkpointing restart NetCDF file. ! DO ng=1,Ngrids IF (LrstGST) THEN CALL get_gst (ng, iADM) ido(ng)=-2 laup2(1)=.FALSE. ! cnorm laup2(2)=.FALSE. ! getv0 laup2(3)=.FALSE. ! initv laup2(4)=.FALSE. ! update laup2(5)=.TRUE. ! ushift ELSE CALL def_gst (ng, iADM) END IF IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END DO #endif ! RETURN END SUBROUTINE ROMS_initialize ! SUBROUTINE ROMS_run (RunInterval) ! !======================================================================= ! ! ! This routine computes the eigenvectors of the stochastic optimal ! ! matrix, defined as: ! ! ! ! S = INTEGRAL [ R'(t,To) X R(t,To) dt ] from To to T ! ! ! ! where T is the forecast interval from the initial state at To, ! ! R(t,To) is the forward tangent propagator of the linearized ! ! dynamical model, R'(t,To) is the adjoint of R(t,To), and the ! ! matrix X defines the norm of interest. Here, X is a seminorm of ! ! the chosen functional. ! ! ! ! The above operator is symmetric and the ARPACK library is used ! ! to select eigenvectors and eigenvalues: ! ! ! ! Lehoucq, R.B., D.C. Sorensen, and C. Yang, 1997: ARPACK user's ! ! guide: solution of large scale eigenvalue problems with ! ! implicit restarted Arnoldi Methods, Rice University, 140p. ! ! ! !======================================================================= ! ! Imported variable declarations ! real(dp), intent(in) :: RunInterval ! seconds ! ! Local variable declarations. ! logical :: ITERATE #ifdef CHECKPOINTING logical :: LwrtGST #endif ! integer :: Fcount, Is, Ie integer :: i, iter, ng, tile integer :: NconvRitz(Ngrids) ! real(r8) :: Enorm ! TYPE (T_GST), allocatable :: ad_state(:) TYPE (T_GST), allocatable :: state(:) ! character (len=55) :: string character (len=*), parameter :: MyFile = & & __FILE__//", ROMS_run" ! !----------------------------------------------------------------------- ! Implicit Restarted Arnoldi Method (IRAM) for the computation of ! optimal perturbation Ritz eigenfunctions. !----------------------------------------------------------------------- ! ! Allocate nested grid pointers for state vectors. ! IF (.not.allocated(ad_state)) THEN allocate ( ad_state(Ngrids) ) END IF IF (.not.allocated(state)) THEN allocate ( state(Ngrids) ) END IF ! ! Iterate until either convergence or maximum iterations has been ! exceeded. ! iter=0 ITERATE=.TRUE. #ifdef CHECKPOINTING LwrtGST=.TRUE. #endif ! ITER_LOOP : DO WHILE (ITERATE) iter=iter+1 ! ! Reverse communication interface. ! DO ng=1,Ngrids #ifdef PROFILE CALL wclock_on (ng, iADM, 38, __LINE__, MyFile) #endif #ifdef DISTRIBUTE CALL pdsaupd (OCN_COMM_WORLD, & & ido(ng), bmat, Nsize(ng), which, NEV, & & Ritz_tol, & & STORAGE(ng)%resid(Nstr(ng)), NCV, & & STORAGE(ng)%Bvec(Nstr(ng),1), Nsize(ng), & & iparam(1,ng), ipntr(1,ng), & & STORAGE(ng)%SworkD, & & SworkL(1,ng), LworkL, info(ng)) #else CALL dsaupd (ido(ng), bmat, Nsize(ng), which, NEV, & & Ritz_tol, & & STORAGE(ng)%resid, NCV, & & STORAGE(ng)%Bvec, Nsize(ng), & & iparam(1,ng), ipntr(1,ng), & & STORAGE(ng)%SworkD, & & SworkL(1,ng), LworkL, info(ng)) #endif Nconv(ng)=iaup2(4) #ifdef PROFILE CALL wclock_off (ng, iADM, 38, __LINE__, MyFile) #endif #ifdef CHECKPOINTING ! ! If appropriate, write out check point data into GST restart NetCDF ! file. Notice that the restart data is always saved if MaxIterGST ! is reached without convergence. It is also saved when convergence ! is achieved (ido=99). ! IF ((MOD(iter,nGST).eq.0).or.(iter.ge.MaxIterGST).or. & & (ANY(ido.eq.99))) THEN CALL wrt_gst (ng, iADM) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN END IF #endif END DO ! ! Terminate computations if maximum number of iterations is reached. ! This will faciliate splitting the analysis in several computational ! cycles using the restart option. ! IF ((iter.ge.MaxIterGST).and.ANY(ido.ne.99)) THEN ITERATE=.FALSE. EXIT ITER_LOOP END IF ! ! Perform matrix-vector operation: R'(t,0)XR(0,t)u ! IF (ANY(ABS(ido).eq.1)) THEN DO ng=1,Ngrids Fcount=ADM(ng)%load ADM(ng)%Nrec(Fcount)=0 Fcount=TLM(ng)%load TLM(ng)%Nrec(Fcount)=0 ADM(ng)%Rindex=0 TLM(ng)%Rindex=0 END DO ! ! Set state vectors to process by the propagator via pointer ! equivalence. ! DO ng=1,Ngrids IF (ASSOCIATED(state(ng)%vector)) THEN nullify (state(ng)%vector) END IF Is=ipntr(1,ng) Ie=Is+Nsize(ng)-1 state(ng)%vector => STORAGE(ng)%SworkD(Is:Ie) IF (ASSOCIATED(ad_state(ng)%vector)) THEN nullify (ad_state(ng)%vector) END IF Is=ipntr(2,ng) Ie=Is+Nsize(ng)-1 ad_state(ng)%vector => STORAGE(ng)%SworkD(Is:Ie) END DO CALL propagator_so_semi (RunInterval, state, ad_state) IF (FoundError(exit_flag, NoError, __LINE__, MyFile)) RETURN ELSE IF (ANY(info.ne.0)) THEN DO ng=1,Ngrids IF (info(ng).ne.0) THEN IF (Master) THEN CALL IRAM_error (info(ng), string) WRITE (stdout,10) 'DSAUPD', TRIM(string), & & ', info = ', info(ng) END IF RETURN END IF END DO ELSE ! ! Compute Ritz vectors (the only choice left is IDO=99). They are ! generated in ARPACK in decreasing magnitude of its eigenvalue. ! The most significant is first. ! DO ng=1,Ngrids NconvRitz(ng)=iparam(5,ng) IF (Master) THEN WRITE (stdout,20) 'Number of converged Ritz values:', & & iparam(5,ng) WRITE (stdout,20) 'Number of Arnoldi iterations:', & & iparam(3,ng) END IF #ifdef PROFILE CALL wclock_on (ng, iADM, 38, __LINE__, MyFile) #endif #ifdef DISTRIBUTE CALL pdseupd (OCN_COMM_WORLD, & & Lrvec, howmany, pick(1,ng), & & RvalueR(1,ng), & & STORAGE(ng)%Rvector(Nstr(ng),1), & & Nsize(ng), sigmaR, & & bmat, Nsize(ng), which, NEV, Ritz_tol, & & STORAGE(ng)%resid(Nstr(ng)), NCV, & & STORAGE(ng)%Bvec(Nstr(ng),1), Nsize(ng), & & iparam(1,ng), ipntr(1,ng), & & STORAGE(ng)%SworkD, & & SworkL(1,ng), LworkL, info(ng)) #else CALL dseupd (Lrvec, howmany, pick(1,ng), & & RvalueR(1,ng), & & STORAGE(ng)%Rvector, Nsize(ng), & & sigmaR, bmat, Nsize, which, NEV, Ritz_tol, & & STORAGE(ng)%resid, NCV, & & STORAGE(ng)%Bvec, Nsize(ng), & & iparam(1,ng), ipntr(1,ng), & & STORAGE(ng)%SworkD, & & SworkL(1,ng), LworkL, info(ng)) #endif #ifdef PROFILE CALL wclock_off (ng, iADM, 38, __LINE__, MyFile) #endif END DO IF (ANY(info.ne.0)) THEN DO ng=1,Ngrids IF (info(ng).ne.0) THEN IF (Master) THEN CALL IRAM_error (info(ng), string) WRITE (stdout,10) 'DSEUPD', TRIM(string), & & ', info = ', info(ng) END IF RETURN END IF END DO ELSE ! ! Close existing adjoint NetCDF file. ! SourceFile=MyFile DO ng=1,Ngrids CALL close_file (ng, iADM, ADM(ng), ADM(ng)%name) IF (FoundError(exit_flag, NoError, & & __LINE__, MyFile)) RETURN END DO ! ! Clear forcing arrays. ! DO ng=1,Ngrids DO tile=first_tile(ng),last_tile(ng),+1 CALL initialize_forces (ng, tile, 0) END DO END DO ! ! Activate writing of each eigenvector into tangent history NetCDF ! file. The "ocean_time" is the eigenvector number. ! Nrun=0 DO i=1,MAXVAL(NconvRitz) DO ng=1,Ngrids IF ((i.eq.1).or.LmultiGST) THEN Fcount=ADM(ng)%load ADM(ng)%Nrec(Fcount)=0 ADM(ng)%Rindex=0 LwrtADJ(ng)=.TRUE. LwrtState2d(ng)=.TRUE. END IF IF (LmultiGST) THEN LdefADJ(ng)=.TRUE. WRITE (ADM(ng)%name,30) TRIM(ADM(ng)%base), i CALL ad_def_his (ng, LdefADJ(ng)) IF (FoundError(exit_flag, NoError, & & __LINE__, MyFile)) RETURN ELSE IF (.not.LmultiGST.and.(i.eq.1)) THEN LdefADJ(ng)=.TRUE. ADM(ng)%name=TRIM(ADM(ng)%base)//'_ritz.nc' CALL ad_def_his (ng, LdefADJ(ng)) IF (FoundError(exit_flag, NoError, & & __LINE__, MyFile)) RETURN END IF END DO ! ! Compute and write the Ritz eigenvectors. ! DO ng=1,Ngrids Is=Nstr(ng) Ie=Nend(ng) IF (ASSOCIATED(state(ng)%vector)) THEN nullify (state(ng)%vector) END IF IF (ASSOCIATED(ad_state(ng)%vector)) THEN nullify (ad_state(ng)%vector) END IF state(ng)%vector => STORAGE(ng)%Rvector(Is:Ie,i) ad_state(ng)%vector => SworkR(Is:Ie) END DO CALL propagator_so_semi (RunInterval, state, ad_state) IF (FoundError(exit_flag, NoError, & & __LINE__, MyFile)) RETURN ! ! Unpack surface forcing eigenvectors from Rvector and write into ! nonlinear history file. ! DO ng=1,Ngrids DO tile=first_tile(ng),last_tile(ng),+1 CALL ad_unpack (ng, tile, Nstr(ng), Nend(ng), & & state(ng)%vector) END DO #ifdef DISTRIBUTE CALL ad_wrt_his (ng, MyRank) #else CALL ad_wrt_his (ng, -1) #endif IF (FoundError(exit_flag, NoError, & & __LINE__, MyFile)) RETURN END DO ! ! Compute Euclidean norm. ! DO ng=1,Ngrids CALL r_norm2 (ng, iADM, Nstr(ng), Nend(ng), & & -RvalueR(i,ng), & & state(ng)%vector, & & ad_state(ng)%vector, Enorm) norm(i,ng)=Enorm IF (Master) THEN WRITE (stdout,40) i, norm(i,ng), RvalueR(i,ng), i END IF END DO ! ! Write out Ritz eigenvalues and Ritz eigenvector Euclidean norm ! (residual) to nonlinear history NetCDF file. ! SourceFile=MyFile DO ng=1,Ngrids SELECT CASE (ADM(ng)%IOtype) CASE (io_nf90) CALL netcdf_put_fvar (ng, iADM, & & ADM(ng)%name, & & 'Ritz_rvalue', & & RvalueR(i:,ng), & & start = (/ADM(ng)%Rindex/), & & total = (/1/), & & ncid = ADM(ng)%ncid) IF (FoundError(exit_flag, NoError, & & __LINE__, MyFile)) RETURN ! CALL netcdf_put_fvar (ng, iADM, ADM(ng)%name, & & 'Ritz_norm', & & norm(i:,ng), & & start = (/ADM(ng)%Rindex/), & & total = (/1/), & & ncid = ADM(ng)%ncid) IF (FoundError(exit_flag, NoError, & & __LINE__, MyFile)) RETURN IF ((i.eq.1).or.LmultiGST) THEN CALL netcdf_put_fvar (ng, iADM, & & ADM(ng)%name, & & 'SO_trace', & & TRnorm(ng), & & start = (/0/), & & total = (/0/), & & ncid = ADM(ng)%ncid) IF (FoundError(exit_flag, NoError, & & __LINE__, MyFile)) RETURN END IF #if defined PIO_LIB && defined DISTRIBUTE CASE (io_pio) CALL pio_netcdf_put_fvar (ng, iADM, & & ADM(ng)%name, & & 'Ritz_rvalue', & & RvalueR(i:,ng), & & start = (/ADM(ng)%Rindex/), & & total = (/1/), & & pioFile = ADM(ng)%pioFile) IF (FoundError(exit_flag, NoError, & & __LINE__, MyFile)) RETURN ! CALL pio_netcdf_put_fvar (ng, iADM, & & ADM(ng)%name, & & 'Ritz_norm', & & norm(i:,ng), & & start = (/ADM(ng)%Rindex/), & & total = (/1/), & & pioFile = ADM(ng)%pioFile) IF (FoundError(exit_flag, NoError, & & __LINE__, MyFile)) RETURN IF ((i.eq.1).or.LmultiGST) THEN CALL pio_netcdf_put_fvar (ng, iADM, & & ADM(ng)%name, & & 'SO_trace', & & TRnorm(ng), & & start = (/0/), & & total = (/0/), & & pioFile = ADM(ng)%pioFile) IF (FoundError(exit_flag, NoError, & & __LINE__, MyFile)) RETURN END IF #endif END SELECT ! IF (LmultiGST) THEN CALL close_file (ng, iADM, ADM(ng), ADM(ng)%name) IF (FoundError(exit_flag, NoError, & & __LINE__, MyFile)) RETURN END IF END DO END DO END IF END IF ITERATE=.FALSE. END IF END DO ITER_LOOP ! 10 FORMAT (/,1x,'Error in ',a,1x,a,a,1x,i5,/) 20 FORMAT (/,a,1x,i2,/) 30 FORMAT (a,'_',i3.3,'.nc') 40 FORMAT (1x,i4.4,'-th residual',1p,e14.6,0p, & & ' Ritz value',1pe14.6,0p,2x,i4.4) ! RETURN END SUBROUTINE ROMS_run ! SUBROUTINE ROMS_finalize ! !======================================================================= ! ! ! This routine terminates ROMS/TOMS nonlinear and adjoint models ! ! execution. ! ! ! !======================================================================= ! ! Local variable declarations. ! integer :: Fcount, ng, thread ! character (len=*), parameter :: MyFile = & & __FILE__//", ROMS_finalize" ! !----------------------------------------------------------------------- ! If blowing-up, save latest model state into RESTART NetCDF file. !----------------------------------------------------------------------- ! ! If cycling restart records, write solution into the next record. ! IF (exit_flag.eq.1) THEN DO ng=1,Ngrids IF (LwrtRST(ng)) THEN IF (Master) WRITE (stdout,10) 10 FORMAT (/,' Blowing-up: Saving latest model state into ', & & ' RESTART file',/) Fcount=RST(ng)%load IF (LcycleRST(ng).and.(RST(ng)%Nrec(Fcount).ge.2)) THEN RST(ng)%Rindex=2 LcycleRST(ng)=.FALSE. END IF blowup=exit_flag exit_flag=NoError #ifdef DISTRIBUTE CALL wrt_rst (ng, MyRank) #else CALL wrt_rst (ng, -1) #endif END IF END DO END IF ! !----------------------------------------------------------------------- ! Stop model and time profiling clocks, report memory requirements, and ! close output NetCDF files. !----------------------------------------------------------------------- ! ! Stop time clocks. ! IF (Master) THEN WRITE (stdout,20) 20 FORMAT (/,'Elapsed wall CPU time for each process (seconds):',/) END IF ! DO ng=1,Ngrids DO thread=THREAD_RANGE CALL wclock_off (ng, iADM, 0, __LINE__, MyFile) END DO END DO ! ! Report dynamic memory and automatic memory requirements. ! CALL memory ! ! Close IO files. ! DO ng=1,Ngrids CALL close_inp (ng, iADM) END DO CALL close_out ! RETURN END SUBROUTINE ROMS_finalize ! SUBROUTINE IRAM_error (info, string) ! !======================================================================= ! ! ! This routine decodes internal error messages from the Implicit ! ! Restarted Arnoldi Method (IRAM) for the computation of optimal ! ! perturbation Ritz eigenfunctions. ! ! ! !======================================================================= ! ! imported variable declarations. ! integer, intent(in) :: info ! character (len=*), intent(out) :: string ! !----------------------------------------------------------------------- ! Decode error message from IRAM. !----------------------------------------------------------------------- ! IF (info.eq.0) THEN string='Normal exit ' ELSE IF (info.eq.1) THEN string='Maximum number of iterations taken ' ELSE IF (info.eq.3) THEN string='No shifts could be applied during an IRAM cycle ' ELSE IF (info.eq.-1) THEN string='Nstate must be positive ' ELSE IF (info.eq.-2) THEN string='NEV must be positive ' ELSE IF (info.eq.-3) THEN string='NCV must be greater NEV and less than or equal Nstate ' ELSE IF (info.eq.-4) THEN string='Maximum number of iterations must be greater than zero ' ELSE IF (info.eq.-5) THEN string='WHICH must be one of LM, SM, LA, SA or BE ' ELSE IF (info.eq.-6) THEN string='BMAT must be one of I or G ' ELSE IF (info.eq.-7) THEN string='Length of private work array SworkL is not sufficient ' ELSE IF (info.eq.-8) THEN string='Error in DSTEQR in the eigenvalue calculation ' ELSE IF (info.eq.-9) THEN string='Starting vector is zero ' ELSE IF (info.eq.-10) THEN string='IPARAM(7) must be 1, 2, 3, 4, 5 ' ELSE IF (info.eq.-11) THEN string='IPARAM(7) = 1 and BMAT = G are incompatable ' ELSE IF (info.eq.-12) THEN string='IPARAM(1) must be equal to 0 or 1 ' ELSE IF (info.eq.-13) THEN string='NEV and WHICH = BE are incompatable ' ELSE IF (info.eq.-14) THEN string='Did not find any eigenvalues to sufficient accuaracy ' ELSE IF (info.eq.-15) THEN string='HOWMANY must be one of A or S if RVEC = .TRUE. ' ELSE IF (info.eq.-16) THEN string='HOWMANY = S not yet implemented ' ELSE IF (info.eq.-17) THEN string='Different count of converge Ritz values in DSEUPD ' ELSE IF (info.eq.-9999) THEN string='Could not build and Arnoldi factorization ' END IF ! RETURN END SUBROUTINE IRAM_error END MODULE roms_kernel_mod