ROMS
Loading...
Searching...
No Matches
tlcheck_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 Model Linearization Test: !
11! !
12! This driver is used to check the "linearization" of the tangent !
13! linear model using a structure similar to the gradient test. !
14! !
15! The subroutines in this driver control the initialization, time- !
16! stepping, and finalization of ROMS model following ESMF !
17! conventions: !
18! !
19! ROMS_initialize !
20! ROMS_run !
21! ROMS_finalize !
22! !
23!=======================================================================
24!
25 USE mod_param
26 USE mod_parallel
27 USE mod_arrays
28 USE mod_fourdvar
29 USE mod_iounits
30 USE mod_ncparam
31 USE mod_scalars
32 USE mod_stepping
33!
35 USE def_mod_mod, ONLY : def_mod
37 USE get_state_mod, ONLY : get_state
38 USE inp_par_mod, ONLY : inp_par
39#ifdef MCT_LIB
40# ifdef ATM_COUPLING
41 USE mct_coupler_mod, ONLY : initialize_ocn2atm_coupling
42# endif
43# ifdef WAV_COUPLING
44 USE mct_coupler_mod, ONLY : initialize_ocn2wav_coupling
45# endif
46#endif
48 USE strings_mod, ONLY : founderror
49 USE wrt_rst_mod, ONLY : wrt_rst
50!
51 implicit none
52!
53 PUBLIC :: roms_initialize
54 PUBLIC :: roms_run
55 PUBLIC :: roms_finalize
56!
57 CONTAINS
58!
59 SUBROUTINE roms_initialize (first, mpiCOMM)
60!
61!=======================================================================
62! !
63! This routine allocates and initializes ROMS state variables !
64! and internal and external parameters. !
65! !
66!=======================================================================
67!
68! Imported variable declarations.
69!
70 logical, intent(inout) :: first
71!
72 integer, intent(in), optional :: mpiCOMM
73!
74! Local variable declarations.
75!
76 logical :: allocate_vars = .true.
77!
78#ifdef DISTRIBUTE
79 integer :: MyError, MySize
80#endif
81 integer :: chunk_size, ng, thread
82!
83 character (len=*), parameter :: MyFile = &
84 & __FILE__//", ROMS_initialize"
85
86#ifdef DISTRIBUTE
87!
88!-----------------------------------------------------------------------
89! Set distribute-memory (mpi) world communictor.
90!-----------------------------------------------------------------------
91!
92 IF (PRESENT(mpicomm)) THEN
93 ocn_comm_world=mpicomm
94 ELSE
95 ocn_comm_world=mpi_comm_world
96 END IF
97 CALL mpi_comm_rank (ocn_comm_world, myrank, myerror)
98 CALL mpi_comm_size (ocn_comm_world, mysize, myerror)
99#endif
100!
101!-----------------------------------------------------------------------
102! On first pass, initialize model parameters a variables for all
103! nested/composed grids. Notice that the logical switch "first"
104! is used to allow multiple calls to this routine during ensemble
105! configurations.
106!-----------------------------------------------------------------------
107!
108 IF (first) THEN
109 first=.false.
110!
111! Initialize parallel control switches. These scalars switches are
112! independent from standard input parameters.
113!
115!
116! Set the ROMS standard output unit to write verbose execution info.
117! Notice that the default standard out unit in Fortran is 6.
118!
119! In some applications like coupling or disjointed mpi-communications,
120! it is advantageous to write standard output to a specific filename
121! instead of the default Fortran standard output unit 6. If that is
122! the case, it opens such formatted file for writing.
123!
124 IF (set_stdoutunit) THEN
126 set_stdoutunit=.false.
127 END IF
128!
129! Read in model tunable parameters from standard input. Allocate and
130! initialize variables in several modules after the number of nested
131! grids and dimension parameters are known.
132!
133 CALL inp_par (inlm)
134 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
135!
136! Set domain decomposition tile partition range. This range is
137! computed only once since the "first_tile" and "last_tile" values
138! are private for each parallel thread/node.
139!
140#if defined _OPENMP
142#elif defined DISTRIBUTE
144#else
145 mythread=0
146#endif
147 DO ng=1,ngrids
148 chunk_size=(ntilex(ng)*ntilee(ng)+numthreads-1)/numthreads
149 first_tile(ng)=mythread*chunk_size
150 last_tile(ng)=first_tile(ng)+chunk_size-1
151 END DO
152!
153! Initialize internal wall clocks. Notice that the timings does not
154! includes processing standard input because several parameters are
155! needed to allocate clock variables.
156!
157 IF (master) THEN
158 WRITE (stdout,10)
159 10 FORMAT (/,' Process Information:',/)
160 END IF
161!
162 DO ng=1,ngrids
163 DO thread=thread_range
164 CALL wclock_on (ng, inlm, 0, __line__, myfile)
165 END DO
166 END DO
167!
168! Allocate and initialize modules variables.
169!
170 CALL roms_allocate_arrays (allocate_vars)
172 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
173
174 END IF
175
176#if defined MCT_LIB && (defined ATM_COUPLING || defined WAV_COUPLING)
177!
178!-----------------------------------------------------------------------
179! Initialize coupling streams between model(s).
180!-----------------------------------------------------------------------
181!
182 DO ng=1,ngrids
183# ifdef ATM_COUPLING
184 CALL initialize_ocn2atm_coupling (ng, myrank)
185# endif
186# ifdef WAV_COUPLING
187 CALL initialize_ocn2wav_coupling (ng, myrank)
188# endif
189 END DO
190#endif
191!
192 RETURN
193 END SUBROUTINE roms_initialize
194!
195 SUBROUTINE roms_run (RunInterval)
196!
197!=======================================================================
198! !
199! This routine time-steps ROMS nonlinear, tangent linear and !
200! adjoint models. !
201! !
202!=======================================================================
203!
204! Imported variable declarations
205!
206 real(dp), intent(in) :: RunInterval ! seconds
207!
208! Local variable declarations.
209!
210 integer :: IniRec, Ipass, i, ng, status
211 integer :: tile
212!
213 real(r8) :: gp, hp, p
214!
215 character (len=*), parameter :: MyFile = &
216 & __FILE__//", ROMS_run"
217!
218!-----------------------------------------------------------------------
219! Run tangent linear model test.
220!-----------------------------------------------------------------------
221!
222! Currently, the tangent linear model test cannot be run over
223! nested grids.
224!
225 IF (ngrids.gt.1) THEN
226 WRITE (stdout,10) 'Nested grids are not allowed, Ngrids = ', &
227 ngrids
228 stop
229 END IF
230!
231! Initialize relevant parameters.
232!
233 DO ng=1,ngrids
234 lold(ng)=1
235 lnew(ng)=1
236 ntlm(ng)=nhis(ng) ! to allow IO comparison
237#ifdef FORWARD_FLUXES
238 lreadblk(ng)=.false.
239#endif
240 lreadfwd(ng)=.false.
241 END DO
242 nrun=1
243 ipass=1
244 erstr=1
245 erend=maxval(nstatevar)+1
246 inirec=1
247 ig1count=0
248 ig2count=0
249 lcycletlm=.false.
250!
251! Initialize nonlinear model with first guess initial conditions.
252!
253 CALL initial
254 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
255!
256! Define output 4D-Var NetCDF file (DAV struc) containing all
257! processed data at observation locations.
258!
259 DO ng=1,ngrids
260 ldefmod(ng)=.true.
261 CALL def_mod (ng)
262 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
263 END DO
264!
265! Run nonlinear model. Extract and store nonlinear model values at
266! observation locations.
267!
268 DO ng=1,ngrids
269 IF (master) THEN
270 WRITE (stdout,20) 'NL', ng, ntstart(ng), ntend(ng)
271 END IF
272 wrtnlmod(ng)=.true.
273 wrttlmod(ng)=.false.
274 END DO
275!
276#ifdef SOLVE3D
277 CALL main3d (runinterval)
278#else
279 CALL main2d (runinterval)
280#endif
281 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
282
283 DO ng=1,ngrids
284 wrtnlmod(ng)=.false.
285 wrttlmod(ng)=.true.
286 END DO
287!
288! Set structure for the nonlinear forward trajectory to be processed
289! by the tangent linear and adjoint models. Also, set switches to
290! process the FWD structure in routine "check_multifile". Notice that
291! it is possible to split solution into multiple NetCDF files to reduce
292! their size.
293!
294 CALL edit_multifile ('HIS2FWD')
295 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
296 DO ng=1,ngrids
297 lreadfwd(ng)=.true.
298 END DO
299
300#ifdef FORWARD_FLUXES
301!
302! Set structure for the nonlinear surface fluxes to be processed by
303! by the tangent linear and adjoint models. Also, set switches to
304! process the BLK structure in routine "check_multifile". Notice that
305! it is possible to split solution into multiple NetCDF files to reduce
306! their size.
307!
308! The switch LreadFRC is deactivated because all the atmospheric
309! forcing, including shortwave radiation, is read from the NLM bulk
310! flux formulation or is assigned during ESM coupling (BLK structure).
311! There is no need for processing forcing from the FRC structure files.
312!
313 CALL edit_multifile ('QCK2BLK')
314 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
315 DO ng=1,ngrids
316 lreadblk(ng)=.true.
317 lreadfrc(ng)=.false.
318 END DO
319#endif
320!
321! Save and Report cost function between nonlinear model and
322! observations.
323!
324 DO ng=1,ngrids
325 DO i=0,nobsvar(ng)
326 fourdvar(ng)%CostFunOld(i)=fourdvar(ng)%CostFun(i)
327 END DO
328 IF (master) THEN
329 WRITE (stdout,40) fourdvar(ng)%CostFunOld(0)
330 DO i=1,nobsvar(ng)
331 IF (fourdvar(ng)%CostFunOld(i).gt.0.0_r8) THEN
332 IF (i.eq.1) THEN
333 WRITE (stdout,50) fourdvar(ng)%CostFunOld(i), &
334 & trim(obsname(i))
335 ELSE
336 WRITE (stdout,60) fourdvar(ng)%CostFunOld(i), &
337 & trim(obsname(i))
338 END IF
339 END IF
340 END DO
341 END IF
342 END DO
343!
344! Initialize the adjoint model from rest.
345!
346 DO ng=1,ngrids
347!! CALL ad_initial (ng, .TRUE.)
348 CALL ad_initial (ng)
349 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
350 END DO
351!
352! Time-step adjoint model: Compute model state gradient, GRAD(J).
353! Force the adjoint model with the adjoint misfit between nonlinear
354! model and observations.
355!
356 DO ng=1,ngrids
357 IF (master) THEN
358 WRITE (stdout,20) 'AD', ng, ntstart(ng), ntend(ng)
359 END IF
360 END DO
361!
362#ifdef SOLVE3D
363 CALL ad_main3d (runinterval)
364#else
365 CALL ad_main2d (runinterval)
366#endif
367 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
368!
369!-----------------------------------------------------------------------
370! Perturb each tangent linear state variable using the steepest decent
371! direction by grad(J). If ndefTLM < 0, suppress IO for both nonlinear
372! and tangent linear models in the outer and inner loops.
373!-----------------------------------------------------------------------
374!
375! Load adjoint solution.
376!
377 DO ng=1,ngrids
378 CALL get_state (ng, iadm, 3, adm(ng), adm(ng)%Rindex, lnew(ng))
379 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
380 END DO
381!
382! Compute adjoint solution dot product for scaling purposes.
383!
384 DO ng=1,ngrids
385 DO tile=first_tile(ng),last_tile(ng),+1
386 CALL ad_dotproduct (ng, tile, lnew(ng))
387 END DO
388 END DO
389
390#ifdef SOLVE3D
391!
392! OUTER LOOP: First, Perturb all state variables (zeta, u, v, t) at
393! ========== once. Then, perturb one state variable at the time.
394! Notice, that ubar and vbar are not perturbed. They are
395! computed by vertically integarting u and v.
396!
397#else
398!
399! OUTER LOOP: First, perturb all state variables (zeta, ubar, vbar) at
400! ========== once. Then, perturb one state variable at the time.
401!
402#endif
403 outer_loop : DO outer=erstr,erend
404
405 DO ng=1,ngrids
406 CALL get_state (ng, inlm, 1, fwd(ng), inirec, lnew(ng))
407 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
408 END DO
409!
410! INNER LOOP: scale perturbation amplitude by selecting "p" scalar,
411! ========== such that:
412! p = 10 ** REAL(-inner,r8)
413!
414 inner_loop : DO inner=1,ninner
415!
416! Add a perturbation to the nonlinear state initial conditions
417! according with the outer and inner loop iterations. The added
418! term is a function of the steepest descent direction defined
419! by grad(J) times the perturbation amplitude "p".
420!
421!! CALL initial (.FALSE.)
422 CALL initial
423 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
424
425 DO ng=1,ngrids
426 WRITE (his(ng)%name,70) trim(his(ng)%base), nrun
427
428 IF (ndeftlm(ng).lt.0) THEN
429 ldefhis(ng)=.false. ! suppress IO
430 lwrthis(ng)=.false.
431 ELSE
432 ldefhis(ng)=.true.
433 END IF
434 wrtnlmod(ng)=.true.
435 wrttlmod(ng)=.false.
436 END DO
437!
438! Time-step nonlinear model: compute perturbed nonlinear state.
439!
440 DO ng=1,ngrids
441 IF (master) THEN
442 WRITE (stdout,20) 'NL', ng, ntstart(ng), ntend(ng)
443 END IF
444 END DO
445!
446#ifdef SOLVE3D
447 CALL main3d (runinterval)
448#else
449 CALL main2d (runinterval)
450#endif
451 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
452!
453! Get current nonlinear model trajectory.
454!
455 DO ng=1,ngrids
456 fwd(ng)%name=trim(his(ng)%head)//'.nc'
457 CALL get_state (ng, inlm, 1, fwd(ng), inirec, lnew(ng))
458 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
459 END DO
460!
461! Initialize tangent linear with the steepest decent direction
462! (adjoint state, GRAD(J)) times the perturbation amplitude "p".
463!
464 DO ng=1,ngrids
465!! CALL tl_initial (ng, .FALSE.)
466 CALL tl_initial (ng)
467 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
468
469 WRITE (tlm(ng)%name,70) trim(tlm(ng)%base), nrun
470
471 IF (ndeftlm(ng).lt.0) THEN
472 ldeftlm(ng)=.false. ! suppress IO
473 lwrttlm(ng)=.false.
474 ELSE
475 ldeftlm(ng)=.true.
476 END IF
477 END DO
478!
479! Time-step tangent linear model: Compute misfit cost function
480! between model (nonlinear + tangent linear) and observations.
481!
482 DO ng=1,ngrids
483 IF (master) THEN
484 WRITE (stdout,20) 'TL', ng, ntstart(ng), ntend(ng)
485 END IF
486 END DO
487!
488#ifdef SOLVE3D
489 CALL tl_main3d (runinterval)
490#else
491 CALL tl_main2d (runinterval)
492#endif
493 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
494!
495! Close current tangent linear model history file.
496!
497 sourcefile=myfile
498 DO ng=1,ngrids
499 CALL close_file (ng, itlm, tlm(ng))
500 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
501 END DO
502!
503! Advance model run counter.
504!
505 nrun=nrun+1
506
507 END DO inner_loop
508
509 END DO outer_loop
510!
511! Report dot products.
512!
513 DO ng=1,ngrids
514 IF (master) THEN
515 WRITE (stdout,80) &
516 & 'TLM Test - Dot Products Summary: p, g1, g2, (g1-g2)/g1'
517 inner=1
518 DO i=1,min(ig1count,ig2count)
519 p=10.0_r8**real(-inner,r8)
520 IF (mod(i,1+ntimes(ng)/ntlm(ng)).eq.0) inner=inner+1
521 WRITE (stdout,90) i, p, g1(i), g2(i), (g1(i)-g2(i))/g1(i)
522 IF ((mod(i,1+ntimes(ng)/ntlm(ng)).eq.0).and. &
523 & (mod(i,ninner*(1+ntimes(ng)/ntlm(ng))).ne.0)) THEN
524 WRITE (stdout,100)
525 ELSE IF (mod(i,ninner*(1+ntimes(ng)/ntlm(ng))).eq.0) THEN
526 inner=1
527 WRITE (stdout,110)
528 END IF
529 END DO
530 END IF
531 END DO
532!
533 10 FORMAT (/,a,i3,/)
534 20 FORMAT (/,1x,a,1x,'ROMS : started time-stepping:', &
535 & ' (Grid: ',i2.2,' TimeSteps: ',i8.8,' - ',i8.8,')',/)
536 40 FORMAT (/,' Nonlinear Model Cost Function = ',1p,e21.14)
537 50 FORMAT (' --------------- ','cost function = ',1p,e21.14,2x,a)
538 60 FORMAT (17x,'cost function = ',1p,e21.14,2x,a)
539 70 FORMAT (a,'_',i3.3,'.nc')
540 80 FORMAT (/,a,/)
541 90 FORMAT (i4,2x,1pe8.1,3(1x,1p,e20.12,0p))
542100 FORMAT (77('.'))
543110 FORMAT (77('-'))
544!
545 RETURN
546 END SUBROUTINE roms_run
547!
548 SUBROUTINE roms_finalize
549!
550!=======================================================================
551! !
552! This routine terminates ROMS nonlinear, tangent linear, and !
553! adjoint models execution. !
554! !
555!=======================================================================
556!
557! Local variable declarations.
558!
559 integer :: Fcount, ng, thread
560!
561 character (len=*), parameter :: MyFile = &
562 & __FILE__//", ROMS_finalize"
563!
564!-----------------------------------------------------------------------
565! If blowing-up, save latest model state into RESTART NetCDF file.
566!-----------------------------------------------------------------------
567!
568! If cycling restart records, write solution into the next record.
569!
570 IF (exit_flag.eq.1) THEN
571 DO ng=1,ngrids
572 IF (lwrtrst(ng)) THEN
573 IF (master) WRITE (stdout,10)
574 10 FORMAT (/,' Blowing-up: Saving latest model state into ', &
575 & ' RESTART file',/)
576 fcount=rst(ng)%load
577 IF (lcyclerst(ng).and.(rst(ng)%Nrec(fcount).ge.2)) THEN
578 rst(ng)%Rindex=2
579 lcyclerst(ng)=.false.
580 END IF
583#ifdef DISTRIBUTE
584 CALL wrt_rst (ng, myrank)
585#else
586 CALL wrt_rst (ng, -1)
587#endif
588 END IF
589 END DO
590 END IF
591!
592!-----------------------------------------------------------------------
593! Stop model and time profiling clocks, report memory requirements, and
594! close output NetCDF files.
595!-----------------------------------------------------------------------
596!
597! Stop time clocks.
598!
599 IF (master) THEN
600 WRITE (stdout,20)
601 20 FORMAT (/,'Elapsed wall CPU time for each process (seconds):',/)
602 END IF
603!
604 DO ng=1,ngrids
605 DO thread=thread_range
606 CALL wclock_off (ng, inlm, 0, __line__, myfile)
607 END DO
608 END DO
609!
610! Report dynamic memory and automatic memory requirements.
611!
612 CALL memory
613!
614! Close IO files.
615!
616 DO ng=1,ngrids
617 CALL close_inp (ng, inlm)
618 END DO
619 CALL close_out
620!
621 RETURN
622 END SUBROUTINE roms_finalize
623
624 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 initial
Definition initial.F:3
subroutine main2d
Definition main2d.F:746
subroutine main3d(runinterval)
Definition main3d.F:4
subroutine memory
Definition memory.F:3
integer function my_threadnum()
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_mod(ng)
Definition def_mod.F:49
subroutine, public ad_dotproduct(ng, tile, linp)
subroutine, public get_state(ng, model, msg, s, inirec, tindex)
Definition get_state.F:90
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
integer ig2count
type(t_fourdvar), dimension(:), allocatable fourdvar
real(r8), dimension(1000) g2
integer, dimension(:), allocatable nobsvar
real(r8), dimension(1000) g1
logical, dimension(:), allocatable wrttlmod
logical, dimension(:), allocatable wrtnlmod
character(len=40), dimension(:), allocatable obsname
integer, dimension(:), allocatable nstatevar
integer ig1count
type(t_io), dimension(:), allocatable his
type(t_io), dimension(:), allocatable adm
type(t_io), dimension(:), allocatable tlm
type(t_io), dimension(:), allocatable fwd
type(t_io), dimension(:), allocatable rst
integer stdout
character(len=256) sourcefile
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, parameter inlm
Definition mod_param.F:662
integer, dimension(:), allocatable ntilex
Definition mod_param.F:685
integer, parameter iadm
Definition mod_param.F:665
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 ninner
integer, dimension(:), allocatable ntimes
integer, dimension(:), allocatable ntlm
integer blowup
integer erend
logical, dimension(:), allocatable lreadfrc
integer, dimension(:), allocatable ndeftlm
logical, dimension(:), allocatable ldefhis
integer, dimension(:), allocatable ntend
logical, dimension(:), allocatable lcycletlm
logical, dimension(:), allocatable ldefmod
integer exit_flag
integer, dimension(:), allocatable nhis
integer erstr
logical, dimension(:), allocatable lwrthis
logical, dimension(:), allocatable lwrttlm
logical, dimension(:), allocatable lwrtrst
logical, dimension(:), allocatable ldeftlm
integer, dimension(:), allocatable ntstart
integer nrun
logical, dimension(:), allocatable lreadfwd
integer inner
integer noerror
logical, dimension(:), allocatable lcyclerst
integer outer
logical, dimension(:), allocatable lreadblk
integer, dimension(:), allocatable lold
integer, dimension(:), allocatable lnew
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