ROMS
Loading...
Searching...
No Matches
rp_main2d.F
Go to the documentation of this file.
1#include "cppdefs.h"
2#if defined TL_IOMS && !defined SOLVE3D
3 SUBROUTINE rp_main2d (RunInterval)
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This routine is the main driver for representers tangent linear !
13! ROMS when configure as shallow water (barotropic ) ocean !
14! model only. It advances advances forward the representer model !
15! for all nested grids, if any, by the specified time interval !
16! (seconds), RunInterval. !
17! !
18!=======================================================================
19!
20 USE mod_param
21 USE mod_parallel
22# if defined MODEL_COUPLING && defined MCT_LIB
23 USE mod_coupler
24# endif
25 USE mod_iounits
26 USE mod_scalars
27 USE mod_stepping
28!
29 USE dateclock_mod, ONLY : time_string
30# ifdef TIDE_GENERATING_FORCES
31 USE equilibrium_tide_mod, ONLY : equilibrium_tide
32# endif
33# if defined ATM_COUPLING_NOT_YET && defined MCT_LIB
34 USE mct_coupler_mod, ONLY : ocn2atm_coupling
35# endif
36# if defined WAV_COUPLING_NOT_YET && defined MCT_LIB
37 USE mct_coupler_mod, ONLY : ocn2wav_coupling
38# endif
39 USE strings_mod, ONLY : founderror
40 USE rp_diag_mod, ONLY : rp_diag
41# ifdef defined ADJUST_WSTRESS
43# endif
45# ifdef ADJUST_BOUNDARY
47# endif
48# ifdef NEARSHORE_MELLOR_NOT_YET
49!! USE rp_radiation_stress_mod, ONLY : rp_radiation_stress
50# endif
51# if defined SSH_TIDES_NOT_YET || defined UV_TIDES_NOT_YET
52!! USE rp_set_tides_mod, ONLY : rp_set_tides
53# endif
54 USE rp_set_vbc_mod, ONLY : rp_set_vbc
55 USE rp_step2d_mod, ONLY : rp_step2d
56# ifdef FLOATS_NOT_YET
57!! USE rp_step_floats_mod, ONLY : rp_step_floats
58# endif
59# ifdef WEAK_CONSTRAINT
60 USE tl_forcing_mod, ONLY : tl_forcing
61# endif
62# ifdef RP_AVERAGES
63 USE rp_set_avg_mod, ONLY : tl_set_avg
64# endif
65# if defined PROPAGATOR || \
66 (defined masking && (defined read_water || defined write_water))
67 USE wpoints_mod, ONLY : wpoints
68# endif
69!
70 implicit none
71!
72! Imported variable declarations.
73!
74 real(dp), intent(in) :: RunInterval
75!
76! Local variable declarations.
77!
78 integer :: ng, tile
79 integer :: next_indx1
80# ifdef FLOATS_NOT_YET
81 integer :: Lend, Lstr, chunk_size
82# endif
83!
84 real(r8) :: MaxDT, my_StepTime
85!
86 character (len=*), parameter :: MyFile = &
87 & __FILE__
88!
89!=======================================================================
90! Time-step tangent linear vertically integrated equations.
91!=======================================================================
92!
93 my_steptime=0.0_r8
94 maxdt=maxval(dt)
95
96 step_loop : DO WHILE (my_steptime.le.(runinterval+0.5_r8*maxdt))
97
98 my_steptime=my_steptime+maxdt
99!
100! Set time clock.
101!
102 DO ng=1,ngrids
103 iic(ng)=iic(ng)+1
104!$OMP MASTER
105 time(ng)=time(ng)+dt(ng)
106 tdays(ng)=time(ng)*sec2day
107 CALL time_string (time(ng), time_code(ng))
108!$OMP END MASTER
109 END DO
110!$OMP BARRIER
111!
112!-----------------------------------------------------------------------
113! Read in required data, if any, from input NetCDF files.
114!-----------------------------------------------------------------------
115!
116 DO ng=1,ngrids
117!$OMP MASTER
118 CALL rp_get_data (ng)
119!$OMP END MASTER
120!$OMP BARRIER
121 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
122 END DO
123!
124!-----------------------------------------------------------------------
125! If applicable, process input data: time interpolate between data
126! snapshots.
127!-----------------------------------------------------------------------
128!
129 DO ng=1,ngrids
130 DO tile=first_tile(ng),last_tile(ng),+1
131 CALL rp_set_data (ng, tile)
132 END DO
133!$OMP BARRIER
134 END DO
135 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
136
137# ifdef WEAK_CONSTRAINT
138!
139!-----------------------------------------------------------------------
140! If appropriate, add convolved adjoint solution impulse forcing to
141! the representer model solution. Notice that the forcing is only
142! needed after finishing all inner loops. The forcing is continuous.
143! That is, it is time interpolated at every time-step from available
144! snapshots (FrequentImpulse=TRUE).
145!-----------------------------------------------------------------------
146!
147 DO ng=1,ngrids
148 IF (frequentimpulse(ng)) THEN
149 DO tile=first_tile(ng),last_tile(ng),+1
150 CALL tl_forcing (ng, tile, kstp(ng), nstp(ng))
151 END DO
152!$OMP BARRIER
153 END IF
154 END DO
155# endif
156!
157!-----------------------------------------------------------------------
158! If not a restart, initialize all time levels and compute other
159! initial fields.
160!-----------------------------------------------------------------------
161!
162 DO ng=1,ngrids
163 IF (iic(ng).eq.ntstart(ng)) THEN
164!
165! Initialize free-surface.
166!
167 DO tile=first_tile(ng),last_tile(ng),+1
168 CALL rp_ini_zeta (ng, tile, irpm)
169 END DO
170!$OMP BARRIER
171!
172! Initialize other state variables.
173!
174 DO tile=last_tile(ng),first_tile(ng),-1
175 CALL rp_ini_fields (ng, tile, irpm)
176 END DO
177!$OMP BARRIER
178 END IF
179 END DO
180!
181!-----------------------------------------------------------------------
182! Compute and report diagnostics. If appropriate, accumulate time-
183! averaged output data which needs a irreversible loop in shared-memory
184! jobs.
185!-----------------------------------------------------------------------
186!
187 DO ng=1,ngrids
188 DO tile=first_tile(ng),last_tile(ng),+1 ! irreversible
189# ifdef RP_AVERAGES
190 CALL tl_set_avg (ng, tile)
191# endif
192# ifdef DIAGNOSTICS
193!! CALL rp_set_diags (ng, tile)
194# endif
195# ifdef TIDE_GENERATING_FORCES
196 CALL equilibrium_tide (ng, tile, irpm)
197# endif
198 CALL rp_diag (ng, tile)
199 END DO
200!$OMP BARRIER
201 END DO
202
203# if defined ATM_COUPLING_NOT_YET && defined MCT_LIB
204!
205!-----------------------------------------------------------------------
206! Couple to atmospheric model every CoupleSteps(Iatmos) timesteps: get
207! air/sea fluxes.
208!-----------------------------------------------------------------------
209!
210 DO ng=1,ngrids
211 IF ((iic(ng).ne.ntstart(ng)).and. &
212 & mod(iic(ng)-1,couplesteps(iatmos,ng)).eq.0) THEN
213 DO tile=last_tile(ng),first_tile(ng),-1
214 CALL ocn2atm_coupling (ng, tile)
215 END DO
216!$OMP BARRIER
217 END IF
218 END DO
219# endif
220
221# if defined WAV_COUPLING_NOT_YET && defined MCT_LIB
222!
223!-----------------------------------------------------------------------
224! Couple to waves model every CoupleSteps(Iwaves) timesteps: get
225! waves/sea fluxes.
226!-----------------------------------------------------------------------
227!
228 DO ng=1,ngrids
229 IF ((iic(ng).ne.ntstart(ng)).and. &
230 & mod(iic(ng)-1,couplesteps(iwaves,ng)).eq.0) THEN
231 DO tile=first_tile(ng),last_tile(ng),+1
232 CALL ocn2wav_coupling (ng, tile)
233 END DO
234!$OMP BARRIER
235 END IF
236 END DO
237# endif
238
239# ifdef NEARSHORE_MELLOR_NOT_YET
240!
241!-----------------------------------------------------------------------
242! Compute radiation stress terms.
243!-----------------------------------------------------------------------
244!
245 DO ng=1,ngrids
246 DO tile=last_tile(ng),first_tile(ng),-1
247 CALL rp_radiation_stress (ng, tile)
248 END DO
249!$OMP BARRIER
250 END DO
251# endif
252!
253!-----------------------------------------------------------------------
254! Set vertical boundary conditions. Process tidal forcing.
255!-----------------------------------------------------------------------
256!
257 DO ng=1,ngrids
258 DO tile=first_tile(ng),last_tile(ng),+1
259 CALL rp_set_vbc (ng, tile)
260# if defined SSH_TIDES_NOT_YET || defined UV_TIDES_NOT_YET
261 CALL rp_set_tides (ng, tile)
262# endif
263 END DO
264!$OMP BARRIER
265 END DO
266
267# ifdef ADJUST_BOUNDARY
268!
269!-----------------------------------------------------------------------
270! Interpolate open boundary increments and adjust open boundaries.
271! Skip the last output timestep.
272!-----------------------------------------------------------------------
273!
274 DO ng=1,ngrids
275 IF (iic(ng).lt.(ntend(ng)+1)) THEN
276 DO tile=first_tile(ng),last_tile(ng),+1
277 CALL rp_obc_adjust (ng, tile, lbinp(ng))
278 END DO
279!$OMP BARRIER
280 END IF
281 END DO
282# endif
283
284# ifdef ADJUST_WSTRESS
285!
286!-----------------------------------------------------------------------
287! Interpolate surface forcing increments and adjust surface forcing.
288! Skip the last output timestep.
289!-----------------------------------------------------------------------
290!
291 DO ng=1,ngrids
292 IF (iic(ng).lt.(ntend(ng)+1)) THEN
293 DO tile=first_tile(ng),last_tile(ng),+1
294 CALL rp_frc_adjust (ng, tile, lfinp(ng))
295 END DO
296!$OMP BARRIER
297 END IF
298 END DO
299# endif
300!
301!-----------------------------------------------------------------------
302! If appropriate, write out fields into output NetCDF files. Notice
303! that IO data is written in delayed and serial mode. Exit if last
304! time step.
305!-----------------------------------------------------------------------
306!
307 DO ng=1,ngrids
308!$OMP MASTER
309 CALL rp_output (ng)
310!$OMP END MASTER
311!$OMP BARRIER
312 IF ((founderror(exit_flag, noerror, __line__, myfile)).or. &
313 & ((iic(ng).eq.(ntend(ng)+1)).and.(ng.eq.ngrids))) RETURN
314 END DO
315!
316!-----------------------------------------------------------------------
317! Solve the vertically integrated primitive equations for the
318! free-surface and momentum components.
319!-----------------------------------------------------------------------
320!
321! Set time indices for predictor step. The PREDICTOR_2D_STEP switch
322! it is assumed to be false before the first time-step.
323!
324 DO ng=1,ngrids
325 iif(ng)=1
326 nfast(ng)=1
327 next_indx1=3-indx1(ng)
328 IF (.not.predictor_2d_step(ng)) THEN
329 predictor_2d_step(ng)=.true.
330 IF (first_2d_step) THEN
331 kstp(ng)=indx1(ng)
332 ELSE
333 kstp(ng)=3-indx1(ng)
334 END IF
335 knew(ng)=3
336 krhs(ng)=indx1(ng)
337 END IF
338!
339! Predictor step - Advance barotropic equations using 2D time-step
340! ============== predictor scheme.
341!
342 DO tile=last_tile(ng),first_tile(ng),-1
343 CALL rp_step2d (ng, tile)
344 END DO
345!$OMP BARRIER
346 END DO
347!
348! Set time indices for corrector step.
349!
350 DO ng=1,ngrids
351 IF (predictor_2d_step(ng)) THEN
352 predictor_2d_step(ng)=.false.
353 knew(ng)=next_indx1
354 kstp(ng)=3-knew(ng)
355 krhs(ng)=3
356 IF (iif(ng).lt.(nfast(ng)+1)) indx1(ng)=next_indx1
357 END IF
358!
359! Corrector step - Apply 2D time-step corrector scheme. Notice that
360! ============== there is not need for a corrector step during the
361! auxiliary (nfast+1) time-step.
362!
363 IF (iif(ng).lt.(nfast(ng)+1)) THEN
364 DO tile=first_tile(ng),last_tile(ng),+1
365 CALL rp_step2d (ng, tile)
366 END DO
367!$OMP BARRIER
368 END IF
369 END DO
370
371# ifdef FLOATS_NOT_YET
372!
373!-----------------------------------------------------------------------
374! Compute Lagrangian drifters trajectories: Split all the drifters
375! between all the computational threads, except in distributed-memory
376! and serial configurations. In distributed-memory, the parallel node
377! containing the drifter is selected internally since the state
378! variables do not have a global scope.
379!-----------------------------------------------------------------------
380!
381 DO ng=1,ngrids
382 IF (lfloats(ng)) THEN
383# ifdef _OPENMP
384 chunk_size=(nfloats(ng)+numthreads-1)/numthreads
385 lstr=1+mythread*chunk_size
386 lend=min(nfloats(ng),lstr+chunk_size-1)
387# else
388 lstr=1
389 lend=nfloats(ng)
390# endif
391 CALL rp_step_floats (ng, lstr, lend)
392!$OMP BARRIER
393!
394! Shift floats time indices.
395!
396 nfp1(ng)=mod(nfp1(ng)+1,nft+1)
397 nf(ng)=mod(nf(ng)+1,nft+1)
398 nfm1(ng)=mod(nfm1(ng)+1,nft+1)
399 nfm2(ng)=mod(nfm2(ng)+1,nft+1)
400 nfm3(ng)=mod(nfm3(ng)+1,nft+1)
401 END IF
402 END DO
403# endif
404 END DO step_loop
405
406 RETURN
407 END SUBROUTINE rp_main2d
408#else
409 SUBROUTINE rp_main2d
410 RETURN
411 END SUBROUTINE rp_main2d
412#endif
subroutine, public time_string(mytime, date_string)
Definition dateclock.F:1272
integer, dimension(:,:), allocatable couplesteps
integer iatmos
integer iwaves
integer numthreads
integer, dimension(:), allocatable first_tile
integer mythread
integer, dimension(:), allocatable last_tile
integer, dimension(:), allocatable nfloats
Definition mod_param.F:543
integer, parameter irpm
Definition mod_param.F:664
integer, parameter nft
Definition mod_param.F:539
integer ngrids
Definition mod_param.F:113
logical, dimension(:), allocatable lfloats
integer, dimension(:), allocatable iic
real(dp), dimension(:), allocatable dt
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
real(dp), dimension(:), allocatable time
integer, dimension(:), allocatable ntstart
integer, dimension(:), allocatable iif
integer noerror
integer, dimension(:), allocatable nfm2
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 krhs
integer, dimension(:), allocatable nstp
subroutine, public rp_diag(ng, tile)
Definition rp_diag.F:26
subroutine, public rp_frc_adjust(ng, tile, linp)
subroutine, public rp_ini_fields(ng, tile, model)
subroutine, public rp_ini_zeta(ng, tile, model)
subroutine, public rp_obc_adjust(ng, tile, linp)
subroutine, public rp_set_vbc(ng, tile)
Definition rp_set_vbc.F:33
subroutine, public rp_step2d(ng, tile)
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52
subroutine, public tl_forcing(ng, tile, kfrc, nfrc)
Definition tl_forcing.F:26
subroutine rp_get_data(ng)
Definition rp_get_data.F:4
subroutine rp_main2d
Definition rp_main2d.F:410
subroutine rp_output(ng)
Definition rp_output.F:4
subroutine rp_set_data(ng, tile)
Definition rp_set_data.F:4
subroutine wpoints(ng, tile, model)
Definition wpoints.F:40