ROMS
Loading...
Searching...
No Matches
so_roms.h
Go to the documentation of this file.
1 MODULE roms_kernel_mod
2!
3!git $Id$
4!================================================== Hernan G. Arango ===
5! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md !
8!=======================================================================
9! !
10! ROMS Stochastic Optimal Driver: !
11! !
12! These routines control the initialization, time-stepping, and !
13! finalization of ROMS model following ESMF conventions: !
14! !
15! ROMS_initialize !
16! ROMS_run !
17! ROMS_finalize !
18! !
19! Reference: !
20! !
21! Moore, A.M. et al., 2004: A comprehensive ocean prediction and !
22! analysis system based on the tangent linear and adjoint of a !
23! regional ocean model, Ocean Modelling, 7, 227-258. !
24! !
25!=======================================================================
26!
27 USE mod_param
28 USE mod_parallel
29 USE mod_arrays
30 USE mod_iounits
31 USE mod_ncparam
32 USE mod_netcdf
33#if defined PIO_LIB && defined DISTRIBUTE
35#endif
36 USE mod_scalars
37 USE mod_stepping
38 USE mod_storage
39!
41!
43#ifdef CHECKPOINTING
44 USE def_gst_mod, ONLY : def_gst
45 USE get_gst_mod, ONLY : get_gst
46#endif
47 USE inp_par_mod, ONLY : inp_par
48#ifdef MCT_LIB
49# ifdef ATM_COUPLING
50 USE mct_coupler_mod, ONLY : initialize_ocn2atm_coupling
51# endif
52# ifdef WAV_COUPLING
53 USE mct_coupler_mod, ONLY : initialize_ocn2wav_coupling
54# endif
55#endif
56 USE packing_mod, ONLY : r_norm2
58 USE strings_mod, ONLY : founderror
59#ifdef CHECKPOINTING
60 USE wrt_gst_mod, ONLY : wrt_gst
61#endif
62 USE wrt_rst_mod, ONLY : wrt_rst
63!
64 implicit none
65!
66 PRIVATE :: iram_error
67 PUBLIC :: roms_initialize
68 PUBLIC :: roms_run
69 PUBLIC :: roms_finalize
70!
71 CONTAINS
72!
73 SUBROUTINE roms_initialize (first, mpiCOMM)
74!
75!=======================================================================
76! !
77! This routine allocates and initializes ROMS state variables !
78! and internal and external parameters. !
79! !
80!=======================================================================
81!
82! Imported variable declarations.
83!
84 logical, intent(inout) :: first
85!
86 integer, intent(in), optional :: mpiCOMM
87!
88! Local variable declarations.
89!
90 logical :: allocate_vars = .true.
91
92#ifdef DISTRIBUTE
93 integer :: MyError, MySize
94#endif
95 integer :: chunk_size, ng, thread
96#ifdef _OPENMP
97 integer :: my_threadnum
98#endif
99!
100 character (len=*), parameter :: MyFile = &
101 & __FILE__//", ROMS_initialize"
102
103#ifdef DISTRIBUTE
104!
105!-----------------------------------------------------------------------
106! Set distribute-memory (mpi) world communictor.
107!-----------------------------------------------------------------------
108!
109 IF (PRESENT(mpicomm)) THEN
110 ocn_comm_world=mpicomm
111 ELSE
112 ocn_comm_world=mpi_comm_world
113 END IF
114 CALL mpi_comm_rank (ocn_comm_world, myrank, myerror)
115 CALL mpi_comm_size (ocn_comm_world, mysize, myerror)
116#endif
117!
118!-----------------------------------------------------------------------
119! On first pass, initialize model parameters a variables for all
120! nested/composed grids. Notice that the logical switch "first"
121! is used to allow multiple calls to this routine during ensemble
122! configurations.
123!-----------------------------------------------------------------------
124!
125 IF (first) THEN
126 first=.false.
127!
128! Initialize parallel control switches. These scalars switches are
129! independent from standard input parameters.
130!
132!
133! Set the ROMS standard output unit to write verbose execution info.
134! Notice that the default standard out unit in Fortran is 6.
135!
136! In some applications like coupling or disjointed mpi-communications,
137! it is advantageous to write standard output to a specific filename
138! instead of the default Fortran standard output unit 6. If that is
139! the case, it opens such formatted file for writing.
140!
141 IF (set_stdoutunit) THEN
143 set_stdoutunit=.false.
144 END IF
145!
146! Read in model tunable parameters from standard input. Allocate and
147! initialize variables in several modules after the number of nested
148! grids and dimension parameters are known.
149!
150 CALL inp_par (itlm)
151 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
152!
153! Set domain decomposition tile partition range. This range is
154! computed only once since the "first_tile" and "last_tile" values
155! are private for each parallel thread/node.
156!
157#if defined _OPENMP
158 mythread=my_threadnum()
159#elif defined DISTRIBUTE
161#else
162 mythread=0
163#endif
164 DO ng=1,ngrids
165 chunk_size=(ntilex(ng)*ntilee(ng)+numthreads-1)/numthreads
166 first_tile(ng)=mythread*chunk_size
167 last_tile(ng)=first_tile(ng)+chunk_size-1
168 END DO
169!
170! Initialize internal wall clocks. Notice that the timings does not
171! includes processing standard input because several parameters are
172! needed to allocate clock variables.
173!
174 IF (master) THEN
175 WRITE (stdout,10)
176 10 FORMAT (/,' Process Information:',/)
177 END IF
178!
179 DO ng=1,ngrids
180 DO thread=thread_range
181 CALL wclock_on (ng, itlm, 0, __line__, myfile)
182 END DO
183 END DO
184!
185! Allocate and initialize modules variables.
186!
187 CALL roms_allocate_arrays (allocate_vars)
189 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
190
191 END IF
192
193#if defined MCT_LIB && (defined ATM_COUPLING || defined WAV_COUPLING)
194!
195!-----------------------------------------------------------------------
196! Initialize coupling streams between model(s).
197!-----------------------------------------------------------------------
198!
199 DO ng=1,ngrids
200# ifdef ATM_COUPLING
201 CALL initialize_ocn2atm_coupling (ng, myrank)
202# endif
203# ifdef WAV_COUPLING
204 CALL initialize_ocn2wav_coupling (ng, myrank)
205# endif
206 END DO
207#endif
208!
209!-----------------------------------------------------------------------
210! Initialize tangent linear for all grids first in order to compute
211! the size of the state vector, Nstate. This size is computed in
212! routine "wpoints".
213!-----------------------------------------------------------------------
214
215#ifdef FORWARD_FLUXES
216!
217! Set the BLK structure to contain the nonlinear model surface fluxes
218! needed by the tangent linear and adjoint models. Also, set switches
219! to process that structure in routine "check_multifile". Notice that
220! it is possible to split the solution into multiple NetCDF files to
221! reduce their size.
222!
223! The switch LreadFRC is deactivated because all the atmospheric
224! forcing, including shortwave radiation, is read from the NLM
225! surface fluxes or is assigned during ESM coupling. Such fluxes
226! are available from the QCK structure. There is no need for reading
227! and processing from the FRC structure input forcing-files.
228!
229 CALL edit_multifile ('QCK2BLK')
230 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
231 DO ng=1,ngrids
232 lreadblk(ng)=.true.
233 lreadfrc(ng)=.false.
234 END DO
235#endif
236!
237! Initialize perturbation tangent linear model.
238!
239 DO ng=1,ngrids
240 lreadfwd(ng)=.true.
241 CALL tl_initial (ng)
242 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
243 END DO
244!
245! Allocate arrays associated with Generalized Stability Theory (GST)
246! analysis.
247!
249!
250! Initialize various IO flags.
251!
252 nrun=0
253 DO ng=1,ngrids
254 lwrtadj(ng)=.true.
255 lwrthis(ng)=.true.
256 lwrttlm(ng)=.true.
257 lwrtper(ng)=.false.
258 lcycletlm(ng)=.false.
259 lcycleadj(ng)=.false.
260 nadj(ng)=ntimes(ng)/nintervals
261 ntlm(ng)=ntimes(ng)
262 END DO
263!
264! Initialize ARPACK parameters.
265!
266 lrvec=.true. ! Compute Ritz vectors
267 bmat='I' ! standard eigenvalue problem
268 which='LM' ! compute NEV largest eigenvalues
269 howmany='A' ! compute NEV Ritz vectors
270 DO ng=1,ngrids
271 ido(ng)=0 ! reverse communication flag
272 info(ng)=0 ! random initial residual vector
273 iparam(1,ng)=1 ! exact shifts
274 iparam(3,ng)=maxitergst ! maximum number of Arnoldi iterations
275 iparam(4,ng)=1 ! block size in the recurrence
276 iparam(7,ng)=1 ! type of eigenproblem being solved
277 END DO
278!
279! ARPACK debugging parameters.
280!
281 logfil=stdout ! output logical unit
282 ndigit=-3 ! number of decimal digits
283 msaupd=1 ! iterations, timings, Ritz
284 msaup2=1 ! norms, Ritz values
285 msaitr=0
286 mseigt=0
287 msapps=0
288 msgets=0
289 mseupd=0
290!
291! Determine size of the eigenproblem (Nsize) and size of work space
292! array SworkL (LworkL).
293!
294 DO ng=1,ngrids
295 nconv(ng)=0
296 nsize(ng)=nend(ng)-nstr(ng)+1
297 END DO
298
299#ifdef CHECKPOINTING
300!
301! If restart, read in checkpointing data GST restart NetCDF file.
302! Otherwise, create checkpointing restart NetCDF file.
303!
304 DO ng=1,ngrids
305 IF (lrstgst) THEN
306 CALL get_gst (ng, itlm)
307 ido(ng)=-2
308 laup2(1)=.false. ! cnorm
309 laup2(2)=.false. ! getv0
310 laup2(3)=.false. ! initv
311 laup2(4)=.false. ! update
312 laup2(5)=.true. ! ushift
313 ELSE
314 CALL def_gst (ng, itlm)
315 END IF
316 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
317 END DO
318#endif
319!
320 RETURN
321 END SUBROUTINE roms_initialize
322!
323 SUBROUTINE roms_run (RunInterval)
324!
325!=======================================================================
326! !
327! This routine computes the singular vectors of R(0,t) by a single !
328! integration of a perturbation "u" forward in time with the !
329! tangent linear model over [0,t], multiplication of the result !
330! by "X", followed by an integration of the result backwards in !
331! time with the adjoint model over [t,0]. This is equivalmet to !
332! the matrix-vector operation: !
333! !
334! transpose[R(t,0)] X R(0,t) u !
335! !
336! The above operator is symmetric and the ARPACK library is used !
337! to select eigenvectors and eigenvalues: !
338! !
339! Lehoucq, R.B., D.C. Sorensen, and C. Yang, 1997: ARPACK user's !
340! guide: solution of large scale eigenvalue problems with !
341! implicit restarted Arnoldi Methods, Rice University, 140p. !
342! !
343!=======================================================================
344!
345! Imported variable declarations
346!
347 real(dp), intent(in) :: RunInterval ! seconds
348!
349! Local variable declarations.
350!
351 logical :: ITERATE
352#ifdef CHECKPOINTING
353 logical :: LwrtGST
354#endif
355!
356 integer :: Is, Ie
357 integer :: i, iter, ng, tile
358 integer :: NconvRitz(Ngrids)
359!
360 real(r8) :: Enorm
361!
362 TYPE (T_GST), allocatable :: ad_state(:)
363 TYPE (T_GST), allocatable :: state(:)
364!
365 character (len=55) :: string
366
367 character (len=*), parameter :: MyFile = &
368 & __FILE__//", ROMS_run"
369!
370!-----------------------------------------------------------------------
371! Implicit Restarted Arnoldi Method (IRAM) for the computation of
372! optimal perturbation Ritz eigenfunctions.
373!-----------------------------------------------------------------------
374!
375! Allocate nested grid pointers for state vectors.
376!
377 IF (.not.allocated(ad_state)) THEN
378 allocate ( ad_state(ngrids) )
379 END IF
380 IF (.not.allocated(state)) THEN
381 allocate ( state(ngrids) )
382 END IF
383!
384! Iterate until either convergence or maximum iterations has been
385! exceeded.
386!
387 iter=0
388 iterate=.true.
389#ifdef CHECKPOINTING
390 lwrtgst=.true.
391#endif
392!
393 iter_loop : DO WHILE (iterate)
394 iter=iter+1
395!
396! Reverse communication interface.
397!
398 DO ng=1,ngrids
399#ifdef PROFILE
400 CALL wclock_on (ng, itlm, 38, __line__, myfile)
401#endif
402#ifdef DISTRIBUTE
403 CALL pdsaupd (ocn_comm_world, &
404 & ido(ng), bmat, nsize(ng), which, nev, &
405 & ritz_tol, &
406 & storage(ng)%resid(nstr(ng)), ncv, &
407 & storage(ng)%Bvec(nstr(ng),1), nsize(ng), &
408 & iparam(1,ng), ipntr(1,ng), &
409 & storage(ng)%SworkD, &
410 & sworkl(1,ng), lworkl, info(ng))
411#else
412 CALL dsaupd (ido(ng), bmat, nsize(ng), which, nev, &
413 & ritz_tol, &
414 & storage(ng)%resid, ncv, &
415 & storage(ng)%Bvec, nsize(ng), &
416 & iparam(1,ng), ipntr(1,ng), &
417 & storage(ng)%SworkD, &
418 & sworkl(1,ng), lworkl, info(ng))
419#endif
420 nconv(ng)=iaup2(4)
421#ifdef PROFILE
422 CALL wclock_off (ng, itlm, 38, __line__, myfile)
423#endif
424#ifdef CHECKPOINTING
425!
426! If appropriate, write out check point data into GST restart NetCDF
427! file. Notice that the restart data is always saved if MaxIterGST
428! is reached without convergence. It is also saved when convergence
429! is achieved (ido=99).
430!
431 IF ((mod(iter,ngst).eq.0).or.(iter.ge.maxitergst).or. &
432 & (any(ido.eq.99))) THEN
433 CALL wrt_gst (ng, itlm)
434 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
435 END IF
436#endif
437 END DO
438!
439! Terminate computations if maximum number of iterations is reached.
440! This will faciliate splitting the analysis in several computational
441! cycles using the restart option.
442!
443 IF ((iter.ge.maxitergst).and.any(ido.ne.99)) THEN
444 iterate=.false.
445 EXIT iter_loop
446 END IF
447!
448! Perform matrix-vector operation: R`(t,0)XR(0,t)u
449!
450 IF (any(abs(ido).eq.1)) THEN
451!
452! Set state vectors to process by the propagator via pointer
453! equivalence.
454!
455 DO ng=1,ngrids
456 IF (ASSOCIATED(state(ng)%vector)) THEN
457 nullify (state(ng)%vector)
458 END IF
459 is=ipntr(1,ng)
460 ie=is+nsize(ng)-1
461 state(ng)%vector => storage(ng)%SworkD(is:ie)
462
463 IF (ASSOCIATED(ad_state(ng)%vector)) THEN
464 nullify (ad_state(ng)%vector)
465 END IF
466 is=ipntr(2,ng)
467 ie=is+nsize(ng)-1
468 ad_state(ng)%vector => storage(ng)%SworkD(is:ie)
469 END DO
470
471 CALL propagator_so (runinterval, iter, state, ad_state)
472 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
473 ELSE
474 IF (any(info.ne.0)) THEN
475 DO ng=1,ngrids
476 IF (info(ng).ne.0) THEN
477 IF (master) THEN
478 CALL iram_error (info(ng), string)
479 WRITE (stdout,10) 'DSAUPD', trim(string), &
480 & ', info = ', info(ng)
481 END IF
482 RETURN
483 END IF
484 END DO
485 ELSE
486!
487! Compute Ritz vectors (the only choice left is IDO=99). They are
488! generated in ARPACK in decreasing magnitude of its eigenvalue.
489! The most significant is first.
490!
491 DO ng=1,ngrids
492 nconvritz(ng)=iparam(5,ng)
493 IF (master) THEN
494 WRITE (stdout,20) 'Number of converged Ritz values:', &
495 & iparam(5,ng)
496 WRITE (stdout,20) 'Number of Arnoldi iterations:', &
497 & iparam(3,ng)
498 END IF
499#ifdef PROFILE
500 CALL wclock_on (ng, itlm, 38, __line__, myfile)
501#endif
502#ifdef DISTRIBUTE
503 CALL pdseupd (ocn_comm_world, &
504 & lrvec, howmany, pick(1,ng), &
505 & rvaluer(1,ng), &
506 & storage(ng)%Rvector(nstr(ng),1), &
507 & nsize(ng), sigmar, &
508 & bmat, nsize(ng), which, nev, ritz_tol, &
509 & storage(ng)%resid(nstr(ng)), ncv, &
510 & storage(ng)%Bvec(nstr(ng),1), nsize(ng), &
511 & iparam(1,ng), ipntr(1,ng), &
512 & storage(ng)%SworkD, &
513 & sworkl(1,ng), lworkl, info(ng))
514#else
515 CALL dseupd (lrvec, howmany, pick(1,ng), &
516 & rvaluer(1,ng), &
517 & storage(ng)%Rvector, nsize(ng), &
519 & storage(ng)%resid, ncv, &
520 & storage(ng)%Bvec, nsize(ng), &
521 & iparam(1,ng), ipntr(1,ng), &
522 & storage(ng)%SworkD, &
523 & sworkl(1,ng), lworkl, info(ng))
524#endif
525#ifdef PROFILE
526 CALL wclock_off (ng, itlm, 38, __line__, myfile)
527#endif
528 END DO
529
530 IF (any(info.ne.0)) THEN
531 DO ng=1,ngrids
532 IF (info(ng).ne.0) THEN
533 IF (master) THEN
534 CALL iram_error (info(ng), string)
535 WRITE (stdout,10) 'DSEUPD', trim(string), &
536 & ', info = ', info(ng)
537 END IF
538 RETURN
539 END IF
540 END DO
541 ELSE
542!
543! Activate writing of each eigenvector into single or multiple tangent
544! linear history NetCDF files.
545!
546 nrun=0
547
548 DO i=1,maxval(nconvritz)
549 DO ng=1,ngrids
550 IF (lmultigst) THEN
551 WRITE (tlm(ng)%name,30) trim(tlm(ng)%base), i
552 END IF
553 END DO
554!
555! Compute and write Ritz eigenvectors.
556!
557 DO ng=1,ngrids
558 is=nstr(ng)
559 ie=nend(ng)
560 IF (ASSOCIATED(state(ng)%vector)) THEN
561 nullify (state(ng)%vector)
562 END IF
563
564 IF (ASSOCIATED(ad_state(ng)%vector)) THEN
565 nullify (ad_state(ng)%vector)
566 END IF
567 state(ng)%vector => storage(ng)%Rvector(is:ie,i)
568 ad_state(ng)%vector => sworkr(is:ie)
569 END DO
570
571 CALL propagator_so (runinterval, -i, state, ad_state)
573 & __line__, myfile)) RETURN
574!
575 DO ng=1,ngrids
576 CALL r_norm2 (ng, itlm, nstr(ng), nend(ng), &
577 & -rvaluer(i,ng), &
578 & state(ng)%vector, &
579 & ad_state(ng)%vector, enorm)
580 norm(i,ng)=enorm
581 IF (master) THEN
582 WRITE (stdout,40) i, norm(i,ng), rvaluer(i,ng), i
583 END IF
584 END DO
585!
586! Write out Ritz eigenvalues and Ritz eigenvector Euclidean norm
587! (residual) to NetCDF file(s). Notice that we write the same value
588! twice in the TLM file for the initial and final perturbation of
589! the eigenvector.
590!
591 sourcefile=myfile
592 DO ng=1,ngrids
593 SELECT CASE (tlm(ng)%IOtype)
594 CASE (io_nf90)
595 CALL netcdf_put_fvar (ng, itlm, &
596 & tlm(ng)%name, &
597 & 'Ritz_rvalue', &
598 & rvaluer(i:,ng), &
599 & start = (/tlm(ng)%Rindex/), &
600 & total = (/1/), &
601 & ncid = tlm(ng)%ncid)
602
604 & __line__, myfile)) RETURN
605!
606 CALL netcdf_put_fvar (ng, itlm, &
607 & tlm(ng)%name, &
608 & 'Ritz_norm', &
609 & norm(i:,ng), &
610 & start = (/tlm(ng)%Rindex/), &
611 & total = (/1/), &
612 & ncid = tlm(ng)%ncid)
613
615 & __line__, myfile)) RETURN
616
617#if defined PIO_LIB && defined DISTRIBUTE
618 CASE (io_pio)
619 CALL pio_netcdf_put_fvar (ng, itlm, &
620 & tlm(ng)%name, &
621 & 'Ritz_rvalue', &
622 & rvaluer(i:,ng), &
623 & start = (/tlm(ng)%Rindex/), &
624 & total = (/1/), &
625 & piofile = tlm(ng)%pioFile)
626
628 & __line__, myfile)) RETURN
629!
630 CALL pio_netcdf_put_fvar (ng, itlm, &
631 & tlm(ng)%name, &
632 & 'Ritz_norm', &
633 & norm(i:,ng), &
634 & start = (/tlm(ng)%Rindex/), &
635 & total = (/1/), &
636 & piofile = tlm(ng)%pioFile)
637
639 & __line__, myfile)) RETURN
640#endif
641 END SELECT
642!
643 IF (lmultigst) THEN
644 CALL close_file (ng, itlm, tlm(ng), tlm(ng)%name)
646 & __line__, myfile)) RETURN
647 END IF
648!
649 SELECT CASE (adm(ng)%IOtype)
650 CASE (io_nf90)
651 CALL netcdf_put_fvar (ng, iadm, &
652 & adm(ng)%name, &
653 & 'Ritz_rvalue', &
654 & rvaluer(i:,ng), &
655 & start = (/adm(ng)%Rindex/), &
656 & total = (/1/), &
657 & ncid = adm(ng)%ncid)
658
660 & __line__, myfile)) RETURN
661!
662 CALL netcdf_put_fvar (ng, iadm, &
663 & adm(ng)%name, &
664 & 'Ritz_norm', &
665 & norm(i:,ng), &
666 & start = (/adm(ng)%Rindex/), &
667 & total = (/1/), &
668 & ncid = adm(ng)%ncid)
669
671 & __line__, myfile)) RETURN
672
673#if defined PIO_LIB && defined DISTRIBUTE
674 CASE (io_pio)
675 CALL pio_netcdf_put_fvar (ng, iadm, &
676 & adm(ng)%name, &
677 & 'Ritz_rvalue', &
678 & rvaluer(i:,ng), &
679 & start = (/adm(ng)%Rindex/), &
680 & total = (/1/), &
681 & piofile = adm(ng)%pioFile)
682
684 & __line__, myfile)) RETURN
685!
686 CALL pio_netcdf_put_fvar (ng, iadm, &
687 & adm(ng)%name, &
688 & 'Ritz_norm', &
689 & norm(i:,ng), &
690 & start = (/adm(ng)%Rindex/), &
691 & total = (/1/), &
692 & piofile = adm(ng)%pioFile)
693
695 & __line__, myfile)) RETURN
696#endif
697 END SELECT
698 END DO
699 END DO
700 END IF
701 END IF
702 iterate=.false.
703 END IF
704
705 END DO iter_loop
706!
707 10 FORMAT (/,1x,'Error in ',a,1x,a,a,1x,i5,/)
708 20 FORMAT (/,a,1x,i2,/)
709 30 FORMAT (a,'_',i3.3,'.nc')
710 40 FORMAT (1x,i4.4,'-th residual',1p,e14.6,0p, &
711 & ' Ritz value',1pe14.6,0p,2x,i4.4)
712!
713 RETURN
714 END SUBROUTINE roms_run
715!
716 SUBROUTINE roms_finalize
717!
718!=======================================================================
719! !
720! This routine terminates ROMS nonlinear and adjoint models !
721! execution. !
722! !
723!=======================================================================
724!
725! Local variable declarations.
726!
727 integer :: Fcount, ng, thread
728!
729 character (len=*), parameter :: MyFile = &
730 & __FILE__//", ROMS_finalize"
731!
732!-----------------------------------------------------------------------
733! If blowing-up, save latest model state into RESTART NetCDF file.
734!-----------------------------------------------------------------------
735!
736! If cycling restart records, write solution into the next record.
737!
738 IF (exit_flag.eq.1) THEN
739 DO ng=1,ngrids
740 IF (lwrtrst(ng)) THEN
741 IF (master) WRITE (stdout,10)
742 10 FORMAT (/,' Blowing-up: Saving latest model state into ', &
743 & ' RESTART file',/)
744 fcount=rst(ng)%load
745 IF (lcyclerst(ng).and.(rst(ng)%Nrec(fcount).ge.2)) THEN
746 rst(ng)%Rindex=2
747 lcyclerst(ng)=.false.
748 END IF
751#ifdef DISTRIBUTE
752 CALL wrt_rst (ng, myrank)
753#else
754 CALL wrt_rst (ng, -1)
755#endif
756 END IF
757 END DO
758 END IF
759!
760!-----------------------------------------------------------------------
761! Stop model and time profiling clocks, report memory requirements, and
762! close output NetCDF files.
763!-----------------------------------------------------------------------
764!
765! Stop time clocks.
766!
767 IF (master) THEN
768 WRITE (stdout,20)
769 20 FORMAT (/,'Elapsed wall CPU time for each process (seconds):',/)
770 END IF
771!
772 DO ng=1,ngrids
773 DO thread=thread_range
774 CALL wclock_off (ng, itlm, 0, __line__, myfile)
775 END DO
776 END DO
777!
778! Report dynamic memory and automatic memory requirements.
779!
780 CALL memory
781!
782! Close IO files.
783!
784 DO ng=1,ngrids
785 CALL close_inp (ng, itlm)
786 END DO
787 CALL close_out
788!
789 RETURN
790 END SUBROUTINE roms_finalize
791!
792 SUBROUTINE iram_error (info, string)
793!
794!=======================================================================
795! !
796! This routine decodes internal error messages from the Implicit !
797! Restarted Arnoldi Method (IRAM) for the computation of optimal !
798! perturbation Ritz eigenfunctions. !
799! !
800!=======================================================================
801!
802! imported variable declarations.
803!
804 integer, intent(in) :: info
805!
806 character (len=*), intent(out) :: string
807!
808!-----------------------------------------------------------------------
809! Decode error message from IRAM.
810!-----------------------------------------------------------------------
811!
812 IF (info.eq.0) THEN
813 string='Normal exit '
814 ELSE IF (info.eq.1) THEN
815 string='Maximum number of iterations taken '
816 ELSE IF (info.eq.3) THEN
817 string='No shifts could be applied during an IRAM cycle '
818 ELSE IF (info.eq.-1) THEN
819 string='Nstate must be positive '
820 ELSE IF (info.eq.-2) THEN
821 string='NEV must be positive '
822 ELSE IF (info.eq.-3) THEN
823 string='NCV must be greater NEV and less than or equal Nstate '
824 ELSE IF (info.eq.-4) THEN
825 string='Maximum number of iterations must be greater than zero '
826 ELSE IF (info.eq.-5) THEN
827 string='WHICH must be one of LM, SM, LA, SA or BE '
828 ELSE IF (info.eq.-6) THEN
829 string='BMAT must be one of I or G '
830 ELSE IF (info.eq.-7) THEN
831 string='Length of private work array SworkL is not sufficient '
832 ELSE IF (info.eq.-8) THEN
833 string='Error in DSTEQR in the eigenvalue calculation '
834 ELSE IF (info.eq.-9) THEN
835 string='Starting vector is zero '
836 ELSE IF (info.eq.-10) THEN
837 string='IPARAM(7) must be 1, 2, 3, 4, 5 '
838 ELSE IF (info.eq.-11) THEN
839 string='IPARAM(7) = 1 and BMAT = G are incompatable '
840 ELSE IF (info.eq.-12) THEN
841 string='IPARAM(1) must be equal to 0 or 1 '
842 ELSE IF (info.eq.-13) THEN
843 string='NEV and WHICH = BE are incompatable '
844 ELSE IF (info.eq.-14) THEN
845 string='Did not find any eigenvalues to sufficient accuaracy '
846 ELSE IF (info.eq.-15) THEN
847 string='HOWMANY must be one of A or S if RVEC = .TRUE. '
848 ELSE IF (info.eq.-16) THEN
849 string='HOWMANY = S not yet implemented '
850 ELSE IF (info.eq.-17) THEN
851 string='Different count of converge Ritz values in DSEUPD '
852 ELSE IF (info.eq.-9999) THEN
853 string='Could not build and Arnoldi factorization '
854 END IF
855!
856 RETURN
857 END SUBROUTINE iram_error
858
859 END MODULE roms_kernel_mod
subroutine edit_multifile(task)
subroutine memory
Definition memory.F:3
subroutine, public close_out
Definition close_io.F:175
subroutine, public close_file(ng, model, s, ncname, lupdate)
Definition close_io.F:43
subroutine, public close_inp(ng, model)
Definition close_io.F:92
subroutine, public def_gst(ng, model)
Definition def_gst.F:45
subroutine, public get_gst(ng, model)
Definition get_gst.F:38
subroutine, public inp_par(model)
Definition inp_par.F:56
subroutine, public roms_initialize_arrays
Definition mod_arrays.F:351
subroutine, public roms_allocate_arrays(allocate_vars)
Definition mod_arrays.F:114
type(t_io), dimension(:), allocatable adm
type(t_io), dimension(:), allocatable tlm
type(t_io), dimension(:), allocatable rst
integer stdout
character(len=256) sourcefile
integer, parameter io_nf90
Definition mod_ncparam.F:95
integer, parameter io_pio
Definition mod_ncparam.F:96
subroutine, public initialize_parallel
integer numthreads
integer, dimension(:), allocatable first_tile
integer mythread
logical master
integer, dimension(:), allocatable last_tile
integer ocn_comm_world
integer, dimension(:), allocatable nstr
Definition mod_param.F:646
integer, dimension(:), allocatable ntilex
Definition mod_param.F:685
integer, parameter iadm
Definition mod_param.F:665
integer, dimension(:), allocatable nsize
Definition mod_param.F:648
integer ngrids
Definition mod_param.F:113
integer, dimension(:), allocatable ntilee
Definition mod_param.F:686
integer, parameter itlm
Definition mod_param.F:663
integer, dimension(:), allocatable nend
Definition mod_param.F:647
integer, dimension(:), allocatable ntimes
integer, dimension(:), allocatable ntlm
logical lmultigst
integer, dimension(:), allocatable nconv
integer blowup
logical, dimension(:), allocatable lreadfrc
logical, dimension(:), allocatable lcycleadj
logical, dimension(:), allocatable lwrtadj
logical lrstgst
logical, dimension(:), allocatable lcycletlm
integer exit_flag
real(dp) ritz_tol
logical, dimension(:), allocatable lwrtper
integer nintervals
integer maxitergst
logical, dimension(:), allocatable lwrthis
logical, dimension(:), allocatable lwrttlm
logical, dimension(:), allocatable lwrtrst
integer nrun
integer, dimension(:), allocatable nadj
logical, dimension(:), allocatable lreadfwd
integer noerror
logical, dimension(:), allocatable lcyclerst
logical, dimension(:), allocatable lreadblk
integer msaup2
integer, dimension(8) iaup2
integer, dimension(:), allocatable ido
real(r8), dimension(:,:), allocatable sworkl
integer, dimension(:,:), allocatable ipntr
integer ncv
character(len=1) howmany
logical, dimension(5) laup2
logical lrvec
Definition mod_storage.F:95
integer msaupd
character(len=1) bmat
integer nev
type(t_storage), dimension(:), allocatable storage
Definition mod_storage.F:91
real(r8), dimension(:,:), allocatable norm
integer, dimension(:), allocatable info
subroutine, public allocate_storage
character(len=2) which
integer lworkl
integer, dimension(:,:), allocatable iparam
integer mseigt
integer msapps
integer ndigit
integer msgets
logical, dimension(:,:), allocatable pick
Definition mod_storage.F:96
real(r8) sigmar
integer logfil
real(r8), dimension(:), pointer sworkr
real(r8), dimension(:,:), allocatable rvaluer
integer msaitr
integer mseupd
subroutine, public propagator_so(runinterval, iter, state, ad_state)
subroutine, public roms_finalize
Definition ad_roms.h:283
subroutine, public roms_run(runinterval)
Definition ad_roms.h:239
subroutine, private iram_error(info, icall, string)
Definition afte_roms.h:862
subroutine, public roms_initialize(first, mpicomm)
Definition ad_roms.h:52
integer function, public stdout_unit(mymaster)
Definition stdout_mod.F:48
logical, save set_stdoutunit
Definition stdout_mod.F:41
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52
subroutine, public wrt_gst(ng, model)
Definition wrt_gst.F:39
subroutine, public wrt_rst(ng, tile)
Definition wrt_rst.F:63
subroutine r_norm2(ng, model, mstr, mend, evalue, evector, state, norm2)
Definition packing.F:175
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3
subroutine tl_initial(ng)
Definition tl_initial.F:4