ROMS
Loading...
Searching...
No Matches
main3d.F File Reference
#include "cppdefs.h"
Include dependency graph for main3d.F:

Go to the source code of this file.

Functions/Subroutines

subroutine main3d (runinterval)
 

Function/Subroutine Documentation

◆ main3d()

subroutine main3d ( real(dp), intent(in) runinterval)

Definition at line 3 of file main3d.F.

4!
5!git $Id$
6!=======================================================================
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md Hernan G. Arango !
10!========================================== Alexander F. Shchepetkin ===
11! !
12! This routine is the main driver for ROMS nonlinear model (NLM) !
13! when configurated as a full 3D baroclinic ocean model. It advances !
14! forward the NLM for all nested grids, if any, for the specified !
15! time interval (seconds), RunInterval. !
16! !
17# if defined STEP2D_FB_LF_AM3
18! Numerical 2D time-stepping kernel: FB AB3-AM4 !
19# elif defined STEP2D_FB_LF_AM3
20! Numerical 2D time-stepping kernel: FB LF-AM3 !
21# else
22! Numerical 2D time-stepping kernel: LF-AM3 (Legacy scheme) !
23# endif
24! !
25!=======================================================================
26!
27 USE mod_param
28 USE mod_parallel
29# if defined MODEL_COUPLING && defined MCT_LIB
30 USE mod_coupler
31# endif
32 USE mod_iounits
33# ifdef NESTING
34 USE mod_nesting
35# endif
36 USE mod_scalars
37 USE mod_stepping
38!
39# ifdef ANA_VMIX
40 USE analytical_mod, ONLY : ana_vmix
41# endif
42# ifdef BIOLOGY
43 USE biology_mod, ONLY : biology
44# endif
45# ifdef BBL_MODEL
46 USE bbl_mod, ONLY : bblm
47# endif
48# ifdef BULK_FLUXES
49 USE bulk_flux_mod, ONLY : bulk_flux
50# endif
51# ifdef BVF_MIXING
52 USE bvf_mix_mod, ONLY : bvf_mix
53# endif
54 USE dateclock_mod, ONLY : time_string
55 USE diag_mod, ONLY : diag
56# ifdef TLM_CHECK
58# endif
59# ifdef TIDE_GENERATING_FORCES
60 USE equilibrium_tide_mod, ONLY : equilibrium_tide
61# endif
62# if defined NLM_OUTER || \
63 defined rbl4dvar || \
64 defined rbl4dvar_ana_sensitivity || \
65 defined rbl4dvar_fct_sensitivity || \
66 defined sp4dvar
67 USE forcing_mod, ONLY : forcing
68# endif
69# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
71# endif
72# ifdef GLS_MIXING
75# endif
76# if defined NLM_OUTER || defined RBL4DVAR || \
77 defined rbl4dvar_ana_sensitivity || defined rbl4dvar_fct_sensitivity
78 USE frc_iau_mod, ONLY : frc_iau
79# endif
80# if defined DIFF_3DCOEF || defined VISC_3DCOEF
81 USE hmixing_mod, ONLY : hmixing
82# endif
83# if defined ICE_MODEL && defined ALBEDO && defined SHORTWAVE
84 USE ice_albedo_mod, ONLY : ice_albedo
85# endif
86# ifdef LMD_MIXING
87 USE lmd_vmix_mod, ONLY : lmd_vmix
88# endif
89# if defined ATM_COUPLING && defined MCT_LIB
90 USE mct_coupler_mod, ONLY : ocn2atm_coupling
91# endif
92# if defined WAV_COUPLING && defined MCT_LIB
93 USE mct_coupler_mod, ONLY : ocn2wav_coupling
94# endif
95# ifdef MY25_MIXING
98# endif
99# ifdef NESTING
100 USE nesting_mod, ONLY : nesting
101# ifndef ONE_WAY
102 USE nesting_mod, ONLY : do_twoway
103# endif
104# endif
105# if defined ADJUST_BOUNDARY
107# endif
108 USE omega_mod, ONLY : omega
110# ifndef TS_FIXED
111 USE rho_eos_mod, ONLY : rho_eos
112# endif
113 USE rhs3d_mod, ONLY : rhs3d
114# ifdef ICE_MODEL
115 USE seaice_mod, ONLY : seaice
116# endif
117# ifdef SEDIMENT
118 USE sediment_mod, ONLY : sediment
119# endif
120# ifdef AVERAGES
121 USE set_avg_mod, ONLY : set_avg
122# endif
123 USE set_depth_mod, ONLY : set_depth
125# if defined SSH_TIDES || defined UV_TIDES
126 USE set_tides_mod, ONLY : set_tides
127# endif
128 USE set_vbc_mod, ONLY : set_vbc
129# if !(defined STEP2D_FB_AB3_AM4 || defined STEP2D_FB_LF_AM3)
130 USE set_zeta_mod, ONLY : set_zeta
131# endif
132 USE step2d_mod, ONLY : step2d
133# ifndef TS_FIXED
134 USE step3d_t_mod, ONLY : step3d_t
135# endif
136 USE step3d_uv_mod, ONLY : step3d_uv
137# ifdef FLOATS
138 USE step_floats_mod, ONLY : step_floats
139# endif
140 USE strings_mod, ONLY : founderror
141# ifdef WEC_VF
142# if defined WDISS_THORGUZA || defined WDISS_CHURTHOR
143 USE wec_dissip_mod, ONLY : wec_dissip
144# endif
145# ifdef WEC_ROLLER
146 USE wec_roller_mod, ONLY : wec_roller
147# endif
148 USE wec_stokes_mod, ONLY : wec_stokes
149 USE wec_vf_mod, ONLY : wec_vf
150# ifdef WAVE_MIXING
151 USE wec_wave_mix_mod, ONLY : wec_wave_mix
152# endif
153# endif
154# ifdef WEC
155 USE wec_wvelocity_mod, ONLY : wec_wvelocity
156# endif
157 USE wvelocity_mod, ONLY : wvelocity
158!
159 implicit none
160!
161! Imported variable declarations.
162!
163 real(dp), intent(in) :: RunInterval
164!
165! Local variable declarations.
166!
167 logical :: DoNestLayer, Time_Step
168!
169 integer :: Nsteps, Rsteps
170 integer :: ig, il, istep, ng, nl, tile
171 integer :: my_iif, next_indx1
172# ifdef FLOATS
173 integer :: Lend, Lstr, chunk_size
174# endif
175!
176 character (len=*), parameter :: MyFile = &
177 & __FILE__
178!
179!=======================================================================
180! Time-step nonlinear 3D primitive equations by the specified time.
181!=======================================================================
182!
183! Time-step the 3D kernel for the specified time interval (seconds),
184! RunInterval.
185!
186 time_step=.true.
187 donestlayer=.true.
188!
189 kernel_loop : DO WHILE (time_step)
190!
191! In nesting applications, the number of nesting layers (NestLayers) is
192! used to facilitate refinement grids and composite/refinament grids
193! combinations. Otherwise, the solution it is looped once for a single
194! grid application (NestLayers = 1).
195!
196 nl=0
197#ifdef NESTING
198 twowayinterval(1:ngrids)=0.0_r8
199#endif
200!
201 nest_layer : DO WHILE (donestlayer)
202!
203! Determine number of time steps to compute in each nested grid layer
204! based on the specified time interval (seconds), RunInterval. Non
205! nesting applications have NestLayers=1. Notice that RunInterval is
206! set in the calling driver. Its value may span the full period of the
207! simulation, a multi-model coupling interval (RunInterval > ifac*dt),
208! or just a single step (RunInterval=0).
209!
210 CALL ntimesteps (inlm, runinterval, nl, nsteps, rsteps)
211 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
212 IF ((nl.le.0).or.(nl.gt.nestlayers)) EXIT
213!
214! Time-step governing equations for Nsteps.
215!
216 step_loop : DO istep=1,nsteps
217!
218! Set time indices and time clock.
219!
220 DO ig=1,gridsinlayer(nl)
221 ng=gridnumber(ig,nl)
222 nstp(ng)=1+mod(iic(ng)-ntstart(ng),2)
223 nnew(ng)=3-nstp(ng)
224 nrhs(ng)=nstp(ng)
225# ifdef JEDI
226 jic(ng)=jic(ng)+1
227 time4jedi(ng)=time4jedi(ng)+dt(ng)
228# endif
229 tdays(ng)=time(ng)*sec2day
230 IF (step_counter(ng).eq.rsteps) time_step=.false.
231 END DO
232!
233!-----------------------------------------------------------------------
234! Read in required data, if any, from input NetCDF files.
235!-----------------------------------------------------------------------
236!
237 DO ig=1,gridsinlayer(nl)
238 ng=gridnumber(ig,nl)
239!$OMP MASTER
240 IF (processinputdata(ng)) THEN
241 CALL get_data (ng)
242 END IF
243!$OMP END MASTER
244!$OMP BARRIER
246 & __line__, myfile)) RETURN
247 END DO
248!
249!-----------------------------------------------------------------------
250! If applicable, process input data: time interpolate between data
251! snapshots.
252!-----------------------------------------------------------------------
253!
254 DO ig=1,gridsinlayer(nl)
255 ng=gridnumber(ig,nl)
256 IF (processinputdata(ng)) THEN
257 DO tile=first_tile(ng),last_tile(ng),+1
258 CALL set_data (ng, tile)
259 END DO
260 END IF
261!$OMP BARRIER
262 END DO
263 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
264
265# if defined NLM_OUTER || \
266 defined rbl4dvar || \
267 defined rbl4dvar_ana_sensitivity || \
268 defined rbl4dvar_fct_sensitivity
269!
270!-----------------------------------------------------------------------
271! If appropriate, add convolved adjoint solution impulse forcing to
272! the nonlinear model solution. Notice that the forcing is only needed
273! after finishing all the inner loops. The forcing is continuous.
274! That is, it is time interpolated at every time-step from available
275! snapshots (FrequentImpulse=TRUE).
276!-----------------------------------------------------------------------
277!
278 DO ig=1,gridsinlayer(nl)
279 ng=gridnumber(ig,nl)
280# ifdef RPCG
281 IF ((iic(ng).gt.int(timeiau(ng)/dt(ng))).or. &
282 & timeiau(ng).eq.0.0_dp) THEN
283 iauswitch(ng)=.false.
284 END IF
285 IF (frequentimpulse(ng).or.iauswitch(ng)) THEN
286# endif
287# if defined WEAK_NOINTERP
288 IF ((iic(ng).gt.1).and.(iic(ng).ne.ntend(ng)+1).and. &
289 & (mod(iic(ng)-1,nadj(ng)).eq.0)) THEN
290 IF (master) THEN
291 WRITE (stdout,*) ' FORCING NLM at iic = ',iic(ng)
292 END IF
293# endif
294# ifdef RPCG
295 IF (frequentimpulse(ng).and.iauswitch(ng)) THEN
296 DO tile=first_tile(ng),last_tile(ng),+1
297 CALL frc_iau (ng, tile, 1)
298 END DO
299 END IF
300 DO tile=first_tile(ng),last_tile(ng),+1
301 CALL forcing (ng, tile, kstp(ng), nstp(ng))
302 CALL set_depth (ng, tile, inlm)
303 END DO
304!$OMP BARRIER
305# else
306 IF (frequentimpulse(ng)) THEN
307 DO tile=first_tile(ng),last_tile(ng),+1
308 CALL forcing (ng, tile, kstp(ng), nstp(ng))
309 CALL set_depth (ng, tile, inlm)
310 END DO
311!$OMP BARRIER
312 END IF
313# endif
314# if defined WEAK_NOINTERP
315 END IF
316# endif
317# ifdef RPCG
318 END IF
319# endif
320 END DO
321# endif
322
323# ifndef JEDI
324!
325!-----------------------------------------------------------------------
326! On the first timestep, it computes the initial depths and level
327! thicknesses from the initial free-surface field. Additionally, it
328! initializes the nonlinear state variables for all time levels and
329! applies lateral boundary conditions.
330!-----------------------------------------------------------------------
331!
332 DO ig=1,gridsinlayer(nl)
333 ng=gridnumber(ig,nl)
334 IF (iic(ng).eq.ntstart(ng)) THEN
335 CALL post_initial (ng, inlm)
336 END IF
337 END DO
338# endif
339!
340!-----------------------------------------------------------------------
341! Compute horizontal mass fluxes (Hz*u/n and Hz*v/m), density related
342! quatities and report global diagnostics.
343!-----------------------------------------------------------------------
344!
345 DO ig=1,gridsinlayer(nl)
346 ng=gridnumber(ig,nl)
347 DO tile=first_tile(ng),last_tile(ng),+1
348 CALL set_massflux (ng, tile, inlm)
349# ifndef TS_FIXED
350 CALL rho_eos (ng, tile, inlm)
351# endif
352# ifdef TIDE_GENERATING_FORCES
353 CALL equilibrium_tide (ng, tile, inlm)
354# endif
355 CALL diag (ng, tile)
356# ifdef TLM_CHECK
357 CALL nl_dotproduct (ng, tile, lnew(ng))
358# endif
359 END DO
360!$OMP BARRIER
361 END DO
362 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
363
364# if defined ATM_COUPLING && defined MCT_LIB
365!
366!-----------------------------------------------------------------------
367! Couple ocean to atmosphere model every "CoupleSteps(Iatmos)"
368! timesteps: get air/sea fluxes.
369!-----------------------------------------------------------------------
370!
371 DO ig=1,gridsinlayer(nl)
372 ng=gridnumber(ig,nl)
373 IF ((iic(ng).ne.ntstart(ng)).and. &
374 & mod(iic(ng)-1,couplesteps(iatmos,ng)).eq.0) THEN
375 DO tile=last_tile(ng),first_tile(ng),-1
376 CALL ocn2atm_coupling (ng, tile)
377 END DO
378!$OMP BARRIER
379 END IF
380 END DO
381# endif
382
383# if defined WAV_COUPLING && defined MCT_LIB
384!
385!-----------------------------------------------------------------------
386! Couple to ocean to waves model every "CoupleSteps(Iwaves)"
387! timesteps: get waves/ocean fluxes.
388!-----------------------------------------------------------------------
389!
390 DO ig=1,gridsinlayer(nl)
391 ng=gridnumber(ig,nl)
392 IF ((iic(ng).ne.ntstart(ng)).and. &
393 & mod(iic(ng)-1,couplesteps(iwaves,ng)).eq.0) THEN
394 DO tile=first_tile(ng),last_tile(ng),+1
395 CALL ocn2wav_coupling (ng, tile)
396 END DO
397!$OMP BARRIER
398 END IF
399 END DO
400# endif
401
402# ifdef WEC_VF
403!
404!-----------------------------------------------------------------------
405! Compute wave dissipation, roller, and stokes terms.
406!-----------------------------------------------------------------------
407!
408 DO ig=1,gridsinlayer(nl)
409 ng=gridnumber(ig,nl)
410 DO tile=first_tile(ng),last_tile(ng),+1
411# if defined WDISS_THORGUZA || defined WDISS_CHURTHOR
412 CALL wec_dissip (ng, tile)
413# endif
414# ifdef WEC_ROLLER
415 CALL wec_roller (ng, tile)
416# endif
417 CALL wec_stokes (ng, tile)
418 END DO
419!$OMP BARRIER
420 END DO
421 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
422# endif
423!
424!-----------------------------------------------------------------------
425! Set fields for vertical boundary conditions. Process tidal forcing,
426! if any.
427!-----------------------------------------------------------------------
428!
429 DO ig=1,gridsinlayer(nl)
430 ng=gridnumber(ig,nl)
431 DO tile=first_tile(ng),last_tile(ng),+1
432# ifdef BULK_FLUXES
433# if defined ICE_MODEL && defined ALBEDO && defined SHORTWAVE
434 CALL ice_albedo (ng, tile, inlm)
435# endif
436# if defined FOUR_DVAR && defined PRIOR_BULK_FLUXES
437 IF (nrun.eq.1) CALL bulk_flux (ng, tile)
438# else
439 CALL bulk_flux (ng, tile)
440# endif
441# endif
442# ifdef BBL_MODEL
443 CALL bblm (ng, tile)
444# endif
445 CALL set_vbc (ng, tile)
446# if defined SSH_TIDES || defined UV_TIDES
447 CALL set_tides (ng, tile)
448# endif
449 END DO
450!$OMP BARRIER
451 END DO
452
453# ifdef NESTING
454!
455! If composite or mosaic grids, process additional points in the
456! contact zone between connected grids for bottom stress variables.
457!
458 DO ig=1,gridsinlayer(nl)
459 ng=gridnumber(ig,nl)
460 IF (any(compositegrid(:,ng))) THEN
461 CALL nesting (ng, inlm, nbstr)
462 END IF
463 END DO
464# endif
465
466# if defined ICE_MODEL
467!
468!-----------------------------------------------------------------------
469! Timestep the ice model.
470!-----------------------------------------------------------------------
471!
472 CALL seaice (inlm)
473# endif
474
475# ifdef ADJUST_BOUNDARY
476!
477!-----------------------------------------------------------------------
478! Interpolate open boundary increments and adjust open boundary.
479! Load open boundary into storage arrays. Skip the last output
480! timestep.
481!-----------------------------------------------------------------------
482!
483 DO ig=1,gridsinlayer(nl)
484 ng=gridnumber(ig,nl)
485 IF (iic(ng).lt.(ntend(ng)+1)) THEN
486 DO tile=first_tile(ng),last_tile(ng),+1
487 CALL obc_adjust (ng, tile, lbinp(ng))
488 CALL load_obc (ng, tile, lbout(ng))
489 END DO
490!$OMP BARRIER
491 END IF
492 END DO
493# endif
494
495# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
496!
497!-----------------------------------------------------------------------
498! Interpolate surface forcing increments and adjust surface forcing.
499! Load surface forcing into storage arrays. Skip the last output
500! timestep.
501!-----------------------------------------------------------------------
502!
503 DO ig=1,gridsinlayer(nl)
504 ng=gridnumber(ig,nl)
505 IF (iic(ng).lt.(ntend(ng)+1)) THEN
506 DO tile=first_tile(ng),last_tile(ng),+1
507 CALL frc_adjust (ng, tile, lfinp(ng))
508 CALL load_frc (ng, tile, lfout(ng))
509 END DO
510!$OMP BARRIER
511 END IF
512 END DO
513# endif
514!
515!-----------------------------------------------------------------------
516! Compute time-dependent vertical/horizontal mixing coefficients for
517! momentum and tracers. Compute S-coordinate vertical velocity,
518! diagnostically from horizontal mass divergence.
519!-----------------------------------------------------------------------
520!
521 DO ig=1,gridsinlayer(nl)
522 ng=gridnumber(ig,nl)
523 DO tile=last_tile(ng),first_tile(ng),-1
524# if defined ANA_VMIX
525 CALL ana_vmix (ng, tile, inlm)
526# elif defined LMD_MIXING
527 CALL lmd_vmix (ng, tile)
528# elif defined BVF_MIXING
529 CALL bvf_mix (ng, tile)
530# endif
531# if defined DIFF_3DCOEF || defined VISC_3DCOEF
532 CALL hmixing (ng, tile)
533# endif
534 CALL omega (ng, tile, inlm)
535 CALL wvelocity (ng, tile, nstp(ng))
536# ifdef WEC
537 CALL wec_wvelocity (ng, tile, nstp(ng))
538# endif
539 END DO
540!$OMP BARRIER
541 END DO
542
543# if !(defined STEP2D_FB_AB3_AM4 || defined STEP2D_FB_LF_AM3) || \
544 defined diagnostics || defined averages
545!
546!-----------------------------------------------------------------------
547! Set free-surface to it time-averaged value. If applicable,
548! accumulate time-averaged output data which needs a irreversible
549! loop in shared-memory jobs.
550!-----------------------------------------------------------------------
551!
552 DO ig=1,gridsinlayer(nl)
553 ng=gridnumber(ig,nl)
554 DO tile=first_tile(ng),last_tile(ng),+1 ! irreversible
555# if !(defined STEP2D_FB_AB3_AM4 || defined STEP2D_FB_LF_AM3)
556 CALL set_zeta (ng, tile)
557# endif
558# ifdef DIAGNOSTICS
559 CALL set_diags (ng, tile)
560# endif
561# ifdef AVERAGES
562 CALL set_avg (ng, tile)
563# endif
564 END DO
565!$OMP BARRIER
566 END DO
567# endif
568
569# ifdef NESTING
570!
571! If composite or mosaic grids, process additional points in the
572! contact zone between connected grids for 3D kernel free-surface.
573!
574 DO ig=1,gridsinlayer(nl)
575 ng=gridnumber(ig,nl)
576 IF (any(compositegrid(:,ng))) THEN
577 CALL nesting (ng, inlm, nzeta)
578 END IF
579 END DO
580# endif
581!
582!-----------------------------------------------------------------------
583! If appropriate, write out fields into output NetCDF files. Notice
584! that IO data is written in delayed and serial mode. Exit if last
585! time step.
586!-----------------------------------------------------------------------
587!
588 DO ig=1,gridsinlayer(nl)
589 ng=gridnumber(ig,nl)
590!$OMP MASTER
591 CALL output (ng)
592!$OMP END MASTER
593!$OMP BARRIER
594 IF ((founderror(exit_flag, noerror, __line__, myfile)).or.&
595 & ((iic(ng).eq.(ntend(ng)+1)).and.(ng.eq.ngrids))) THEN
596 RETURN
597 END IF
598 END DO
599
600# ifdef NESTING
601!
602!-----------------------------------------------------------------------
603! If refinement grid, interpolate (space, time) state variables
604! contact points from donor grid extracted data.
605# ifdef NESTING_DEBUG
606!
607! Also, fill BRY_CONTACT(:,:)%Mflux to check for mass conservation
608! between coarse and fine grids. This is only done for diagnostic
609! purposes. Also, debugging is possible with very verbose output
610! to fort.300 is allowed by activating uppercase(nesting_debug).
611# endif
612!-----------------------------------------------------------------------
613!
614 DO ig=1,gridsinlayer(nl)
615 ng=gridnumber(ig,nl)
616 IF (refinedgrid(ng).and.(refinescale(ng).gt.0)) THEN
617 CALL nesting (ng, inlm, nputd)
618# ifdef NESTING_DEBUG
619 CALL nesting (ng, inlm, nmflx)
620# endif
621 END IF
622 END DO
623# endif
624!
625!-----------------------------------------------------------------------
626! Compute right-hand-side terms for 3D equations.
627!-----------------------------------------------------------------------
628!
629 DO ig=1,gridsinlayer(nl)
630 ng=gridnumber(ig,nl)
631 DO tile=last_tile(ng),first_tile(ng),-1
632 CALL rhs3d (ng, tile)
633# ifdef MY25_MIXING
634 CALL my25_prestep (ng, tile)
635# elif defined GLS_MIXING
636 CALL gls_prestep (ng, tile)
637# endif
638 END DO
639!$OMP BARRIER
640 END DO
641
642# ifdef NESTING
643!
644! If composite or mosaic grids, process additional points in the
645! contact zone between connected grids for right-hand-side terms
646! (tracers).
647!
648 DO ig=1,gridsinlayer(nl)
649 ng=gridnumber(ig,nl)
650 IF (any(compositegrid(:,ng))) THEN
651 CALL nesting (ng, inlm, nrhst)
652 END IF
653 END DO
654# endif
655
656# ifdef STEP2D_FB_AB3_AM4
657!
658!-----------------------------------------------------------------------
659! Solve nonlinear vertically integrated primitive equations for
660! free-surface and momentum components using a generalized Forward-
661! Backward, 3rd-order Adams-Bashforth / 4th-order Adams-Moulton
662! (FB AB3-AM4) time stepping scheme (Shchepetkin and McWilliams,
663! 2009).
664!-----------------------------------------------------------------------
665!
666 loop_2d : DO my_iif=1,maxval(nfast)
667
668 DO ig=1,gridsinlayer(nl)
669 ng=gridnumber(ig,nl)
670 IF (my_iif.le.nfast(ng)) THEN
671 iif(ng)=my_iif
672 kstp(ng)=knew(ng)
673 knew(ng)=kstp(ng)+1
674 IF (knew(ng).gt.4) knew(ng)=1
675
676 IF (mod(knew(ng),2).eq.0) THEN ! zig-zag
677 DO tile=first_tile(ng),last_tile(ng),+1 ! processing
678 CALL step2d (ng, tile) ! sequence
679 END DO
680!$OMP BARRIER
681 ELSE
682 DO tile=last_tile(ng),first_tile(ng),-1
683 CALL step2d (ng, tile)
684 END DO
685!$OMP BARRIER
686 END IF
687 END IF
688 END DO
689
690# ifdef NESTING
691!
692! If composite or mosaic grids, process additional points in the
693! contact zone between connected grids for the state variables
694! associated with the 2D engine CORRECTOR STEP section (KNEW INDEX).
695!
696 DO ig=1,gridsinlayer(nl)
697 ng=gridnumber(ig,nl)
698 IF (any(compositegrid(:,ng))) THEN
699 CALL nesting (ng, inlm, n2dcs)
700 END IF
701 END DO
702# endif
703 END DO loop_2d
704
705# else
706
707# ifdef STEP2D_FB_LF_AM3
708!
709!-----------------------------------------------------------------------
710! Solve the vertically integrated primitive equations for the
711! free-surface and barotropic momentum components using a predictor-
712! corrector LeapFrog / 3rd-order Adams-Moulton with a Forward-Backward
713! feeback (FB LF-AM3) stepping scheme (Shchepetkin and McWilliams,
714! 2009).
715!-----------------------------------------------------------------------
716!
717 loop_2d : DO my_iif=1,maxval(nfast)
718!
719! Predictor LF substep with FB-feedback.
720!
721 DO ig=1,gridsinlayer(nl)
722 ng=gridnumber(ig,nl)
723 IF (my_iif.le.nfast(ng)) THEN
724 iif(ng)=my_iif
725 kstp(ng)=next_kstp(ng)
726 knew(ng)=3
727
728 DO tile=last_tile(ng),first_tile(ng),-1
729 CALL step2d (ng, tile)
730 END DO
731!$OMP BARRIER
732 END IF
733 END DO
734
735# ifdef NESTING
736!
737! If composite or mosaic grids, process additional points in the
738! contact zone between connected grids for the state variables
739! associated with the 2D engine PREDICTOR STEP section.
740# ifdef NESTING_DEBUG
741! If refinement, check mass flux conservation between coarse and
742! fine grids during debugging.
743! Warning: very verbose output to fort.300 ascii file to check
744! mass flux conservation.
745# endif
746!
747 DO ig=1,gridsinlayer(nl)
748 ng=gridnumber(ig,nl)
749 IF (any(compositegrid(:,ng))) THEN
750 CALL nesting (ng, inlm, n2dps)
751 END IF
752# ifdef NESTING_DEBUG
753 IF (refinedgrid(ng).and.(refinescale(ng).gt.0)) THEN
754 CALL nesting (ng, inlm, nmflx)
755 END IF
756# endif
757 END DO
758# endif
759!
760! Corrector AM3 substep with FB-feedback.
761!
762 DO ig=1,gridsinlayer(nl)
763 ng=gridnumber(ig,nl)
764 IF (my_iif.le.nfast(ng)) THEN
765 knew(ng)=3-kstp(ng)
766 next_kstp(ng)=knew(ng)
767
768 DO tile=first_tile(ng),last_tile(ng),+1
769 CALL step2d (ng, tile)
770 END DO
771!$OMP BARRIER
772 END IF
773 END DO
774
775# ifdef NESTING
776!
777! If composite or mosaic grids, process additional points in the
778! contact zone between connected grids for the state variables
779! associated with the 2D engine CORRECTOR STEP section (KNEW INDEX).
780# ifdef NESTING_DEBUG
781! If refinement, check mass flux conservation between coarse and
782! fine grids during debugging.
783! Warning: very verbose output to fort.300 ascii file to check
784! mass flux conservation.
785# endif
786!
787 DO ig=1,gridsinlayer(nl)
788 ng=gridnumber(ig,nl)
789 IF (any(compositegrid(:,ng))) THEN
790 CALL nesting (ng, inlm, n2dcs)
791 END IF
792# ifdef NESTING_DEBUG
793 IF (refinedgrid(ng).and.(refinescale(ng).gt.0)) THEN
794 CALL nesting (ng, inlm, nmflx)
795 END IF
796# endif
797 END DO
798# endif
799 END DO loop_2d
800
801# else
802!
803!-----------------------------------------------------------------------
804! Solve the vertically integrated primitive equations for the
805! free-surface and barotropic momentum components using a predictor-
806! corrector LeapFrog with 3rd-order Adams-Moulton (LF-AM3) time
807! stepping scheme.
808!-----------------------------------------------------------------------
809!
810 loop_2d : DO my_iif=1,maxval(nfast)+1
811!
812! Set time indices for predictor step. The PREDICTOR_2D_STEP switch
813! it is assumed to be false before the first time-step.
814!
815 DO ig=1,gridsinlayer(nl)
816 ng=gridnumber(ig,nl)
817 next_indx1=3-indx1(ng)
818 IF (.not.predictor_2d_step(ng).and. &
819 & my_iif.le.(nfast(ng)+1)) THEN
820 predictor_2d_step(ng)=.true.
821 iif(ng)=my_iif
822 IF (first_2d_step) THEN
823 kstp(ng)=indx1(ng)
824 ELSE
825 kstp(ng)=3-indx1(ng)
826 END IF
827 knew(ng)=3
828 krhs(ng)=indx1(ng)
829 END IF
830!
831! Predictor step - Advance barotropic equations using 2D time-step
832! ============== predictor scheme. No actual time-stepping is
833! performed during the auxiliary (nfast+1) time-step. It is needed
834! to finalize the fast-time averaging of 2D fields, if any, and
835! compute the new time-evolving depths.
836!
837 IF (my_iif.le.(nfast(ng)+1)) THEN
838 DO tile=last_tile(ng),first_tile(ng),-1
839 CALL step2d (ng, tile)
840 END DO
841!$OMP BARRIER
842 END IF
843 END DO
844
845# ifdef NESTING
846!
847! If composite or mosaic grids, process additional points in the
848! contact zone between connected grids for the state variables
849! associated with the 2D engine PREDICTOR STEP section.
850# ifdef NESTING_DEBUG
851! If refinement, check mass flux conservation between coarse and
852! fine grids during debugging.
853! Warning: very verbose output to fort.300 ascii file to check
854! mass flux conservation.
855# endif
856!
857 DO ig=1,gridsinlayer(nl)
858 ng=gridnumber(ig,nl)
859 IF (any(compositegrid(:,ng))) THEN
860 CALL nesting (ng, inlm, n2dps)
861 END IF
862# ifdef NESTING_DEBUG
863 IF (refinedgrid(ng).and.(refinescale(ng).gt.0)) THEN
864 CALL nesting (ng, inlm, nmflx)
865 END IF
866# endif
867 END DO
868# endif
869!
870! Set time indices for corrector step.
871!
872 DO ig=1,gridsinlayer(nl)
873 ng=gridnumber(ig,nl)
874 IF (predictor_2d_step(ng)) THEN
875 predictor_2d_step(ng)=.false.
876 knew(ng)=next_indx1
877 kstp(ng)=3-knew(ng)
878 krhs(ng)=3
879 IF (iif(ng).lt.(nfast(ng)+1)) indx1(ng)=next_indx1
880 END IF
881!
882! Corrector step - Apply 2D time-step corrector scheme. Notice that
883! ============== there is not need for a corrector step during the
884! auxiliary (nfast+1) time-step.
885!
886 IF (iif(ng).lt.(nfast(ng)+1)) THEN
887 DO tile=first_tile(ng),last_tile(ng),+1
888 CALL step2d (ng, tile)
889 END DO
890!$OMP BARRIER
891 END IF
892 END DO
893
894# ifdef NESTING
895!
896! If composite or mosaic grids, process additional points in the
897! contact zone between connected grids for the state variables
898! associated with the 2D engine CORRECTOR STEP section (KNEW INDEX).
899# ifdef NESTING_DEBUG
900! If refinement, check mass flux conservation between coarse and
901! fine grids during debugging.
902! Warning: very verbose output to fort.300 ascii file to check
903! mass flux conservation.
904# endif
905!
906 DO ig=1,gridsinlayer(nl)
907 ng=gridnumber(ig,nl)
908 IF (any(compositegrid(:,ng))) THEN
909 CALL nesting (ng, inlm, n2dcs)
910 END IF
911# ifdef NESTING_DEBUG
912 IF (refinedgrid(ng).and.(refinescale(ng).gt.0)) THEN
913 CALL nesting (ng, inlm, nmflx)
914 END IF
915# endif
916 END DO
917# endif
918 END DO loop_2d
919# endif
920
921# endif
922
923# ifdef NESTING
924# if defined MASKING && defined WET_DRY
925!
926!-----------------------------------------------------------------------
927! If nesting and wetting and drying, scale horizontal interpolation
928! weights to account for land/sea masking in contact areas. This needs
929! to be done at very time-step since the Land/Sea masking is time
930! dependent.
931!-----------------------------------------------------------------------
932!
933 DO ig=1,gridsinlayer(nl)
934 ng=gridnumber(ig,nl)
935 CALL nesting (ng, inlm, nmask)
936 END DO
937# endif
938!
939!-----------------------------------------------------------------------
940! If composite or mosaic grids, process additional points in the
941! contact zone between connected grids for the time-averaged
942! momentum fluxes (DU_avg1, DV_avg1) and free-surface (Zt_avg).
943!-----------------------------------------------------------------------
944!
945 DO ig=1,gridsinlayer(nl)
946 ng=gridnumber(ig,nl)
947 IF (any(compositegrid(:,ng))) THEN
948 CALL nesting (ng, inlm, n2dfx)
949 END IF
950 END DO
951# endif
952
953# if !(defined STEP2D_FB_AB3_AM4 || defined STEP2D_FB_LF_AM3)
954!
955!-----------------------------------------------------------------------
956! Recompute depths and thicknesses using the new time filtered
957! free-surface.
958!-----------------------------------------------------------------------
959!
960 DO ig=1,gridsinlayer(nl)
961 ng=gridnumber(ig,nl)
962 DO tile=last_tile(ng),first_tile(ng),-1
963 CALL set_depth (ng, tile, inlm)
964 END DO
965!$OMP BARRIER
966 END DO
967# endif
968
969# ifdef NESTING
970!
971! If nesting, determine vertical indices and vertical interpolation
972! weights in the contact zone using new depth arrays.
973!
974 DO ig=1,gridsinlayer(nl)
975 ng=gridnumber(ig,nl)
976! CALL nesting (ng, iNLM, nzwgt)
977 END DO
978# endif
979!
980!-----------------------------------------------------------------------
981! Time-step 3D momentum equations.
982!-----------------------------------------------------------------------
983!
984! Time-step 3D momentum equations and couple with vertically
985! integrated equations.
986!
987 DO ig=1,gridsinlayer(nl)
988 ng=gridnumber(ig,nl)
989 DO tile=last_tile(ng),first_tile(ng),-1
990 CALL step3d_uv (ng, tile)
991 END DO
992!$OMP BARRIER
993 END DO
994
995# ifdef NESTING
996!
997! If composite or mosaic grids, process additional points in the
998! contact zone between connected grids for 3D momentum (u,v),
999! adjusted 2D momentum (ubar,vbar), and fluxes (Huon, Hvom).
1000!
1001 DO ig=1,gridsinlayer(nl)
1002 ng=gridnumber(ig,nl)
1003 IF (any(compositegrid(:,ng))) THEN
1004 CALL nesting (ng, inlm, n3duv)
1005 END IF
1006 END DO
1007# endif
1008!
1009!-----------------------------------------------------------------------
1010! Time-step vertical mixing turbulent equations and passive tracer
1011! source and sink terms, if applicable.
1012!-----------------------------------------------------------------------
1013!
1014 DO ig=1,gridsinlayer(nl)
1015 ng=gridnumber(ig,nl)
1016 DO tile=first_tile(ng),last_tile(ng),+1
1017 CALL omega (ng, tile, inlm)
1018# ifdef MY25_MIXING
1019 CALL my25_corstep (ng, tile)
1020# elif defined GLS_MIXING
1021 CALL gls_corstep (ng, tile)
1022# endif
1023# if defined WEC_VF && defined WAVE_MIXING
1024 CALL wec_wave_mix (ng, tile)
1025# endif
1026# ifdef BIOLOGY
1027 CALL biology (ng, tile)
1028# endif
1029# ifdef SEDIMENT
1030 CALL sediment (ng, tile)
1031# endif
1032 END DO
1033!$OMP BARRIER
1034 END DO
1035
1036# ifndef TS_FIXED
1037!
1038!-----------------------------------------------------------------------
1039! Time-step tracer equations.
1040!-----------------------------------------------------------------------
1041!
1042 DO ig=1,gridsinlayer(nl)
1043 ng=gridnumber(ig,nl)
1044 DO tile=last_tile(ng),first_tile(ng),-1
1045 CALL step3d_t (ng, tile)
1046 END DO
1047!$OMP BARRIER
1048 END DO
1049
1050# ifdef NESTING
1051!
1052! If composite or mosaic grids, process additional points in the
1053! contact zone between connected grids for Tracer Variables.
1054!
1055 DO ig=1,gridsinlayer(nl)
1056 ng=gridnumber(ig,nl)
1057 IF (any(compositegrid(:,ng))) THEN
1058 CALL nesting (ng, inlm, n3dtv)
1059 END IF
1060 END DO
1061# endif
1062# endif
1063
1064# ifdef NESTING
1065# ifndef ONE_WAY
1066!
1067!-----------------------------------------------------------------------
1068! If refinement grids, perform two-way coupling between fine and
1069! coarse grids. Correct coarse grid tracers values at the refinement
1070! grid with refined accumulated fluxes. Then, replace coarse grid
1071! state variable with averaged refined grid values (two-way nesting).
1072! Update coarse grid depth variables.
1073!
1074! The two-way exchange of infomation between nested grids needs to be
1075! done at the correct time-step and in the right sequence.
1076!-----------------------------------------------------------------------
1077!
1078 DO il=nestlayers,1,-1
1079 DO ig=1,gridsinlayer(il)
1080 ng=gridnumber(ig,il)
1081 IF (do_twoway(inlm, nl, il, ng, istep)) THEN
1082 CALL nesting (ng, inlm, n2way)
1083 END IF
1084 END DO
1085 END DO
1086# endif
1087!
1088!-----------------------------------------------------------------------
1089! If donor to a finer grid, extract data for the external contact
1090! points. This is the latest solution for the coarser grid.
1091!
1092! It is stored in the REFINED structure so it can be used for the
1093! space-time interpolation when "nputD" argument is used above.
1094!-----------------------------------------------------------------------
1095!
1096 DO ig=1,gridsinlayer(nl)
1097 ng=gridnumber(ig,nl)
1098 IF (donortofiner(ng)) THEN
1099 CALL nesting (ng, inlm, ngetd)
1100 END IF
1101 END DO
1102# endif
1103
1104# ifdef FLOATS
1105!
1106!-----------------------------------------------------------------------
1107! Compute Lagrangian drifters trajectories: Split all the drifters
1108! between all the computational threads, except in distributed-memory
1109! and serial configurations. In distributed-memory, the parallel node
1110! containing the drifter is selected internally since the state
1111! variables do not have a global scope.
1112!-----------------------------------------------------------------------
1113!
1114 DO ig=1,gridsinlayer(nl)
1115 ng=gridnumber(ig,nl)
1116 IF (lfloats(ng)) THEN
1117# ifdef _OPENMP
1118 chunk_size=(nfloats(ng)+numthreads-1)/numthreads
1119 lstr=1+mythread*chunk_size
1120 lend=min(nfloats(ng),lstr+chunk_size-1)
1121# else
1122 lstr=1
1123 lend=nfloats(ng)
1124# endif
1125 CALL step_floats (ng, lstr, lend)
1126!$OMP BARRIER
1127!
1128! Shift floats time indices.
1129!
1130 nfp1(ng)=mod(nfp1(ng)+1,nft+1)
1131 nf(ng)=mod(nf(ng)+1,nft+1)
1132 nfm1(ng)=mod(nfm1(ng)+1,nft+1)
1133 nfm2(ng)=mod(nfm2(ng)+1,nft+1)
1134 nfm3(ng)=mod(nfm3(ng)+1,nft+1)
1135 END IF
1136 END DO
1137# endif
1138!
1139!-----------------------------------------------------------------------
1140! Advance time index and time clock.
1141!-----------------------------------------------------------------------
1142!
1143 DO ig=1,gridsinlayer(nl)
1144 ng=gridnumber(ig,nl)
1145 iic(ng)=iic(ng)+1
1146 time(ng)=time(ng)+dt(ng)
1147 step_counter(ng)=step_counter(ng)-1
1148 CALL time_string (time(ng), time_code(ng))
1149 END DO
1150
1151 END DO step_loop
1152
1153 END DO nest_layer
1154
1155 END DO kernel_loop
1156
1157 RETURN
subroutine get_data(ng)
Definition get_data.F:3
subroutine ana_vmix(ng, tile, model)
Definition ana_vmix.h:3
subroutine, public bblm(ng, tile)
Definition mb_bbl.h:55
subroutine, public biology(ng, tile)
Definition ecosim.h:57
subroutine, public bulk_flux(ng, tile)
Definition bulk_flux.F:107
subroutine, public bvf_mix(ng, tile)
Definition bvf_mix.F:28
subroutine, public time_string(mytime, date_string)
Definition dateclock.F:1272
Definition diag.F:2
subroutine, public diag(ng, tile)
Definition diag.F:24
subroutine, public nl_dotproduct(ng, tile, linp)
Definition dotproduct.F:40
subroutine, public forcing(ng, tile, kfrc, nfrc)
Definition forcing.F:30
subroutine, public load_frc(ng, tile, lout)
Definition frc_adjust.F:221
subroutine, public frc_adjust(ng, tile, linp)
Definition frc_adjust.F:36
subroutine, public frc_iau(ng, tile, irec)
Definition frc_iau.F:240
subroutine, public gls_corstep(ng, tile)
Definition gls_corstep.F:38
subroutine, public gls_prestep(ng, tile)
Definition gls_prestep.F:30
subroutine, public lmd_vmix(ng, tile)
Definition lmd_vmix.F:34
integer, dimension(:,:), allocatable couplesteps
integer iatmos
integer iwaves
integer stdout
integer, parameter nrhst
Definition mod_nesting.F:87
integer, parameter ngetd
Definition mod_nesting.F:77
integer, parameter n2dfx
Definition mod_nesting.F:92
integer, parameter nputd
Definition mod_nesting.F:79
integer, parameter nzeta
Definition mod_nesting.F:88
real(r8), dimension(:), allocatable twowayinterval
integer, parameter nmask
Definition mod_nesting.F:78
logical, dimension(:), allocatable donortofiner
integer, parameter n2way
Definition mod_nesting.F:80
integer, parameter n3duv
Definition mod_nesting.F:93
integer, parameter n2dps
Definition mod_nesting.F:90
integer, parameter n3dtv
Definition mod_nesting.F:94
integer, parameter nbstr
Definition mod_nesting.F:86
integer, parameter nmflx
Definition mod_nesting.F:75
integer, parameter n2dcs
Definition mod_nesting.F:91
integer numthreads
integer, dimension(:), allocatable first_tile
integer mythread
logical master
integer, dimension(:), allocatable last_tile
integer, dimension(:), allocatable nfloats
Definition mod_param.F:543
integer, parameter inlm
Definition mod_param.F:662
integer, dimension(:,:), allocatable gridnumber
Definition mod_param.F:127
integer nestlayers
Definition mod_param.F:118
integer, parameter nft
Definition mod_param.F:539
integer ngrids
Definition mod_param.F:113
integer, dimension(:), allocatable gridsinlayer
Definition mod_param.F:122
logical, dimension(:), allocatable lfloats
integer, dimension(:), allocatable iic
real(dp), dimension(:), allocatable dt
real(dp), dimension(:), allocatable timeiau
integer, dimension(:), allocatable next_kstp
logical, dimension(:), allocatable predictor_2d_step
integer, dimension(:), allocatable jic
real(dp), dimension(:), allocatable tdays
logical, dimension(:), allocatable frequentimpulse
integer, dimension(:), allocatable nfast
real(dp), parameter sec2day
integer, dimension(:), allocatable ntend
character(len=22), dimension(:), allocatable time_code
integer exit_flag
logical, dimension(:), allocatable iauswitch
real(dp), dimension(:), allocatable time4jedi
integer, dimension(:), allocatable indx1
logical, dimension(:,:), allocatable compositegrid
real(dp), dimension(:), allocatable time
logical, dimension(:), allocatable refinedgrid
integer, dimension(:), allocatable refinescale
integer, dimension(:), allocatable ntstart
integer nrun
integer, dimension(:), allocatable step_counter
logical, dimension(:), allocatable processinputdata
integer, dimension(:), allocatable nadj
integer, dimension(:), allocatable iif
integer noerror
integer, dimension(:), allocatable nfm2
integer, dimension(:), allocatable lbout
integer, dimension(:), allocatable kstp
integer, dimension(:), allocatable knew
integer, dimension(:), allocatable nfm1
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable lfinp
integer, dimension(:), allocatable lbinp
integer, dimension(:), allocatable nf
integer, dimension(:), allocatable nfm3
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nfp1
integer, dimension(:), allocatable lnew
integer, dimension(:), allocatable krhs
integer, dimension(:), allocatable lfout
integer, dimension(:), allocatable nstp
subroutine, public my25_corstep(ng, tile)
subroutine, public my25_prestep(ng, tile)
subroutine, public nesting(ng, model, isection)
Definition nesting.F:140
subroutine, public obc_adjust(ng, tile, linp)
Definition obc_adjust.F:36
subroutine, public load_obc(ng, tile, lout)
Definition obc_adjust.F:503
subroutine, public omega(ng, tile, model)
Definition omega.F:42
subroutine, public post_initial(ng, model)
subroutine, public rho_eos(ng, tile, model)
Definition rho_eos.F:48
subroutine, public rhs3d(ng, tile)
Definition rhs3d.F:26
subroutine, public sediment(ng, tile)
Definition sediment.F:78
subroutine, public set_avg(ng, tile)
Definition set_avg.F:34
subroutine, public set_depth(ng, tile, model)
Definition set_depth.F:34
subroutine, public set_massflux(ng, tile, model)
subroutine, public set_tides(ng, tile)
Definition set_tides.F:26
subroutine, public set_vbc(ng, tile)
Definition set_vbc.F:28
subroutine, public set_zeta(ng, tile)
Definition set_zeta.F:26
subroutine, public step2d(ng, tile)
Definition step2d_FB.h:75
subroutine, public step3d_t(ng, tile)
Definition step3d_t.F:41
subroutine, public step3d_uv(ng, tile)
Definition step3d_uv.F:50
subroutine, public step_floats(ng, lstr, lend)
Definition step_floats.F:43
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52
subroutine, public wvelocity(ng, tile, ninp)
Definition wvelocity.F:27
subroutine ntimesteps(model, runinterval, nl, nsteps, rsteps)
Definition ntimestep.F:3
subroutine output(ng)
Definition output.F:4
subroutine set_data(ng, tile)
Definition set_data.F:4
subroutine set_diags(ng, tile)
Definition set_diags.F:4

References analytical_mod::ana_vmix(), bbl_mod::bblm(), biology_mod::biology(), bulk_flux_mod::bulk_flux(), bvf_mix_mod::bvf_mix(), mod_scalars::compositegrid, mod_coupler::couplesteps, diag_mod::diag(), mod_nesting::donortofiner, mod_scalars::dt, mod_scalars::exit_flag, mod_parallel::first_tile, forcing_mod::forcing(), strings_mod::founderror(), frc_adjust_mod::frc_adjust(), frc_iau_mod::frc_iau(), mod_scalars::frequentimpulse, get_data(), gls_corstep_mod::gls_corstep(), gls_prestep_mod::gls_prestep(), mod_param::gridnumber, mod_param::gridsinlayer, mod_coupler::iatmos, mod_scalars::iauswitch, mod_scalars::iic, mod_scalars::iif, mod_scalars::indx1, mod_param::inlm, mod_coupler::iwaves, mod_scalars::jic, mod_stepping::knew, mod_stepping::krhs, mod_stepping::kstp, mod_parallel::last_tile, mod_stepping::lbinp, mod_stepping::lbout, mod_stepping::lfinp, mod_scalars::lfloats, mod_stepping::lfout, lmd_vmix_mod::lmd_vmix(), mod_stepping::lnew, frc_adjust_mod::load_frc(), obc_adjust_mod::load_obc(), mod_parallel::master, my25_corstep_mod::my25_corstep(), my25_prestep_mod::my25_prestep(), mod_parallel::mythread, mod_nesting::n2dcs, mod_nesting::n2dfx, mod_nesting::n2dps, mod_nesting::n2way, mod_nesting::n3dtv, mod_nesting::n3duv, mod_scalars::nadj, mod_nesting::nbstr, nesting_mod::nesting(), mod_param::nestlayers, mod_scalars::next_kstp, mod_stepping::nf, mod_scalars::nfast, mod_param::nfloats, mod_stepping::nfm1, mod_stepping::nfm2, mod_stepping::nfm3, mod_stepping::nfp1, mod_param::nft, mod_nesting::ngetd, mod_param::ngrids, dotproduct_mod::nl_dotproduct(), mod_nesting::nmask, mod_nesting::nmflx, mod_stepping::nnew, mod_scalars::noerror, mod_nesting::nputd, mod_stepping::nrhs, mod_nesting::nrhst, mod_scalars::nrun, mod_stepping::nstp, mod_scalars::ntend, ntimesteps(), mod_scalars::ntstart, mod_parallel::numthreads, mod_nesting::nzeta, obc_adjust_mod::obc_adjust(), omega_mod::omega(), output(), post_initial_mod::post_initial(), mod_scalars::predictor_2d_step, mod_scalars::processinputdata, mod_scalars::refinedgrid, mod_scalars::refinescale, rho_eos_mod::rho_eos(), rhs3d_mod::rhs3d(), mod_scalars::sec2day, sediment_mod::sediment(), set_avg_mod::set_avg(), set_data(), set_depth_mod::set_depth(), set_diags(), set_massflux_mod::set_massflux(), set_tides_mod::set_tides(), set_vbc_mod::set_vbc(), set_zeta_mod::set_zeta(), mod_iounits::stdout, step2d_mod::step2d(), step3d_t_mod::step3d_t(), step3d_uv_mod::step3d_uv(), mod_scalars::step_counter, step_floats_mod::step_floats(), mod_scalars::tdays, mod_scalars::time, mod_scalars::time4jedi, mod_scalars::time_code, dateclock_mod::time_string(), mod_scalars::timeiau, mod_nesting::twowayinterval, and wvelocity_mod::wvelocity().

Referenced by rbl4dvar_mod::analysis(), i4dvar_mod::background(), r4dvar_mod::background(), rbl4dvar_mod::background(), i4dvar_mod::posterior_analysis(), and roms_kernel_mod::roms_run().

Here is the call graph for this function:
Here is the caller graph for this function: