ROMS
Loading...
Searching...
No Matches
main2d.F
Go to the documentation of this file.
1#include "cppdefs.h"
2#if defined NONLINEAR && !defined SOLVE3D
3 SUBROUTINE main2d (RunInterval)
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 configured as a 2D barotropic shallow water ocean model. It !
14! advances forward the NLM for all nested grids, if any, by the !
15! specified 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 USE dateclock_mod, ONLY : time_string
40 USE diag_mod, ONLY : diag
41# ifdef TLM_CHECK
43# endif
44# ifdef TIDE_GENERATING_FORCES
45 USE equilibrium_tide_mod, ONLY : equilibrium_tide
46# endif
47# if defined NLM_OUTER || \
48 defined rbl4dvar || \
49 defined rbl4dvar_ana_sensitivity || \
50 defined rbl4dvar_fct_sensitivity || \
51 defined sp4dvar
52 USE forcing_mod, ONLY : forcing
53# endif
54# ifdef ADJUST_WSTRESS
56# endif
58# if defined ATM_COUPLING && defined MCT_LIB
59 USE mct_coupler_mod, ONLY : ocn2atm_coupling
60# endif
61# if defined WAV_COUPLING && defined MCT_LIB
62 USE mct_coupler_mod, ONLY : ocn2wav_coupling
63# endif
64# ifdef NESTING
65 USE nesting_mod, ONLY : nesting
66# ifndef ONE_WAY
67 USE nesting_mod, ONLY : do_twoway
68# endif
69# endif
70# if defined ADJUST_BOUNDARY
72# endif
73# ifdef AVERAGES
74 USE set_avg_mod, ONLY : set_avg
75# endif
76# if defined SSH_TIDES || defined UV_TIDES
77 USE set_tides_mod, ONLY : set_tides
78# endif
79 USE set_vbc_mod, ONLY : set_vbc
80 USE step2d_mod, ONLY : step2d
81# ifdef FLOATS
83# endif
84 USE strings_mod, ONLY : founderror
85!
86 implicit none
87!
88! Imported variable declarations.
89!
90 real(dp), intent(in) :: RunInterval
91!
92! Local variable declarations.
93!
94 logical :: DoNestLayer, Time_Step
95!
96 integer :: Nsteps, Rsteps
97 integer :: ig, il, istep, ng, nl, tile
98 integer :: next_indx1
99# ifdef FLOATS
100 integer :: Lend, Lstr, chunk_size
101# endif
102!
103 character (len=*), parameter :: MyFile = &
104 & __FILE__
105!
106!=======================================================================
107! Time-step nonlinear vertically integrated equations.
108!=======================================================================
109!
110! Time-step the 3D kernel for the specified time interval (seconds),
111! RunInterval.
112!
113 time_step=.true.
114 donestlayer=.true.
115!
116 kernel_loop : DO WHILE (time_step)
117!
118! In nesting applications, the number of nesting layers (NestLayers) is
119! used to facilitate refinement grids and composite/refinament grids
120! combinations. Otherwise, the solution it is looped once for a single
121! grid application (NestLayers = 1).
122!
123 nl=0
124# ifdef NESTING
125 twowayinterval(1:ngrids)=0.0_r8
126# endif
127!
128 nest_layer : DO WHILE (donestlayer)
129!
130! Determine number of time steps to compute in each nested grid layer
131! based on the specified time interval (seconds), RunInterval. Non
132! nesting applications have NestLayers=1. Notice that RunInterval is
133! set in the calling driver. Its value may span the full period of the
134! simulation, a multi-model coupling interval (RunInterval > ifac*dt),
135! or just a single step (RunInterval=0).
136!
137 CALL ntimesteps (inlm, runinterval, nl, nsteps, rsteps)
138 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
139 IF ((nl.le.0).or.(nl.gt.nestlayers)) EXIT
140!
141! Time-step governing equations for Nsteps.
142!
143 step_loop : DO istep=1,nsteps
144!
145! Set time indices and time clock.
146!
147 DO ig=1,gridsinlayer(nl)
148 ng=gridnumber(ig,nl)
149 tdays(ng)=time(ng)*sec2day
150 IF (step_counter(ng).eq.rsteps) time_step=.false.
151 END DO
152!
153!-----------------------------------------------------------------------
154! Read in required data, if any, from input NetCDF files.
155!-----------------------------------------------------------------------
156!
157 DO ig=1,gridsinlayer(nl)
158 ng=gridnumber(ig,nl)
159!$OMP MASTER
160 CALL get_data (ng)
161!$OMP END MASTER
162!$OMP BARRIER
164 & __line__, myfile)) RETURN
165 END DO
166!
167!-----------------------------------------------------------------------
168! If applicable, process input data: time interpolate between data
169! snapshots.
170!-----------------------------------------------------------------------
171!
172 DO ig=1,gridsinlayer(nl)
173 ng=gridnumber(ig,nl)
174 DO tile=first_tile(ng),last_tile(ng),+1
175 CALL set_data (ng, tile)
176 END DO
177!$OMP BARRIER
178 END DO
179 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
180
181
182# if defined NLM_OUTER || \
183 defined rbl4dvar || \
184 defined rbl4dvar_ana_sensitivity || \
185 defined rbl4dvar_fct_sensitivity
186!
187!-----------------------------------------------------------------------
188! If appropriate, add convolved adjoint solution impulse forcing to
189! the nonlinear model solution. Notice that the forcing is only needed
190! after finishing all the inner loops. The forcing is continuous.
191! That is, it is time interpolated at every time-step from available
192! snapshots (FrequentImpulse=TRUE).
193!-----------------------------------------------------------------------
194!
195 DO ig=1,gridsinlayer(nl)
196 ng=gridnumber(ig,nl)
197# if defined WEAK_NOINTERP
198 IF ((iic(ng).gt.1).and.(iic(ng).ne.ntend(ng)+1).and. &
199 & (mod(iic(ng)-1,nadj(ng)).eq.0)) THEN
200 IF (master) THEN
201 WRITE (stdout,*) ' FORCING NLM at iic = ',iic(ng)
202 END IF
203# endif
204 IF (frequentimpulse(ng)) THEN
205 DO tile=first_tile(ng),last_tile(ng),+1
206 CALL forcing (ng, tile, kstp(ng), nstp(ng))
207 END DO
208!$OMP BARRIER
209 END IF
210# ifdef WEAK_NOINTERP
211 END IF
212# endif
213 END DO
214# endif
215!
216!-----------------------------------------------------------------------
217! Initialize all time levels and compute other initial fields.
218!-----------------------------------------------------------------------
219!
220 DO ig=1,gridsinlayer(nl)
221 ng=gridnumber(ig,nl)
222 IF (iic(ng).eq.ntstart(ng)) THEN
223!
224! Initialize free-surface.
225!
226 DO tile=first_tile(ng),last_tile(ng),+1
227 CALL ini_zeta (ng, tile, inlm)
228 END DO
229!$OMP BARRIER
230!
231! Initialize other state variables.
232!
233 DO tile=last_tile(ng),first_tile(ng),-1
234 CALL ini_fields (ng, tile, inlm)
235 END DO
236!$OMP BARRIER
237
238# ifdef NESTING
239!
240! Extract donor grid initial data at contact points and store it in
241! REFINED structure so it can be used for the space-time interpolation.
242!
243 IF (refinedgrid(ng)) THEN
244 CALL nesting (ng, inlm, ngetd)
245 END IF
246# endif
247 END IF
248 END DO
249!
250!-----------------------------------------------------------------------
251! Compute and report diagnostics. If appropriate, accumulate time-
252! averaged output data which needs a irreversible loop in shared-memory
253! jobs.
254!-----------------------------------------------------------------------
255!
256 DO ig=1,gridsinlayer(nl)
257 ng=gridnumber(ig,nl)
258 DO tile=first_tile(ng),last_tile(ng),+1 ! irreversible
259# ifdef AVERAGES
260 CALL set_avg (ng, tile)
261# endif
262# ifdef DIAGNOSTICS
263 CALL set_diags (ng, tile)
264# endif
265# ifdef TIDE_GENERATING_FORCES
266 CALL equilibrium_tide (ng, tile, inlm)
267# endif
268 CALL diag (ng, tile)
269# ifdef TLM_CHECK
270 CALL nl_dotproduct (ng, tile, lnew(ng))
271# endif
272 END DO
273!$OMP BARRIER
274 END DO
275 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
276
277# if defined ATM_COUPLING && defined MCT_LIB
278!
279!-----------------------------------------------------------------------
280! Couple ocean to atmosphere model every "CoupleSteps(Iatmos)"
281! timesteps: get air/sea fluxes.
282!-----------------------------------------------------------------------
283!
284 DO ig=1,gridsinlayer(nl)
285 ng=gridnumber(ig,nl)
286 IF ((iic(ng).ne.ntstart(ng)).and. &
287 & mod(iic(ng),couplesteps(iatmos,ng)).eq.0) THEN
288 DO tile=last_tile(ng),first_tile(ng),-1
289 CALL ocn2atm_coupling (ng, tile)
290 END DO
291!$OMP BARRIER
292 END IF
293 END DO
294# endif
295
296# ifdef ADJUST_BOUNDARY
297!
298!-----------------------------------------------------------------------
299! Interpolate open boundary increments and adjust open boundary.
300! Load open boundary into storage arrays. Skip the last output
301! timestep.
302!-----------------------------------------------------------------------
303!
304 DO ig=1,gridsinlayer(nl)
305 ng=gridnumber(ig,nl)
306 IF (iic(ng).lt.(ntend(ng)+1)) THEN
307 DO tile=first_tile(ng),last_tile(ng),+1
308 CALL obc_adjust (ng, tile, lbinp(ng))
309 CALL load_obc (ng, tile, lbout(ng))
310 END DO
311!$OMP BARRIER
312 END IF
313 END DO
314# endif
315
316# ifdef ADJUST_WSTRESS
317!
318!-----------------------------------------------------------------------
319! Interpolate surface forcing increments and adjust surface forcing.
320! Load surface forcing to storage arrays.
321!-----------------------------------------------------------------------
322!
323 DO ig=1,gridsinlayer(nl)
324 ng=gridnumber(ig,nl)
325 IF (iic(ng).lt.(ntend(ng)+1)) THEN
326 DO tile=first_tile(ng),last_tile(ng),+1
327 CALL frc_adjust (ng, tile, lfinp(ng))
328 CALL load_frc (ng, tile, lfout(ng))
329 END DO
330!$OMP BARRIER
331 END IF
332 END DO
333# endif
334
335# if defined WAV_COUPLING && defined MCT_LIB
336!
337!-----------------------------------------------------------------------
338! Couple ocean to waves model every "CoupleSteps(Iwaves)"
339! timesteps: get waves/sea fluxes.
340!-----------------------------------------------------------------------
341!
342 DO ig=1,gridsinlayer(nl)
343 ng=gridnumber(ig,nl)
344 IF ((iic(ng).ne.ntstart(ng)).and. &
345 & mod(iic(ng)-1,couplesteps(iwaves,ng)).eq.0) THEN
346 DO tile=first_tile(ng),last_tile(ng),+1
347 CALL ocn2wav_coupling (ng, tile)
348 END DO
349!$OMP BARRIER
350 END IF
351 END DO
352# endif
353
354!
355!-----------------------------------------------------------------------
356! Set vertical boundary conditions. Process tidal forcing.
357!-----------------------------------------------------------------------
358!
359 DO ig=1,gridsinlayer(nl)
360 ng=gridnumber(ig,nl)
361 DO tile=first_tile(ng),last_tile(ng),+1
362 CALL set_vbc (ng, tile)
363# if defined SSH_TIDES || defined UV_TIDES
364 CALL set_tides (ng, tile)
365# endif
366 END DO
367!$OMP BARRIER
368 END DO
369
370# ifdef NESTING
371!
372! If composite or mosaic grids, process additional points in the
373! contact zone between connected grids for bottom stress variables.
374!
375 DO ig=1,gridsinlayer(nl)
376 ng=gridnumber(ig,nl)
377 IF (any(compositegrid(:,ng))) THEN
378 CALL nesting (ng, inlm, nbstr)
379 END IF
380 END DO
381# endif
382!
383!-----------------------------------------------------------------------
384! If appropriate, write out fields into output NetCDF files. Notice
385! that IO data is written in delayed and serial mode. Exit if last
386! time step.
387!-----------------------------------------------------------------------
388!
389 DO ig=1,gridsinlayer(nl)
390 ng=gridnumber(ig,nl)
391!$OMP MASTER
392 CALL output (ng)
393!$OMP END MASTER
394!$OMP BARRIER
395 IF ((founderror(exit_flag, noerror, __line__, myfile)).or.&
396 & ((iic(ng).eq.(ntend(ng)+1)).and.(ng.eq.ngrids))) THEN
397 RETURN
398 END IF
399 END DO
400
401# ifdef NESTING
402!
403!-----------------------------------------------------------------------
404! If refinement grid, interpolate (space, time) state variables
405! contact points from donor grid extracted data.
406# ifdef NESTING_DEBUG
407!
408! Also, fill BRY_CONTACT(:,:)%Mflux to check for mass conservation
409! between coarse and fine grids. This is only done for diagnostic
410! purposes. Also, debugging is possible with very verbose output
411! to fort.300 is allowed by activating uppercase(nesting_debug).
412# endif
413!-----------------------------------------------------------------------
414!
415 DO ig=1,gridsinlayer(nl)
416 ng=gridnumber(ig,nl)
417 IF (refinedgrid(ng).and.(refinescale(ng).gt.0)) THEN
418 CALL nesting (ng, inlm, nputd)
419# ifdef NESTING_DEBUG
420 CALL nesting (ng, inlm, nmflx)
421# endif
422 END IF
423 END DO
424# endif
425
426# ifdef STEP2D_FB_AB3_AM4
427!
428!-----------------------------------------------------------------------
429! Solve nonlinear vertically integrated primitive equations for
430! free-surface and momentum components using a generalized Forward-
431! Backward, 3rd-order Adams-Bashforth / 4th-order Adams-Moulton
432! (FB AB3-AM4) time stepping scheme (Shchepetkin and McWilliams,
433! 2009).
434!-----------------------------------------------------------------------
435!
436 DO ig=1,gridsinlayer(nl)
437 ng=gridnumber(ig,nl)
438 iif(ng)=1
439 nfast(ng)=1
440 kstp(ng)=knew(ng)
441 knew(ng)=kstp(ng)+1
442 IF (knew(ng).gt.4) knew(ng)=1
443 IF (mod(knew(ng),2).eq.0) THEN ! zig-zag
444 DO tile=first_tile(ng),last_tile(ng),+1 ! processing
445 CALL step2d (ng, tile) ! sequence
446 END DO
447!$OMP BARRIER
448 ELSE
449 DO tile=last_tile(ng),first_tile(ng),-1
450 CALL step2d (ng, tile)
451 END DO
452!$OMP BARRIER
453 END IF
454 END DO
455
456# else
457
458# ifdef STEP2D_FB_LF_AM3
459!
460!-----------------------------------------------------------------------
461! Solve nonlinear vertically integrated primitive equations for
462! free-surface and momentum components using a predictor-corrector
463! LeapFrog / 3rd-order Adams-Moulton with a Forward-Backward
464! feeback (FB LF-AM3) stepping scheme (Shchepetkin and McWilliams,
465! 2009).
466!-----------------------------------------------------------------------
467!
468! Predictor LF substep with FB-feedback.
469!
470 DO ig=1,gridsinlayer(nl)
471 ng=gridnumber(ig,nl)
472 iif(ng)=1
473 nfast(ng)=1
474 kstp(ng)=next_kstp(ng)
475 knew(ng)=3
476
477 DO tile=last_tile(ng),first_tile(ng),-1
478 CALL step2d (ng, tile)
479 END DO
480!$OMP BARRIER
481 END DO
482
483# ifdef NESTING
484!
485! If composite or mosaic grids, process additional points in the
486! contact zone between connected grids for the state variables
487! associated with the 2D engine PREDICTOR STEP section.
488# ifdef NESTING_DEBUG
489! If refinement, check mass flux conservation between coarse and
490! fine grids during debugging.
491# endif
492!
493 DO ig=1,gridsinlayer(nl)
494 ng=gridnumber(ig,nl)
495 IF (any(compositegrid(:,ng))) THEN
496 CALL nesting (ng, inlm, n2dps)
497 END IF
498# ifdef NESTING_DEBUG
499 IF (refinedgrid(ng).and.(refinescale(ng).gt.0)) THEN
500 CALL nesting (ng, inlm, nmflx)
501 END IF
502# endif
503 END DO
504# endif
505!
506! Corrector AM3 substep with FB-feedback.
507!
508 DO ig=1,gridsinlayer(nl)
509 ng=gridnumber(ig,nl)
510 knew(ng)=3-kstp(ng)
511 next_kstp(ng)=knew(ng)
512
513 DO tile=first_tile(ng),last_tile(ng),+1
514 CALL step2d (ng, tile)
515 END DO
516!$OMP BARRIER
517 END DO
518
519# else
520!
521!-----------------------------------------------------------------------
522! Solve nonlinear vertically integrated primitive equations for
523! free-surface and momentum components using a predictor-corrector
524! LeapFrog with 3rd-order Adams-Moulton (LF-AM3) time stepping scheme
525! (ROMS legacy 2D kernel).
526!-----------------------------------------------------------------------
527!
528! Set time indices for predictor step. The PREDICTOR_2D_STEP switch
529! it is assumed to be false before the first time-step.
530!
531 DO ig=1,gridsinlayer(nl)
532 ng=gridnumber(ig,nl)
533 iif(ng)=1
534 nfast(ng)=1
535 next_indx1=3-indx1(ng)
536 IF (.not.predictor_2d_step(ng)) THEN
537 predictor_2d_step(ng)=.true.
538 IF (first_2d_step) THEN
539 kstp(ng)=indx1(ng)
540 ELSE
541 kstp(ng)=3-indx1(ng)
542 END IF
543 knew(ng)=3
544 krhs(ng)=indx1(ng)
545 END IF
546!
547! Predictor step - Advance barotropic equations using 2D time-step
548! ============== predictor scheme.
549!
550 DO tile=last_tile(ng),first_tile(ng),-1
551 CALL step2d (ng, tile)
552 END DO
553!$OMP BARRIER
554 END DO
555
556# ifdef NESTING
557!
558! If composite or mosaic grids, process additional points in the
559! contact zone between connected grids for the state variables
560! associated with the 2D engine PREDICTOR STEP section.
561# ifdef NESTING_DEBUG
562! If refinement, check mass flux conservation between coarse and
563! fine grids during debugging.
564! Warning: very verbose output to fort.300 ascii file to check
565! mass flux conservation.
566# endif
567!
568 DO ig=1,gridsinlayer(nl)
569 ng=gridnumber(ig,nl)
570 IF (any(compositegrid(:,ng))) THEN
571 CALL nesting (ng, inlm, n2dps)
572 END IF
573# ifdef NESTING_DEBUG
574 IF (refinedgrid(ng).and.(refinescale(ng).gt.0)) THEN
575 CALL nesting (ng, inlm, nmflx)
576 END IF
577# endif
578 END DO
579# endif
580!
581! Set time indices for corrector step.
582!
583 DO ig=1,gridsinlayer(nl)
584 ng=gridnumber(ig,nl)
585 IF (predictor_2d_step(ng)) THEN
586 predictor_2d_step(ng)=.false.
587 knew(ng)=next_indx1
588 kstp(ng)=3-knew(ng)
589 krhs(ng)=3
590 IF (iif(ng).lt.(nfast(ng)+1)) indx1(ng)=next_indx1
591 END IF
592!
593! Corrector step - Apply 2D time-step corrector scheme. Notice that
594! ============== there is not need for a corrector step during the
595! auxiliary (nfast+1) time-step.
596!
597 IF (iif(ng).lt.(nfast(ng)+1)) THEN
598 DO tile=first_tile(ng),last_tile(ng),+1
599 CALL step2d (ng, tile)
600 END DO
601!$OMP BARRIER
602 END IF
603 END DO
604# endif
605
606# endif
607
608# ifdef NESTING
609# if defined MASKING && defined WET_DRY
610!
611!-----------------------------------------------------------------------
612! If nesting and wetting and drying, scale horizontal interpolation
613! weights to account for land/sea masking in contact areas. This needs
614! to be done at very time-step since the Land/Sea masking is time
615! dependent.
616!-----------------------------------------------------------------------
617!
618 DO ig=1,gridsinlayer(nl)
619 ng=gridnumber(ig,nl)
620 CALL nesting (ng, inlm, nmask)
621 END DO
622# endif
623!
624!-----------------------------------------------------------------------
625! If composite or mosaic grids, process additional points in the
626! contact zone between connected grids for the state variables
627! associated with the 2D engine CORRECTOR STEP section (KNEW INDEX).
628# ifdef NESTING_DEBUG
629! If refinement, check mass flux conservation between coarse and
630! fine grids during debugging.
631! Warning: very verbose output to fort.300 ascii file to check
632! mass flux conservation.
633# endif
634!-----------------------------------------------------------------------
635!
636 DO ig=1,gridsinlayer(nl)
637 ng=gridnumber(ig,nl)
638 IF (any(compositegrid(:,ng))) THEN
639 CALL nesting (ng, inlm, n2dcs)
640 END IF
641# ifdef NESTING_DEBUG
642 IF (refinedgrid(ng).and.(refinescale(ng).gt.0)) THEN
643 CALL nesting (ng, inlm, nmflx)
644 END IF
645# endif
646 END DO
647# endif
648
649# ifdef NESTING
650# ifndef ONE_WAY
651!
652!-----------------------------------------------------------------------
653! If refinement grids, perform two-way coupling between fine and
654! coarse grids. Correct coarse grid tracers values at the refinement
655! grid with refined accumulated fluxes. Then, replace coarse grid
656! state variable with averaged refined grid values (two-way nesting).
657! Update coarse grid depth variables.
658!
659! The two-way exchange of infomation between nested grids needs to be
660! done at the correct time-step and in the right sequence.
661!-----------------------------------------------------------------------
662!
663 DO il=nestlayers,1,-1
664 DO ig=1,gridsinlayer(il)
665 ng=gridnumber(ig,il)
666 IF (do_twoway(inlm, nl, il, ng, istep)) THEN
667 CALL nesting (ng, inlm, n2way)
668 END IF
669 END DO
670 END DO
671# endif
672!
673!-----------------------------------------------------------------------
674! If donor to a finer grid, extract data for the external contact
675! points. This is the latest solution for the coarser grid.
676!
677! It is stored in the REFINED structure so it can be used for the
678! space-time interpolation.
679!-----------------------------------------------------------------------
680!
681 DO ig=1,gridsinlayer(nl)
682 ng=gridnumber(ig,nl)
683 IF (donortofiner(ng)) THEN
684 CALL nesting (ng, inlm, ngetd)
685 END IF
686 END DO
687# endif
688
689# ifdef FLOATS
690!
691!-----------------------------------------------------------------------
692! Compute Lagrangian drifters trajectories: Split all the drifters
693! between all the computational threads, except in distributed-memory
694! and serial configurations. In distributed-memory, the parallel node
695! containing the drifter is selected internally since the state
696! variables do not have a global scope.
697!-----------------------------------------------------------------------
698!
699 DO ig=1,gridsinlayer(nl)
700 ng=gridnumber(ig,nl)
701 IF (lfloats(ng)) THEN
702# ifdef _OPENMP
703 chunk_size=(nfloats(ng)+numthreads-1)/numthreads
704 lstr=1+mythread*chunk_size
705 lend=min(nfloats(ng),lstr+chunk_size-1)
706# else
707 lstr=1
708 lend=nfloats(ng)
709# endif
710 CALL step_floats (ng, lstr, lend)
711!$OMP BARRIER
712!
713! Shift floats time indices.
714!
715 nfp1(ng)=mod(nfp1(ng)+1,nft+1)
716 nf(ng)=mod(nf(ng)+1,nft+1)
717 nfm1(ng)=mod(nfm1(ng)+1,nft+1)
718 nfm2(ng)=mod(nfm2(ng)+1,nft+1)
719 nfm3(ng)=mod(nfm3(ng)+1,nft+1)
720 END IF
721 END DO
722# endif
723!
724!-----------------------------------------------------------------------
725! Advance time index and time clock.
726!-----------------------------------------------------------------------
727!
728 DO ig=1,gridsinlayer(nl)
729 ng=gridnumber(ig,nl)
730 iic(ng)=iic(ng)+1
731 time(ng)=time(ng)+dt(ng)
732 step_counter(ng)=step_counter(ng)-1
733 CALL time_string (time(ng), time_code(ng))
734 END DO
735
736 END DO step_loop
737
738 END DO nest_layer
739
740 END DO kernel_loop
741!
742 RETURN
743 END SUBROUTINE main2d
744#else
745 SUBROUTINE main2d
746 RETURN
747 END SUBROUTINE main2d
748#endif
subroutine get_data(ng)
Definition get_data.F:3
subroutine main2d
Definition main2d.F:746
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 ini_zeta(ng, tile, model)
Definition ini_fields.F:709
subroutine, public ini_fields(ng, tile, model)
Definition ini_fields.F:67
integer, dimension(:,:), allocatable couplesteps
integer iatmos
integer iwaves
integer stdout
integer, parameter ngetd
Definition mod_nesting.F:77
integer, parameter nputd
Definition mod_nesting.F:79
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 n2dps
Definition mod_nesting.F:90
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
integer, dimension(:), allocatable next_kstp
logical, dimension(:), allocatable predictor_2d_step
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
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, dimension(:), allocatable step_counter
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 lfinp
integer, dimension(:), allocatable lbinp
integer, dimension(:), allocatable nf
integer, dimension(:), allocatable nfm3
integer, dimension(:), allocatable nfp1
integer, dimension(:), allocatable lnew
integer, dimension(:), allocatable krhs
integer, dimension(:), allocatable lfout
integer, dimension(:), allocatable nstp
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 set_avg(ng, tile)
Definition set_avg.F:34
subroutine, public set_tides(ng, tile)
Definition set_tides.F:26
subroutine, public set_vbc(ng, tile)
Definition set_vbc.F:28
subroutine, public step2d(ng, tile)
Definition step2d_FB.h:75
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 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