ROMS
Loading...
Searching...
No Matches
pert_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 Tangent Linear and Adjoint Models Sanity Test: !
11! !
12! This driver is used to test the tangent linear model and adjoint !
13! models for ROMS. It is used to check whether or not both models !
14! are correct. This driver is usually run in "ensemble" mode where !
15! each member represents a perturbed interior point in the tangent !
16! linear and adjoint models. There is a maximum of Lm*Mm ensemble !
17! members. Each member is ran for few time steps. !
18! !
19! If each interior point is perturbed at one time, the resulting !
20! tangent linear (T) and adjoint (A) M-by-N matrices yield: !
21! !
22! T - tranpose(A) = 0 within round off !
23! !
24! That is, their inner product give a symmetric matrix. Here, M is !
25! the number of state points and N is the number of perturbations. !
26! This option is activated with the INNER_PRODUCT test. !
27! !
28! In realistic applications, it is awkward to perturb all interior !
29! points for each state variable. Alternatively, random check at a !
30! specified point is inexpensive. This option is activate with the !
31! SANITY_CHECK switch. The standard input "User" array is used to !
32! specify such point: !
33! !
34! INT(user(1)) => state tangent variable to perturb !
35! INT(user(2)) => state adjoint variable to perturb !
36! INT(user(3)) => I-index of tangent variable to perturb !
37! INT(user(4)) => I-index of adjoint variable to perturb !
38! INT(user(5)) => J-index of tangent variable to perturb !
39! INT(user(6)) => J-index of adjoint variable to perturb !
40! !
41! In 3D applications: !
42! !
43! INT(user(7)) => J-index of tangent variable to perturb !
44! INT(user(8)) => J-index of adjoint variable to perturb !
45! !
46! The subroutines in this driver control the initialization, time- !
47! stepping, and finalization of ROMS model following ESMF !
48! conventions: !
49! !
50! ROMS_initialize !
51! ROMS_run !
52! ROMS_finalize !
53! !
54!=======================================================================
55!
56 USE mod_param
57 USE mod_parallel
58 USE mod_arrays
59 USE mod_iounits
60 USE mod_ncparam
61 USE mod_ocean
62 USE mod_scalars
63 USE mod_stepping
64!
66#ifdef DISTRIBUTE
67 USE distribute_mod, ONLY : mp_collect
68#endif
69 USE inp_par_mod, ONLY : inp_par
70#ifdef MCT_LIB
71# ifdef ATM_COUPLING
72 USE mct_coupler_mod, ONLY : initialize_ocn2atm_coupling
73# endif
74# ifdef WAV_COUPLING
75 USE mct_coupler_mod, ONLY : initialize_ocn2wav_coupling
76# endif
77#endif
79 USE strings_mod, ONLY : founderror
80 USE wrt_rst_mod, ONLY : wrt_rst
81!
82 implicit none
83!
84 PUBLIC :: roms_initialize
85 PUBLIC :: roms_run
86 PUBLIC :: roms_finalize
87!
88 CONTAINS
89!
90 SUBROUTINE roms_initialize (first, mpiCOMM)
91!
92!=======================================================================
93! !
94! This routine allocates and initializes ROMS state variables !
95! and internal and external parameters. !
96! !
97!=======================================================================
98!
99! Imported variable declarations.
100!
101 logical, intent(inout) :: first
102!
103 integer, intent(in), optional :: mpiCOMM
104!
105! Local variable declarations.
106!
107 logical :: allocate_vars = .true.
108!
109#ifdef DISTRIBUTE
110 integer :: MyError, MySize
111#endif
112 integer :: chunk_size, ng, thread
113#ifdef _OPENMP
114 integer :: my_threadnum
115#endif
116!
117 character (len=*), parameter :: MyFile = &
118 & __FILE__//", ROMS_initialize"
119
120#ifdef DISTRIBUTE
121!
122!-----------------------------------------------------------------------
123! Set distribute-memory (mpi) world communictor.
124!-----------------------------------------------------------------------
125!
126 IF (PRESENT(mpicomm)) THEN
127 ocn_comm_world=mpicomm
128 ELSE
129 ocn_comm_world=mpi_comm_world
130 END IF
131 CALL mpi_comm_rank (ocn_comm_world, myrank, myerror)
132 CALL mpi_comm_size (ocn_comm_world, mysize, myerror)
133#endif
134!
135!-----------------------------------------------------------------------
136! On first pass, initialize model parameters a variables for all
137! nested/composed grids. Notice that the logical switch "first"
138! is used to allow multiple calls to this routine during ensemble
139! configurations.
140!-----------------------------------------------------------------------
141!
142 IF (first) THEN
143 first=.false.
144!
145! Initialize parallel control switches. These scalars switches are
146! independent from standard input parameters.
147!
149!
150! Set the ROMS standard output unit to write verbose execution info.
151! Notice that the default standard out unit in Fortran is 6.
152!
153! In some applications like coupling or disjointed mpi-communications,
154! it is advantageous to write standard output to a specific filename
155! instead of the default Fortran standard output unit 6. If that is
156! the case, it opens such formatted file for writing.
157!
158 IF (set_stdoutunit) THEN
160 set_stdoutunit=.false.
161 END IF
162!
163! Read in model tunable parameters from standard input. Allocate and
164! initialize variables in several modules after the number of nested
165! grids and dimension parameters are known.
166!
167 CALL inp_par (inlm)
168 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
169!
170! Set domain decomposition tile partition range. This range is
171! computed only once since the "first_tile" and "last_tile" values
172! are private for each parallel thread/node.
173!
174#if defined _OPENMP
175 mythread=my_threadnum()
176#elif defined DISTRIBUTE
178#else
179 mythread=0
180#endif
181 DO ng=1,ngrids
182 chunk_size=(ntilex(ng)*ntilee(ng)+numthreads-1)/numthreads
183 first_tile(ng)=mythread*chunk_size
184 last_tile(ng)=first_tile(ng)+chunk_size-1
185 END DO
186!
187! Initialize internal wall clocks. Notice that the timings does not
188! includes processing standard input because several parameters are
189! needed to allocate clock variables.
190!
191 IF (master) THEN
192 WRITE (stdout,10)
193 10 FORMAT (/,' Process Information:',/)
194 END IF
195!
196 DO ng=1,ngrids
197 DO thread=thread_range
198 CALL wclock_on (ng, inlm, 0, __line__, myfile)
199 END DO
200 END DO
201!
202! Allocate and initialize modules variables.
203!
204 CALL roms_allocate_arrays (allocate_vars)
206 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
207
208 END IF
209
210#if defined MCT_LIB && (defined ATM_COUPLING || defined WAV_COUPLING)
211!
212!-----------------------------------------------------------------------
213! Initialize coupling streams between model(s).
214!-----------------------------------------------------------------------
215!
216 DO ng=1,ngrids
217# ifdef ATM_COUPLING
218 CALL initialize_ocn2atm_coupling (ng, myrank)
219# endif
220# ifdef WAV_COUPLING
221 CALL initialize_ocn2wav_coupling (ng, myrank)
222# endif
223 END DO
224#endif
225!
226 RETURN
227 END SUBROUTINE roms_initialize
228!
229 SUBROUTINE roms_run (RunInterval)
230!
231!=======================================================================
232! !
233! This routine time-steps ROMS tangent linear and adjoint !
234! models. !
235! !
236!=======================================================================
237!
238! Imported variable declarations
239!
240 real(dp), intent(in) :: RunInterval ! seconds
241!
242! Local variable declarations.
243!
244 integer :: ng, tile
245#ifdef SANITY_CHECK
246 logical :: BOUNDED_AD, BOUNDED_TL, SAME_VAR
247# ifdef DISTRIBUTE
248 integer :: Istr, Iend, Jstr, Jend
249# endif
250 integer :: IperAD, JperAD, KperAD, ivarAD
251 integer :: IperTL, JperTL, KperTL, ivarTL
252 integer :: i
253!
254 real(r8) :: IniVal = 0.0_r8
255
256 real(r8), dimension(4,Ngrids) :: val
257#endif
258!
259 character (len=*), parameter :: MyFile = &
260 & __FILE__//", ROMS_run"
261!
262!=======================================================================
263! Run model for all nested grids, if any.
264!=======================================================================
265!
266#ifdef INNER_PRODUCT
267!
268! Set end of pertubation loos as the total number of state variables
269! interior points. Is not possible to run the inner product test over
270! nested grids. If nested grids, run each grid separately.
271!
272 IF (ngrids.gt.1) THEN
273 WRITE (stdout,10) 'Nested grids are not allowed, Ngrids = ', &
274 ngrids
275 stop
276 END IF
277 ng=1
278# ifdef SOLVE3D
279 erend=(lm(ng)-1)*mm(ng)+ &
280 & lm(ng)*(mm(ng)-1)+ &
281 & lm(ng)*mm(ng)+ &
282 & (lm(ng)-1)*mm(ng)*n(ng)+ &
283 & lm(ng)*(mm(ng)-1)*n(ng)+ &
284 & lm(ng)*mm(ng)*n(ng)*nat
285# else
286 erend=(lm(ng)-1)*mm(ng)+ &
287 & lm(ng)*(mm(ng)-1)+ &
288 & lm(ng)*mm(ng)
289# endif
290#endif
291#ifdef SANITY_CHECK
292!
293! Set tangent and adjoint variable and random point to perturb.
294!
295 ivartl=int(user(1))
296 ivarad=int(user(2))
297 ipertl=int(user(3))
298 iperad=int(user(4))
299 jpertl=int(user(5))
300 jperad=int(user(6))
301# ifdef SOLVE3D
302 kpertl=int(user(7))
303 kperad=int(user(8))
304 same_var=(ivartl.eq.ivarad).and. &
305 & (ipertl.eq.iperad).and. &
306 & (jpertl.eq.jperad).and. &
307 & (kpertl.eq.kperad)
308# else
309 same_var=(ivartl.eq.ivarad).and. &
310 & (ipertl.eq.iperad).and. &
311 & (jpertl.eq.jperad)
312# endif
313#endif
314!
315! Set relevant IO switches.
316!
317 DO ng=1,ngrids
318 IF (ntlm(ng).gt.0) ldeftlm(ng)=.true.
319 IF (nadj(ng).gt.0) ldefadj(ng)=.true.
320 lreadfwd(ng)=.true.
321 END DO
322 lstiffness=.false.
323
324#ifdef FORWARD_FLUXES
325!
326! Set the BLK structure to contain the nonlinear model surface fluxes
327! needed by the tangent linear and adjoint models. Also, set switches
328! to process that structure in routine "check_multifile". Notice that
329! it is possible to split the solution into multiple NetCDF files to
330! reduce their size.
331!
332! The switch LreadFRC is deactivated because all the atmospheric
333! forcing, including shortwave radiation, is read from the NLM
334! surface fluxes or is assigned during ESM coupling. Such fluxes
335! are available from the QCK structure. There is no need for reading
336! and processing from the FRC structure input forcing-files.
337
338 CALL edit_multifile ('QCK2BLK')
339 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
340 DO ng=1,ngrids
341 lreadblk(ng)=.true.
342 lreadfrc(ng)=.false.
343 END DO
344#endif
345!
346!=======================================================================
347! Perturbation loop.
348!=======================================================================
349!
350 pert_loop : DO nrun=erstr,erend
351!
352!-----------------------------------------------------------------------
353! Time step tangent linear model.
354!-----------------------------------------------------------------------
355!
356 tlmodel=.true.
357 admodel=.false.
358 DO ng=1,ngrids
359 CALL tl_initial (ng)
360 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
361
362 IF (master) THEN
363 WRITE (stdout,10) 'TL', ng, ntstart(ng), ntend(ng)
364 END IF
365 END DO
366!
367#ifdef SOLVE3D
368 CALL tl_main3d (runinterval)
369#else
370 CALL tl_main2d (runinterval)
371#endif
372 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
373
374#ifdef SANITY_CHECK
375!
376! Get check value for tangent linear state variable.
377!
378 DO ng=1,ngrids
379# ifdef DISTRIBUTE
380 istr=bounds(ng)%Istr(myrank)
381 iend=bounds(ng)%Iend(myrank)
382 jstr=bounds(ng)%Jstr(myrank)
383 jend=bounds(ng)%Jend(myrank)
384 bounded_tl=((istr.le.ipertl).and.(ipertl.le.iend)).and. &
385 & ((jstr.le.jpertl).and.(jpertl.le.jend))
386 bounded_ad=((istr.le.iperad).and.(iperad.le.iend)).and. &
387 & ((jstr.le.jperad).and.(jperad.le.jend))
388# else
389 bounded_tl=.true.
390 bounded_ad=.true.
391# endif
392 DO i=1,4
393 val(i,ng)=inival
394 END DO
395 IF (bounded_tl) THEN
396# ifdef SOLVE3D
397 IF (ivartl.eq.isubar) THEN
398 val(1,ng)=ocean(ng)%tl_ubar(ipertl,jpertl,knew(ng))
399 ELSE IF (ivartl.eq.isvbar) THEN
400 val(1,ng)=ocean(ng)%tl_vbar(ipertl,jpertl,knew(ng))
401 ELSE IF (ivartl.eq.isfsur) THEN
402 val(1,ng)=ocean(ng)%tl_zeta(ipertl,jpertl,knew(ng))
403 ELSE IF (ivartl.eq.isuvel) THEN
404 val(1,ng)=ocean(ng)%tl_u(ipertl,jpertl,kpertl,nstp(ng))
405 ELSE IF (ivartl.eq.isvvel) THEN
406 val(1,ng)=ocean(ng)%tl_v(ipertl,jpertl,kpertl,nstp(ng))
407 ELSE
408 DO i=1,nt(ng)
409 IF (ivartl.eq.istvar(i)) THEN
410 val(1,ng)=ocean(ng)%tl_t(ipertl,jpertl,kpertl, &
411 & nstp(ng),i)
412 END IF
413 END DO
414 END IF
415# else
416 IF (ivartl.eq.isubar) THEN
417 val(1,ng)=ocean(ng)%tl_ubar(ipertl,jpertl,knew(ng))
418 ELSE IF (ivartl.eq.isvbar) THEN
419 val(1,ng)=ocean(ng)%tl_vbar(ipertl,jpertl,knew(ng))
420 ELSE IF (ivartl.eq.isfsur) THEN
421 val(1,ng)=ocean(ng)%tl_zeta(ipertl,jpertl,knew(ng))
422 END IF
423# endif
424 END IF
425
426 IF (bounded_ad) THEN
427# ifdef SOLVE3D
428 IF (ivarad.eq.isubar) THEN
429 val(3,ng)=ocean(ng)%tl_ubar(iperad,jperad,knew(ng))
430 ELSE IF (ivarad.eq.isvbar) THEN
431 val(3,ng)=ocean(ng)%tl_vbar(iperad,jperad,knew(ng))
432 ELSE IF (ivarad.eq.isfsur) THEN
433 val(3,ng)=ocean(ng)%tl_zeta(iperad,jperad,knew(ng))
434 ELSE IF (ivarad.eq.isuvel) THEN
435 val(3,ng)=ocean(ng)%tl_u(iperad,jperad,kperad,nstp(ng))
436 ELSE IF (ivarad.eq.isvvel) THEN
437 val(3,ng)=ocean(ng)%tl_v(iperad,jperad,kperad,nstp(ng))
438 ELSE
439 DO i=1,nt(ng)
440 IF (ivarad.eq.istvar(i)) THEN
441 val(3,ng)=ocean(ng)%tl_t(iperad,jperad,kperad, &
442 & nstp(ng),i)
443 END IF
444 END DO
445 END IF
446# else
447 IF (ivarad.eq.isubar) THEN
448 val(3,ng)=ocean(ng)%tl_ubar(iperad,jperad,knew(ng))
449 ELSE IF (ivarad.eq.isvbar) THEN
450 val(3,ng)=ocean(ng)%tl_vbar(iperad,jperad,knew(ng))
451 ELSE IF (ivarad.eq.isfsur) THEN
452 val(3,ng)=ocean(ng)%tl_zeta(iperad,jperad,knew(ng))
453 END IF
454# endif
455 END IF
456 END DO
457#endif
458!
459! Clear all model state arrays arrays.
460!
461 DO ng=1,ngrids
462 DO tile=first_tile(ng),last_tile(ng),+1
463 CALL initialize_ocean (ng, tile, 0)
464 END DO
465 END DO
466!
467!-----------------------------------------------------------------------
468! Time step adjoint model backwards.
469!-----------------------------------------------------------------------
470!
471 tlmodel=.false.
472 admodel=.true.
473 DO ng=1,ngrids
474 CALL ad_initial (ng)
475 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
476
477 IF (master) THEN
478 WRITE (stdout,10) 'AD', ng, ntstart(ng), ntend(ng)
479 END IF
480 END DO
481!
482#ifdef SOLVE3D
483 CALL ad_main3d (runinterval)
484#else
485 CALL ad_main2d (runinterval)
486#endif
487 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
488
489#ifdef SANITY_CHECK
490!
491! Get check value for adjoint state variable.
492!
493 DO ng=1,ngrids
494# ifdef DISTRIBUTE
495 istr=bounds(ng)%Istr(myrank)
496 iend=bounds(ng)%Iend(myrank)
497 jstr=bounds(ng)%Jstr(myrank)
498 jend=bounds(ng)%Jend(myrank)
499 bounded_ad=((istr.le.iperad).and.(iperad.le.iend)).and. &
500 & ((jstr.le.jperad).and.(jperad.le.jend))
501 bounded_tl=((istr.le.ipertl).and.(ipertl.le.iend)).and. &
502 & ((jstr.le.jpertl).and.(jpertl.le.jend))
503# else
504 bounded_ad=.true.
505 bounded_tl=.true.
506# endif
507 IF (bounded_ad) THEN
508# ifdef SOLVE3D
509 IF (ivarad.eq.isubar) THEN
510 val(2,ng)=ocean(ng)%ad_ubar(iperad,jperad,kstp(ng))
511 ELSE IF (ivarad.eq.isvbar) THEN
512 val(2,ng)=ocean(ng)%ad_vbar(iperad,jperad,kstp(ng))
513 ELSE IF (ivarad.eq.isfsur) THEN
514 val(2,ng)=ocean(ng)%ad_zeta(iperad,jperad,kstp(ng))
515 ELSE IF (ivarad.eq.isuvel) THEN
516 val(2,ng)=ocean(ng)%ad_u(iperad,jperad,kperad,nstp(ng))
517 ELSE IF (ivarad.eq.isvvel) THEN
518 val(2,ng)=ocean(ng)%ad_v(iperad,jperad,kperad,nstp(ng))
519 ELSE
520 DO i=1,nt(ng)
521 IF (ivarad.eq.istvar(i)) THEN
522 val(2,ng)=ocean(ng)%ad_t(iperad,jperad,kperad, &
523 & nstp(ng),i)
524 END IF
525 END DO
526 END IF
527# else
528 IF (ivarad.eq.isubar) THEN
529 val(2,ng)=ocean(ng)%ad_ubar(iperad,jperad,kstp(ng))
530 ELSE IF (ivarad.eq.isvbar) THEN
531 val(2,ng)=ocean(ng)%ad_vbar(iperad,jperad,kstp(ng))
532 ELSE IF (ivarad.eq.isfsur) THEN
533 val(2,ng)=ocean(ng)%ad_zeta(iperad,jperad,kstp(ng))
534 END IF
535# endif
536 END IF
537
538 IF (bounded_tl) THEN
539# ifdef SOLVE3D
540 IF (ivartl.eq.isubar) THEN
541 val(4,ng)=ocean(ng)%ad_ubar(ipertl,jpertl,kstp(ng))
542 ELSE IF (ivartl.eq.isvbar) THEN
543 val(4,ng)=ocean(ng)%ad_vbar(ipertl,jpertl,kstp(ng))
544 ELSE IF (ivartl.eq.isfsur) THEN
545 val(4,ng)=ocean(ng)%ad_zeta(ipertl,jpertl,kstp(ng))
546 ELSE IF (ivartl.eq.isuvel) THEN
547 val(4,ng)=ocean(ng)%ad_u(ipertl,jpertl,kpertl,nstp(ng))
548 ELSE IF (ivartl.eq.isvvel) THEN
549 val(4,ng)=ocean(ng)%ad_v(ipertl,jpertl,kpertl,nstp(ng))
550 ELSE
551 DO i=1,nt(ng)
552 IF (ivartl.eq.istvar(i)) THEN
553 val(4,ng)=ocean(ng)%ad_t(ipertl,jpertl,kpertl, &
554 & nstp(ng),i)
555 END IF
556 END DO
557 END IF
558# else
559 IF (ivartl.eq.isubar) THEN
560 val(4,ng)=ocean(ng)%ad_ubar(ipertl,jpertl,kstp(ng))
561 ELSE IF (ivartl.eq.isvbar) THEN
562 val(4,ng)=ocean(ng)%ad_vbar(ipertl,jpertl,kstp(ng))
563 ELSE IF (ivartl.eq.isfsur) THEN
564 val(4,ng)=ocean(ng)%ad_zeta(ipertl,jpertl,kstp(ng))
565 END IF
566# endif
567 END IF
568!
569! Report sanity check values.
570!
571# ifdef DISTRIBUTE
572 CALL mp_collect (ng, itlm, 4, inival, val(:,ng))
573# endif
574 IF (master) THEN
575 IF (same_var) THEN
576 WRITE (stdout,20) 'Perturbing', &
577 & trim(vname(1,idsvar(ivartl)))
578 IF (ivartl.le.3) THEN
579 WRITE (stdout,30) 'Tangent: ', val(1,ng), &
580 & ipertl,jpertl
581 WRITE (stdout,30) 'Adjoint: ', val(2,ng), &
582 & iperad, jperad
583 WRITE (stdout,30) 'Difference: ', val(2,ng)-val(1,ng), &
584 & ipertl,jpertl
585 ELSE
586 WRITE (stdout,40) 'Tangent: ', val(1,ng), &
587 & ipertl,jpertl,kpertl
588 WRITE (stdout,40) 'Adjoint: ', val(2,ng), &
589 & iperad,jperad,kperad
590 WRITE (stdout,40) 'Difference: ', val(2,ng)-val(1,ng), &
591 & ipertl,jpertl,kpertl
592 END IF
593 ELSE
594 IF (ivartl.le.3) THEN
595 WRITE (stdout,50) 'Tangent, Perturbing: ', &
596 & trim(vname(1,idsvar(ivartl))), &
597 & ipertl,jpertl
598 WRITE (stdout,60) trim(vname(1,idsvar(ivartl))), &
599 & val(1,ng),ipertl,jpertl
600 ELSE
601 WRITE (stdout,70) 'Tangent, Perturbing: ', &
602 & trim(vname(1,idsvar(ivartl))), &
603 & ipertl,jpertl,kpertl
604 WRITE (stdout,80) trim(vname(1,idsvar(ivartl))), &
605 & val(1,ng),ipertl,jpertl,kpertl
606 END IF
607 IF (ivarad.le.3) THEN
608 WRITE (stdout,60) trim(vname(1,idsvar(ivarad))), &
609 & val(3,ng),iperad,jperad
610 ELSE
611 WRITE (stdout,80) trim(vname(1,idsvar(ivarad))), &
612 & val(3,ng),iperad,jperad,kperad
613 END IF
614!
615 IF (ivarad.le.3) THEN
616 WRITE (stdout,50) 'Adjoint, Perturbing: ', &
617 & trim(vname(1,idsvar(ivarad))), &
618 & iperad,jperad
619 WRITE (stdout,60) trim(vname(1,idsvar(ivarad))), &
620 & val(2,ng),iperad,jperad
621 ELSE
622 WRITE (stdout,70) 'Adjoint, Perturbing: ', &
623 & trim(vname(1,idsvar(ivarad))), &
624 & iperad,jperad,kperad
625 WRITE (stdout,80) trim(vname(1,idsvar(ivarad))), &
626 & val(2,ng),iperad,jperad,kperad
627 END IF
628 IF (ivartl.le.3) THEN
629 WRITE (stdout,60) trim(vname(1,idsvar(ivartl))), &
630 & val(4,ng),ipertl,jpertl
631 ELSE
632 WRITE (stdout,80) trim(vname(1,idsvar(ivartl))), &
633 & val(4,ng),ipertl,jpertl,kpertl
634 END IF
635 WRITE (stdout,90) val(3,ng)-val(4,ng)
636 END IF
637 END IF
638 END DO
639#endif
640!
641! Clear model state arrays arrays.
642!
643 DO ng=1,ngrids
644 DO tile=first_tile(ng),last_tile(ng),+1
645 CALL initialize_ocean (ng, tile, 0)
646 END DO
647 END DO
648!
649! Close current forward NetCDF file.
650!
651 sourcefile=myfile
652 DO ng=1,ngrids
653 CALL close_file (ng, itlm, fwd(ng), fwd(ng)%name)
654 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
655 END DO
656
657 END DO pert_loop
658!
659 10 FORMAT (/,1x,a,1x,'ROMS: started time-stepping:', &
660 & ' (Grid: ',i2.2,' TimeSteps: ',i8.8,' - ',i8.8,')',/)
661#ifdef SANITY_CHECK
662 20 FORMAT (/,' Sanity Check - ',a,' variable: ',a,t60)
663 30 FORMAT (' Sanity Check - ', a, 1p,e19.12, &
664 & 3x, 'at (i,j) ',2i4)
665 40 FORMAT (' Sanity Check - ', a, 1p,e19.12, &
666 & 3x, 'at (i,j,k) ',3i4)
667 50 FORMAT (/,' Sanity Check - ',a, a, t52, 'at (i,j) ', 2i4)
668 60 FORMAT (' Sanity Check - ', a, ' =', t30, 1p,e19.12, &
669 & t52, 'at (i,j) ',2i4)
670 70 FORMAT (/,' Sanity Check - ',a, a, t52, 'at (i,j,k) ', 3i4)
671 80 FORMAT (' Sanity Check - ', a, ' =', t30, 1p,e19.12, &
672 & t52, 'at (i,j,k) ',3i4)
673 90 FORMAT (/,' Sanity Check - Difference = ', 1p,e19.12)
674#endif
675!
676 RETURN
677 END SUBROUTINE roms_run
678!
679 SUBROUTINE roms_finalize
680!
681!=======================================================================
682! !
683! This routine terminates ROMS tangent linear and adjoint !
684! models execution. !
685! !
686!=======================================================================
687!
688! Local variable declarations.
689!
690 integer :: Fcount, ng, thread
691!
692 character (len=*), parameter :: MyFile = &
693 & __FILE__//", ROMS_finalize"
694!
695!-----------------------------------------------------------------------
696! If blowing-up, save latest model state into RESTART NetCDF file.
697!-----------------------------------------------------------------------
698!
699! If cycling restart records, write solution into the next record.
700!
701 IF (exit_flag.eq.1) THEN
702 DO ng=1,ngrids
703 IF (lwrtrst(ng)) THEN
704 IF (master) WRITE (stdout,10)
705 10 FORMAT (/,' Blowing-up: Saving latest model state into ', &
706 & ' RESTART file',/)
707 fcount=rst(ng)%load
708 IF (lcyclerst(ng).and.(rst(ng)%Nrec(fcount).ge.2)) THEN
709 rst(ng)%Rindex=2
710 lcyclerst(ng)=.false.
711 END IF
714#ifdef DISTRIBUTE
715 CALL wrt_rst (ng, myrank)
716#else
717 CALL wrt_rst (ng, -1)
718#endif
719 END IF
720 END DO
721 END IF
722!
723!-----------------------------------------------------------------------
724! Stop model and time profiling clocks, report memory requirements, and
725! close output NetCDF files.
726!-----------------------------------------------------------------------
727!
728! Stop time clocks.
729!
730 IF (master) THEN
731 WRITE (stdout,20)
732 20 FORMAT (/,'Elapsed wall CPU time for each process (seconds):',/)
733 END IF
734!
735 DO ng=1,ngrids
736 DO thread=thread_range
737 CALL wclock_off (ng, inlm, 0, __line__, myfile)
738 END DO
739 END DO
740!
741! Report dynamic memory and automatic memory requirements.
742!
743 CALL memory
744!
745! Close IO files.
746!
747 DO ng=1,ngrids
748 CALL close_inp (ng, inlm)
749 END DO
750 CALL close_out
751!
752 RETURN
753 END SUBROUTINE roms_finalize
754
755 END MODULE roms_kernel_mod
subroutine ad_initial(ng)
Definition ad_initial.F:4
subroutine ad_main2d
Definition ad_main2d.F:586
subroutine ad_main3d(runinterval)
Definition ad_main3d.F:4
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 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 fwd
type(t_io), dimension(:), allocatable rst
integer stdout
character(len=256) sourcefile
integer isvvel
integer isvbar
integer, dimension(:), allocatable istvar
integer isuvel
integer isfsur
character(len=maxlen), dimension(6, 0:nv) vname
integer, dimension(:), allocatable idsvar
integer isubar
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
subroutine, public initialize_ocean(ng, tile, model)
Definition mod_ocean.F:1526
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 nat
Definition mod_param.F:499
integer, parameter inlm
Definition mod_param.F:662
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
integer, dimension(:), allocatable ntilex
Definition mod_param.F:685
integer, dimension(:), allocatable lm
Definition mod_param.F:455
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 nt
Definition mod_param.F:489
integer, dimension(:), allocatable mm
Definition mod_param.F:456
integer, dimension(:), allocatable ntlm
logical tlmodel
logical lstiffness
integer blowup
real(r8), dimension(:), allocatable user
integer erend
logical, dimension(:), allocatable lreadfrc
logical admodel
logical, dimension(:), allocatable ldefadj
integer, dimension(:), allocatable ntend
integer exit_flag
integer erstr
logical, dimension(:), allocatable lwrtrst
logical, dimension(:), allocatable ldeftlm
integer, dimension(:), allocatable ntstart
integer nrun
integer, dimension(:), allocatable nadj
logical, dimension(:), allocatable lreadfwd
integer noerror
logical, dimension(:), allocatable lcyclerst
logical, dimension(:), allocatable lreadblk
integer, dimension(:), allocatable kstp
integer, dimension(:), allocatable knew
integer, dimension(:), allocatable nstp
subroutine, public roms_finalize
Definition ad_roms.h:283
subroutine, public roms_run(runinterval)
Definition ad_roms.h:239
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_rst(ng, tile)
Definition wrt_rst.F:63
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
subroutine tl_main2d
Definition tl_main2d.F:429
subroutine tl_main3d(runinterval)
Definition tl_main3d.F:4