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