ROMS
Loading...
Searching...
No Matches
set_tides.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#if defined NONLINEAR && (defined SSH_TIDES || defined UV_TIDES)
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group Robert Hetland !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This routine adds tidal elevation (m) and tidal currents (m/s) to !
13! sea surface height and 2D momentum climatologies, respectively. !
14! !
15!=======================================================================
16!
17 implicit none
18!
19 PRIVATE
20 PUBLIC :: set_tides
21!
22 CONTAINS
23!
24!***********************************************************************
25 SUBROUTINE set_tides (ng, tile)
26!***********************************************************************
27!
28 USE mod_param
29 USE mod_grid
30 USE mod_tides
31# ifdef NESTING
32# if defined AVERAGES && defined AVERAGES_DETIDE && \
33 (defined ssh_tides || defined uv_tides)
34 USE mod_scalars
35# endif
36# endif
37 USE mod_stepping
38!
39! Imported variable declarations.
40!
41 integer, intent(in) :: ng, tile
42!
43! Local variable declarations.
44!
45# ifdef NESTING
46# if defined AVERAGES && defined AVERAGES_DETIDE && \
47 (defined ssh_tides || defined uv_tides)
48 integer :: itide, mg
49!
50# endif
51# endif
52
53 character (len=*), parameter :: myfile = &
54 & __FILE__
55!
56# include "tile.h"
57!
58# ifdef PROFILE
59 CALL wclock_on (ng, inlm, 11, __line__, myfile)
60# endif
61 CALL set_tides_tile (ng, tile, &
62 & lbi, ubi, lbj, ubj, &
63 & imins, imaxs, jmins, jmaxs, &
64 & ntc(ng), &
65 & grid(ng) % angler, &
66# ifdef MASKING
67 & grid(ng) % rmask, &
68 & grid(ng) % umask, &
69 & grid(ng) % vmask, &
70# endif
71# ifdef SSH_TIDES
72 & tides(ng) % SSH_Tamp, &
73 & tides(ng) % SSH_Tphase, &
74# endif
75# ifdef UV_TIDES
76 & tides(ng) % UV_Tangle, &
77 & tides(ng) % UV_Tphase, &
78 & tides(ng) % UV_Tmajor, &
79 & tides(ng) % UV_Tminor, &
80# endif
81# if defined AVERAGES && defined AVERAGES_DETIDE && \
82 (defined ssh_tides || defined uv_tides)
83 & tides(ng) % SinOmega, &
84 & tides(ng) % CosOmega, &
85# endif
86 & tides(ng) % Tperiod)
87
88# ifdef NESTING
89# if defined AVERAGES && defined AVERAGES_DETIDE && \
90 (defined ssh_tides || defined uv_tides)
91!
92! If nested grids and detiding, load period and harmonics to other
93! grids.
94!
95 IF (lprocesstides(ng)) THEN
96 DO mg=1,ngrids
97 IF (ng.ne.mg) THEN
98 DO itide=1,ntc(ng)
99 tides(mg)%Tperiod (1:ntc(ng))=tides(ng)%Tperiod (1:ntc(ng))
100 tides(mg)%SinOmega(1:ntc(ng))=tides(ng)%SinOmega(1:ntc(ng))
101 tides(mg)%CosOmega(1:ntc(ng))=tides(ng)%CosOmega(1:ntc(ng))
102 END DO
103 END IF
104 END DO
105 END IF
106# endif
107# endif
108# ifdef PROFILE
109 CALL wclock_off (ng, inlm, 11, __line__, myfile)
110# endif
111!
112 RETURN
113 END SUBROUTINE set_tides
114!
115!***********************************************************************
116 SUBROUTINE set_tides_tile (ng, tile, &
117 & LBi, UBi, LBj, UBj, &
118 & IminS, ImaxS, JminS, JmaxS, &
119 & NTC, &
120 & angler, &
121# ifdef MASKING
122 & rmask, umask, vmask, &
123# endif
124# ifdef SSH_TIDES
125 & SSH_Tamp, SSH_Tphase, &
126# endif
127# ifdef UV_TIDES
128 & UV_Tangle, UV_Tphase, &
129 & UV_Tmajor, UV_Tminor, &
130# endif
131# if defined AVERAGES && defined AVERAGES_DETIDE && \
132 (defined SSH_TIDES || defined UV_TIDES)
133 & sinomega, cosomega, &
134# endif
135 & tperiod)
136!***********************************************************************
137!
138 USE mod_param
139 USE mod_boundary
140 USE mod_clima
141 USE mod_ncparam
142 USE mod_scalars
143!
144# ifdef DISTRIBUTE
145 USE distribute_mod, ONLY : mp_boundary
146# endif
148# ifdef DISTRIBUTE
150# endif
151!
152! Imported variables declarations.
153!
154 integer, intent(in) :: ng, tile
155 integer, intent(in) :: LBi, UBi, LBj, UBj
156 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
157 integer, intent(in) :: NTC
158!
159# ifdef ASSUMED_SHAPE
160 real(r8), intent(in) :: angler(LBi:,LBj:)
161# ifdef MASKING
162 real(r8), intent(in) :: rmask(LBi:,LBj:)
163 real(r8), intent(in) :: umask(LBi:,LBj:)
164 real(r8), intent(in) :: vmask(LBi:,LBj:)
165# endif
166 real(r8), intent(in) :: Tperiod(MTC)
167# ifdef SSH_TIDES
168 real(r8), intent(in) :: SSH_Tamp(LBi:,LBj:,:)
169 real(r8), intent(in) :: SSH_Tphase(LBi:,LBj:,:)
170# endif
171# ifdef UV_TIDES
172 real(r8), intent(in) :: UV_Tangle(LBi:,LBj:,:)
173 real(r8), intent(in) :: UV_Tmajor(LBi:,LBj:,:)
174 real(r8), intent(in) :: UV_Tminor(LBi:,LBj:,:)
175 real(r8), intent(in) :: UV_Tphase(LBi:,LBj:,:)
176# endif
177# if defined AVERAGES && defined AVERAGES_DETIDE && \
178 (defined ssh_tides || defined uv_tides)
179 real(r8), intent(inout) :: SinOmega(:)
180 real(r8), intent(inout) :: CosOmega(:)
181# endif
182# else
183 real(r8), intent(in) :: angler(LBi:UBi,LBj:UBj)
184# ifdef MASKING
185 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
186 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
187 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
188# endif
189 real(r8), intent(in) :: Tperiod(MTC)
190# ifdef SSH_TIDES
191 real(r8), intent(in) :: SSH_Tamp(LBi:UBi,LBj:UBj,MTC)
192 real(r8), intent(in) :: SSH_Tphase(LBi:UBi,LBj:UBj,MTC)
193# endif
194# ifdef UV_TIDES
195 real(r8), intent(in) :: UV_Tangle(LBi:UBi,LBj:UBj,MTC)
196 real(r8), intent(in) :: UV_Tmajor(LBi:UBi,LBj:UBj,MTC)
197 real(r8), intent(in) :: UV_Tminor(LBi:UBi,LBj:UBj,MTC)
198 real(r8), intent(in) :: UV_Tphase(LBi:UBi,LBj:UBj,MTC)
199# endif
200# if defined AVERAGES && defined AVERAGES_DETIDE && \
201 (defined ssh_tides || defined uv_tides)
202 real(r8), intent(inout) :: SinOmega(MTC)
203 real(r8), intent(inout) :: CosOmega(MTC)
204# endif
205# endif
206!
207! Local variables declarations.
208!
209 logical :: update
210
211# ifdef DISTRIBUTE
212 integer :: ILB, IUB, JLB, JUB
213# endif
214 integer :: i, itide, j
215
216 real(r8) :: Cangle, Cphase, Sangle, Sphase
217 real(r8) :: angle, cff, phase, omega, ramp
218 real(r8) :: bry_cor, bry_pgr, bry_str, bry_val
219
220 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Etide
221 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Utide
222 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Uwrk
223 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Vtide
224 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Vwrk
225
226# include "set_bounds.h"
227
228# ifdef DISTRIBUTE
229!
230! Lower and upper bounds for nontiled (global values) boundary arrays.
231!
232 ilb=bounds(ng)%LBi(-1)
233 iub=bounds(ng)%UBi(-1)
234 jlb=bounds(ng)%LBj(-1)
235 jub=bounds(ng)%UBj(-1)
236# endif
237!
238!=======================================================================
239! Process tidal forcing used at the lateral boundaries.
240!=======================================================================
241!
242! In refinement nesting applications, tidal forcing is only necessary
243! in the coarser large scale grid.
244!
245 needed : IF (lprocesstides(ng)) THEN
246!
247! Set time-ramping parameter.
248!
249# ifdef RAMP_TIDES
250 ramp=tanh((tdays(ng)-dstart)/1.0_r8)
251# else
252 ramp=1.0_r8
253# endif
254# if defined AVERAGES && defined AVERAGES_DETIDE && \
255 (defined ssh_tides || defined uv_tides)
256!
257!-----------------------------------------------------------------------
258! Compute harmonic used to detide output fields.
259!-----------------------------------------------------------------------
260!
261 cff=2.0_r8*pi*(time(ng)-tide_start*day2sec)
262 DO itide=1,ntc
263 IF (tperiod(itide).gt.0.0_r8) THEN
264 omega=cff/tperiod(itide)
265 sinomega(itide)=sin(omega)
266 cosomega(itide)=cos(omega)
267 ELSE
268 sinomega(itide)=0.0_r8
269 cosomega(itide)=0.0_r8
270 END IF
271 END DO
272# endif
273# ifdef SSH_TIDES
274!
275!-----------------------------------------------------------------------
276! Add tidal elevation (m) to sea surface height climatology.
277!-----------------------------------------------------------------------
278!
279 etide(:,:)=0.0_r8
280 cff=2.0_r8*pi*(time(ng)-tide_start*day2sec)
281 DO itide=1,ntc
282 IF (tperiod(itide).gt.0.0_r8) THEN
283 omega=cff/tperiod(itide)
284 DO j=jstrr,jendr
285 DO i=istrr,iendr
286 etide(i,j)=etide(i,j)+ &
287 & ramp*ssh_tamp(i,j,itide)* &
288 & cos(omega-ssh_tphase(i,j,itide))
289# ifdef MASKING
290 etide(i,j)=etide(i,j)*rmask(i,j)
291# endif
292 END DO
293 END DO
294 END IF
295 END DO
296
297# ifdef ADD_FSOBC
298!
299! Add sub-tidal forcing and adjust climatology to include tides.
300!
301 IF (lsshclm(ng)) THEN
302 DO j=jstrr,jendr
303 DO i=istrr,iendr
304 clima(ng)%ssh(i,j)=clima(ng)%ssh(i,j)+etide(i,j)
305 END DO
306 END DO
307 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
308 CALL exchange_r2d_tile (ng, tile, &
309 & lbi, ubi, lbj, ubj, &
310 & clima(ng)%ssh)
311 END IF
312# ifdef DISTRIBUTE
313 CALL mp_exchange2d (ng, tile, inlm, 1, &
314 & lbi, ubi, lbj, ubj, &
315 & nghostpoints, &
316 & ewperiodic(ng), nsperiodic(ng), &
317 & clima(ng)%ssh)
318# endif
319 END IF
320# endif
321!
322! If appropriate, load tidal forcing into boundary arrays. The "zeta"
323! boundary arrays are important for the Flather or reduced physics
324! boundary conditions for 2D momentum. To avoid having two boundary
325! points for these arrays, the values of "zeta_west" and "zeta_east"
326! are averaged at u-points. Similarly, the values of "zeta_south"
327! and "zeta_north" is averaged at v-points. Noticed that these
328! arrays are also used for the clamped conditions for free-surface.
329! This averaging is less important for that type ob boundary
330! conditions.
331!
332 IF (lbc(iwest,isfsur,ng)%acquire.or. &
333 & lbc(iwest,isubar,ng)%acquire.or. &
334 & lbc(iwest,isvbar,ng)%acquire) THEN
335 update=.false.
336 IF (domain(ng)%Western_Edge(tile)) THEN
337 DO j=jstrr,jendr
338# ifdef ADD_FSOBC
339 boundary(ng)%zeta_west(j)=boundary(ng)%zeta_west(j)+ &
340 & 0.5_r8*(etide(istr-1,j)+ &
341 & etide(istr ,j))
342# else
343 boundary(ng)%zeta_west(j)=0.5_r8*(etide(istr-1,j)+ &
344 & etide(istr ,j))
345# endif
346 END DO
347 update=.true.
348 END IF
349# ifdef DISTRIBUTE
350 CALL mp_boundary (ng, inlm, jstrr, jendr, &
351 & jlb, jub, 1, 1, update, &
352 & boundary(ng)%zeta_west)
353# endif
354 END IF
355!
356 IF (lbc(ieast,isfsur,ng)%acquire.or. &
357 & lbc(ieast,isubar,ng)%acquire.or. &
358 & lbc(ieast,isvbar,ng)%acquire) THEN
359 update=.false.
360 IF (domain(ng)%Eastern_Edge(tile)) THEN
361 DO j=jstrr,jendr
362# ifdef ADD_FSOBC
363 boundary(ng)%zeta_east(j)=boundary(ng)%zeta_east(j)+ &
364 & 0.5_r8*(etide(iend ,j)+ &
365 & etide(iend+1,j))
366# else
367 boundary(ng)%zeta_east(j)=0.5_r8*(etide(iend ,j)+ &
368 & etide(iend+1,j))
369# endif
370 END DO
371 update=.true.
372 END IF
373# ifdef DISTRIBUTE
374 CALL mp_boundary (ng, inlm, jstrr, jendr, &
375 & jlb, jub, 1, 1, update, &
376 & boundary(ng)%zeta_east)
377# endif
378 END IF
379!
380 IF (lbc(isouth,isfsur,ng)%acquire.or. &
381 & lbc(isouth,isubar,ng)%acquire.or. &
382 & lbc(isouth,isvbar,ng)%acquire) THEN
383 update=.false.
384 IF (domain(ng)%Southern_Edge(tile)) THEN
385 DO i=istrr,iendr
386# ifdef ADD_FSOBC
387 boundary(ng)%zeta_south(i)=boundary(ng)%zeta_south(i)+ &
388 & 0.5_r8*(etide(i,jstr-1)+ &
389 & etide(i,jstr ))
390# else
391 boundary(ng)%zeta_south(i)=0.5_r8*(etide(i,jstr-1)+ &
392 & etide(i,jstr ))
393# endif
394 END DO
395 update=.true.
396 END IF
397# ifdef DISTRIBUTE
398 CALL mp_boundary (ng, inlm, istrr, iendr, &
399 & ilb, iub, 1, 1, update, &
400 & boundary(ng)%zeta_south)
401# endif
402 END IF
403!
404 IF (lbc(inorth,isfsur,ng)%acquire.or. &
405 & lbc(inorth,isubar,ng)%acquire.or. &
406 & lbc(inorth,isvbar,ng)%acquire) THEN
407 update=.false.
408 IF (domain(ng)%Northern_Edge(tile)) THEN
409 DO i=istrr,iendr
410# ifdef ADD_FSOBC
411 boundary(ng)%zeta_north(i)=boundary(ng)%zeta_north(i)+ &
412 & 0.5_r8*(etide(i,jend )+ &
413 & etide(i,jend+1))
414# else
415 boundary(ng)%zeta_north(i)=0.5_r8*(etide(i,jend )+ &
416 & etide(i,jend+1))
417# endif
418 END DO
419 update=.true.
420 END IF
421# ifdef DISTRIBUTE
422 CALL mp_boundary (ng, inlm, istrr, iendr, &
423 & ilb, iub, 1, 1, update, &
424 & boundary(ng)%zeta_north)
425# endif
426 END IF
427# endif
428
429# if defined UV_TIDES
430!
431!-----------------------------------------------------------------------
432! Add tidal currents (m/s) to 2D momentum climatologies.
433!-----------------------------------------------------------------------
434!
435 utide(:,:)=0.0_r8
436 vtide(:,:)=0.0_r8
437 cff=2.0_r8*pi*(time(ng)-tide_start*day2sec)
438 DO itide=1,ntc
439 IF (tperiod(itide).gt.0.0_r8) THEN
440 omega=cff/tperiod(itide)
441 DO j=min(jstrr,jstr-1),jendr
442 DO i=min(istrr,istr-1),iendr
443 angle=uv_tangle(i,j,itide)-angler(i,j)
444 cangle=cos(angle)
445 sangle=sin(angle)
446 phase=omega-uv_tphase(i,j,itide)
447 cphase=cos(phase)
448 sphase=sin(phase)
449 uwrk(i,j)=uv_tmajor(i,j,itide)*cangle*cphase- &
450 & uv_tminor(i,j,itide)*sangle*sphase
451 vwrk(i,j)=uv_tmajor(i,j,itide)*sangle*cphase+ &
452 & uv_tminor(i,j,itide)*cangle*sphase
453 END DO
454 END DO
455 DO j=jstrr,jendr
456 DO i=istr,iendr
457 utide(i,j)=utide(i,j)+ &
458 & ramp*0.5_r8*(uwrk(i-1,j)+uwrk(i,j))
459# ifdef MASKING
460 utide(i,j)=utide(i,j)*umask(i,j)
461# endif
462 END DO
463 END DO
464 DO j=jstr,jendr
465 DO i=istrr,iendr
466 vtide(i,j)=(vtide(i,j)+ &
467 & ramp*0.5_r8*(vwrk(i,j-1)+vwrk(i,j)))
468# ifdef MASKING
469 vtide(i,j)=vtide(i,j)*vmask(i,j)
470# endif
471 END DO
472 END DO
473 END IF
474 END DO
475
476# ifdef ADD_M2OBC
477!
478! Add sub-tidal forcing and adjust climatology to include tides.
479!
480 IF (lm2clm(ng)) THEN
481 DO j=jstrr,jendr
482 DO i=istr,iendr
483 clima(ng)%ubarclm(i,j)=clima(ng)%ubarclm(i,j)+utide(i,j)
484 END DO
485 END DO
486 DO j=jstr,jendr
487 DO i=istrr,iendr
488 clima(ng)%vbarclm(i,j)=clima(ng)%vbarclm(i,j)+vtide(i,j)
489 END DO
490 END DO
491 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
492 CALL exchange_u2d_tile (ng, tile, &
493 & lbi, ubi, lbj, ubj, &
494 & clima(ng)%ubarclm)
495 CALL exchange_v2d_tile (ng, tile, &
496 & lbi, ubi, lbj, ubj, &
497 & clima(ng)%vbarclm)
498 END IF
499# ifdef DISTRIBUTE
500 CALL mp_exchange2d (ng, tile, inlm, 2, &
501 & lbi, ubi, lbj, ubj, &
502 & nghostpoints, &
503 & ewperiodic(ng), nsperiodic(ng), &
504 & clima(ng)%ubarclm, &
505 & clima(ng)%vbarclm)
506# endif
507 END IF
508# endif
509!
510! If appropriate, load tidal forcing into boundary arrays.
511!
512 IF (lbc(iwest,isubar,ng)%acquire.and. &
513 & lbc(iwest,isvbar,ng)%acquire) THEN
514 update=.false.
515 IF (domain(ng)%Western_Edge(tile)) THEN
516 DO j=jstrr,jendr
517# ifdef ADD_M2OBC
518 boundary(ng)%ubar_west(j)=boundary(ng)%ubar_west(j)+ &
519 & utide(istr,j)
520# else
521 boundary(ng)%ubar_west(j)=utide(istr,j)
522# endif
523 END DO
524 DO j=jstr,jendr
525# ifdef ADD_M2OBC
526 boundary(ng)%vbar_west(j)=boundary(ng)%vbar_west(j)+ &
527 & vtide(istr-1,j)
528# else
529 boundary(ng)%vbar_west(j)=vtide(istr-1,j)
530# endif
531 END DO
532 update=.true.
533 END IF
534# ifdef DISTRIBUTE
535 CALL mp_boundary (ng, inlm, jstrr, jendr, &
536 & jlb, jub, 1, 1, update, &
537 & boundary(ng)%ubar_west)
538 CALL mp_boundary (ng, inlm, jstr, jendr, &
539 & jlb, jub, 1, 1, update, &
540 & boundary(ng)%vbar_west)
541# endif
542 END IF
543!
544 IF (lbc(ieast,isubar,ng)%acquire.and. &
545 & lbc(ieast,isvbar,ng)%acquire) THEN
546 update=.false.
547 IF (domain(ng)%Eastern_Edge(tile)) THEN
548 DO j=jstrr,jendr
549# ifdef ADD_M2OBC
550 boundary(ng)%ubar_east(j)=boundary(ng)%ubar_east(j)+ &
551 & utide(iend+1,j)
552# else
553 boundary(ng)%ubar_east(j)=utide(iend+1,j)
554# endif
555 END DO
556 DO j=jstr,jendr
557# ifdef ADD_M2OBC
558 boundary(ng)%vbar_east(j)=boundary(ng)%vbar_east(j)+ &
559 & vtide(iend+1,j)
560# else
561 boundary(ng)%vbar_east(j)=vtide(iend+1,j)
562# endif
563 END DO
564 update=.true.
565 END IF
566# ifdef DISTRIBUTE
567 CALL mp_boundary (ng, inlm, jstrr, jendr, &
568 & jlb, jub, 1, 1, update, &
569 & boundary(ng)%ubar_east)
570 CALL mp_boundary (ng, inlm, jstr, jendr, &
571 & jlb, jub, 1, 1, update, &
572 & boundary(ng)%vbar_east)
573# endif
574 END IF
575!
576 IF (lbc(isouth,isubar,ng)%acquire.and. &
577 & lbc(isouth,isvbar,ng)%acquire) THEN
578 update=.false.
579 IF (domain(ng)%Southern_Edge(tile)) THEN
580 DO i=istr,iendr
581# ifdef ADD_M2OBC
582 boundary(ng)%ubar_south(i)=boundary(ng)%ubar_south(i)+ &
583 & utide(i,jstr-1)
584# else
585 boundary(ng)%ubar_south(i)=utide(i,jstr-1)
586# endif
587 END DO
588 DO i=istrr,iendr
589# ifdef ADD_M2OBC
590 boundary(ng)%vbar_south(i)=boundary(ng)%vbar_south(i)+ &
591 & vtide(i,jstr)
592# else
593 boundary(ng)%vbar_south(i)=vtide(i,jstr)
594# endif
595 END DO
596 update=.true.
597 END IF
598# ifdef DISTRIBUTE
599 CALL mp_boundary (ng, inlm, istr, iendr, &
600 & ilb, iub, 1, 1, update, &
601 & boundary(ng)%ubar_south)
602 CALL mp_boundary (ng, inlm, istrr, iendr, &
603 & ilb, iub, 1, 1, update, &
604 & boundary(ng)%vbar_south)
605# endif
606 END IF
607!
608 IF (lbc(inorth,isubar,ng)%acquire.and. &
609 & lbc(inorth,isvbar,ng)%acquire) THEN
610 update=.false.
611 IF (domain(ng)%Northern_Edge(tile)) THEN
612 DO i=istr,iendr
613# ifdef ADD_M2OBC
614 boundary(ng)%ubar_north(i)=boundary(ng)%ubar_north(i)+ &
615 & utide(i,jend+1)
616# else
617 boundary(ng)%ubar_north(i)=utide(i,jend+1)
618# endif
619 END DO
620 DO i=istrr,iendr
621# ifdef ADD_M2OBC
622 boundary(ng)%vbar_north(i)=boundary(ng)%vbar_north(i)+ &
623 & vtide(i,jend+1)
624# else
625 boundary(ng)%vbar_north(i)=vtide(i,jend+1)
626# endif
627 END DO
628 update=.true.
629 END IF
630# ifdef DISTRIBUTE
631 CALL mp_boundary (ng, inlm, istr, iendr, &
632 & ilb, iub, 1, 1, update, &
633 & boundary(ng)%ubar_north)
634 CALL mp_boundary (ng, inlm, istrr, iendr, &
635 & ilb, iub, 1, 1, update, &
636 & boundary(ng)%vbar_north)
637# endif
638 END IF
639# endif
640 END IF needed
641!
642 RETURN
643 END SUBROUTINE set_tides_tile
644#endif
645 END MODULE set_tides_mod
subroutine mp_boundary(ng, model, imin, imax, lbi, ubi, lbk, ubk, update, a)
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_u2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_v2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
type(t_boundary), dimension(:), allocatable boundary
type(t_clima), dimension(:), allocatable clima
Definition mod_clima.F:153
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer isvbar
integer isfsur
integer isubar
integer, parameter inlm
Definition mod_param.F:662
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
integer nghostpoints
Definition mod_param.F:710
type(t_lbc), dimension(:,:,:), allocatable lbc
Definition mod_param.F:375
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer ngrids
Definition mod_param.F:113
real(dp), parameter day2sec
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
logical, dimension(:), allocatable lsshclm
logical, dimension(:), allocatable lprocesstides
real(dp), dimension(:), allocatable tdays
real(dp) dstart
integer, parameter isouth
logical, dimension(:), allocatable lm2clm
integer, parameter ieast
real(dp), dimension(:), allocatable time
integer, parameter inorth
real(dp) tide_start
real(dp), parameter pi
integer, dimension(:), allocatable ntc
type(t_tides), dimension(:), allocatable tides
Definition mod_tides.F:133
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine set_tides_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, ntc, angler, rmask, umask, vmask, ssh_tamp, ssh_tphase, uv_tangle, uv_tphase, uv_tmajor, uv_tminor, sinomega, cosomega, tperiod)
Definition set_tides.F:136
subroutine, public set_tides(ng, tile)
Definition set_tides.F:26
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