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