ROMS
Loading...
Searching...
No Matches
rp_ini_fields.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#ifdef TL_IOMS
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 initializes other time levels for tangent 2D fields. !
13! It also couples 3D and 2D momentum equations: it initializes 2D !
14! momentum (ubar,vbar) to the vertical integral of initial 3D !
15! momentum (u,v). !
16! !
17!=======================================================================
18!
19 USE mod_param
20 USE mod_grid
21# ifdef SOLVE3D
22 USE mod_coupling
23# endif
24 USE mod_ncparam
25 USE mod_ocean
26 USE mod_scalars
27# if defined SOLVE3D && \
28 defined sediment_not_yet && defined sed_morph_not_yet
29 USE mod_sedbed
30 USE mod_sediment
31# endif
32!
34# ifdef SOLVE3D
36# endif
37# ifdef DISTRIBUTE
39# ifdef SOLVE3D
41# endif
42# endif
43# ifdef SOLVE3D
44 USE rp_t3dbc_mod, ONLY : rp_t3dbc_tile
45 USE rp_u3dbc_mod, ONLY : rp_u3dbc_tile
46 USE rp_v3dbc_mod, ONLY : rp_v3dbc_tile
47# endif
48 USE rp_u2dbc_mod, ONLY : rp_u2dbc_tile
49 USE rp_v2dbc_mod, ONLY : rp_v2dbc_tile
51!
52 implicit none
53!
54 PRIVATE
55 PUBLIC :: rp_ini_fields
56 PUBLIC :: rp_ini_zeta
57# ifdef SOLVE3D
58 PUBLIC :: rp_set_zeta_timeavg
59# endif
60!
61 CONTAINS
62!
63!***********************************************************************
64 SUBROUTINE rp_ini_fields (ng, tile, model)
65!***********************************************************************
66!
67 USE mod_stepping
68!
69! Imported variable declarations.
70!
71 integer, intent(in) :: ng, tile, model
72!
73! Local variable declarations.
74!
75 character (len=*), parameter :: myfile = &
76 & __FILE__
77!
78# include "tile.h"
79!
80# ifdef PROFILE
81 CALL wclock_on (ng, model, 2, __line__, myfile)
82# endif
83 CALL rp_ini_fields_tile (ng, tile, model, &
84 & lbi, ubi, lbj, ubj, &
85 & imins, imaxs, jmins, jmaxs, &
86 & kstp(ng), krhs(ng), knew(ng), &
87# ifdef SOLVE3D
88 & nstp(ng), nnew(ng), &
89# endif
90# ifdef MASKING
91 & grid(ng) % rmask, &
92 & grid(ng) % umask, &
93 & grid(ng) % vmask, &
94# endif
95# ifdef SOLVE3D
96 & grid(ng) % Hz, &
97 & grid(ng) % tl_Hz, &
98 & ocean(ng) % tl_t, &
99 & ocean(ng) % u, &
100 & ocean(ng) % tl_u, &
101 & ocean(ng) % v, &
102 & ocean(ng) % tl_v, &
103# endif
104 & ocean(ng) % ubar, &
105 & ocean(ng) % tl_ubar, &
106 & ocean(ng) % vbar, &
107 & ocean(ng) % tl_vbar, &
108 & ocean(ng) % zeta, &
109 & ocean(ng) % tl_zeta)
110# ifdef PROFILE
111 CALL wclock_off (ng, model, 2, __line__, myfile)
112# endif
113!
114 RETURN
115 END SUBROUTINE rp_ini_fields
116!
117!***********************************************************************
118 SUBROUTINE rp_ini_fields_tile (ng, tile, model, &
119 & LBi, UBi, LBj, UBj, &
120 & IminS, ImaxS, JminS, JmaxS, &
121 & kstp, krhs, knew, &
122# ifdef SOLVE3D
123 & nstp, nnew, &
124# endif
125# ifdef MASKING
126 & rmask, umask, vmask, &
127# endif
128# ifdef SOLVE3D
129 & Hz, tl_Hz, &
130 & tl_t, &
131 & u, tl_u, &
132 & v, tl_v, &
133# endif
134 & ubar, tl_ubar, &
135 & vbar, tl_vbar, &
136 & zeta, tl_zeta)
137!***********************************************************************
138!
139! Imported variable declarations.
140!
141 integer, intent(in) :: ng, tile, model
142 integer, intent(in) :: LBi, UBi, LBj, UBj
143 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
144 integer, intent(in) :: kstp, krhs, knew
145# ifdef SOLVE3D
146 integer, intent(in) :: nstp, nnew
147# endif
148!
149# ifdef ASSUMED_SHAPE
150# ifdef MASKING
151 real(r8), intent(in) :: rmask(LBi:,LBj:)
152 real(r8), intent(in) :: umask(LBi:,LBj:)
153 real(r8), intent(in) :: vmask(LBi:,LBj:)
154# endif
155# ifdef SOLVE3D
156 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
157 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
158 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
159 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
160# endif
161 real(r8), intent(in) :: ubar(LBi:,LBj:,:)
162 real(r8), intent(in) :: vbar(LBi:,LBj:,:)
163 real(r8), intent(in) :: zeta(LBi:,LBj:,:)
164 real(r8), intent(in) :: tl_zeta(LBi:,LBj:,:)
165# ifdef SOLVE3D
166 real(r8), intent(inout) :: tl_t(LBi:,LBj:,:,:,:)
167 real(r8), intent(inout) :: tl_u(LBi:,LBj:,:,:)
168 real(r8), intent(inout) :: tl_v(LBi:,LBj:,:,:)
169# endif
170 real(r8), intent(inout) :: tl_ubar(LBi:,LBj:,:)
171 real(r8), intent(inout) :: tl_vbar(LBi:,LBj:,:)
172
173# else
174
175# ifdef MASKING
176 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
177 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
178 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
179# endif
180# ifdef SOLVE3D
181 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
182 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
183 real(r8), intent(in) :: u(LBi:UBi,LBj:UBj,N(ng),2)
184 real(r8), intent(in) :: v(LBi:UBi,LBj:UBj,N(ng),2)
185# endif
186 real(r8), intent(in) :: ubar(LBi:UBi,LBj:UBj,:)
187 real(r8), intent(in) :: vbar(LBi:UBi,LBj:UBj,:)
188 real(r8), intent(in) :: zeta(LBi:UBi,LBj:UBj,:)
189 real(r8), intent(in) :: tl_zeta(LBi:UBi,LBj:UBj,:)
190# ifdef SOLVE3D
191 real(r8), intent(inout) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
192 real(r8), intent(inout) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
193 real(r8), intent(inout) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
194# endif
195 real(r8), intent(inout) :: tl_ubar(LBi:UBi,LBj:UBj,:)
196 real(r8), intent(inout) :: tl_vbar(LBi:UBi,LBj:UBj,:)
197# endif
198!
199! Local variable declarations.
200!
201 integer :: i, ic, itrc, j, k
202
203 real(r8) :: cff1
204 real(r8) :: tl_cff1, tl_cff2
205
206# ifdef SOLVE3D
207 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: CF
208 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: DC
209
210 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_CF
211 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_DC
212# endif
213
214# include "set_bounds.h"
215
216# ifdef SOLVE3D
217!
218!-----------------------------------------------------------------------
219! Initialize other time levels for 3D momentum.
220!-----------------------------------------------------------------------
221!
222 DO j=jstrb,jendb
223 DO k=1,n(ng)
224 DO i=istrm,iendb
225!^ cff1=u(i,j,k,nstp)
226!^
227 tl_cff1=tl_u(i,j,k,nstp)
228# ifdef MASKING
229!^ cff1=cff1*umask(i,j)
230!^
231 tl_cff1=tl_cff1*umask(i,j)
232# endif
233!^ u(i,j,k,nstp)=cff1
234!^
235 tl_u(i,j,k,nstp)=tl_cff1
236 END DO
237 END DO
238!
239 IF (j.ge.jstrm) THEN
240 DO k=1,n(ng)
241 DO i=istrb,iendb
242!^ cff2=v(i,j,k,nstp)
243!^
244 tl_cff2=tl_v(i,j,k,nstp)
245# ifdef MASKING
246!^ cff2=cff2*vmask(i,j)
247!^
248 tl_cff2=tl_cff2*vmask(i,j)
249# endif
250!^ v(i,j,k,nstp)=cff2
251!^
252 tl_v(i,j,k,nstp)=tl_cff2
253 END DO
254 END DO
255 END IF
256 END DO
257!
258! Apply boundary conditions.
259!
260!^ CALL u3dbc_tile (ng, tile, &
261!^ & LBi, UBi, LBj, UBj, N(ng), &
262!^ & IminS, ImaxS, JminS, JmaxS, &
263!^ & nstp, nstp, &
264!^ & u)
265!^
266 CALL rp_u3dbc_tile (ng, tile, &
267 & lbi, ubi, lbj, ubj, n(ng), &
268 & imins, imaxs, jmins, jmaxs, &
269 & nstp, nstp, &
270 & tl_u)
271!^ CALL v3dbc_tile (ng, tile, &
272!^ & LBi, UBi, LBj, UBj, N(ng), &
273!^ & IminS, ImaxS, JminS, JmaxS, &
274!^ & nstp, nstp, &
275!^ & v)
276!^
277 CALL rp_v3dbc_tile (ng, tile, &
278 & lbi, ubi, lbj, ubj, n(ng), &
279 & imins, imaxs, jmins, jmaxs, &
280 & nstp, nstp, &
281 & tl_v)
282!
283 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
284!^ CALL exchange_u3d_tile (ng, tile, &
285!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
286!^ & u(:,:,:,nstp))
287!^
288 CALL exchange_u3d_tile (ng, tile, &
289 & lbi, ubi, lbj, ubj, 1, n(ng), &
290 & tl_u(:,:,:,nstp))
291!^ CALL exchange_v3d_tile (ng, tile, &
292!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
293!^ & v(:,:,:,nstp))
294!^
295 CALL exchange_v3d_tile (ng, tile, &
296 & lbi, ubi, lbj, ubj, 1, n(ng), &
297 & tl_v(:,:,:,nstp))
298 END IF
299
300# ifdef DISTRIBUTE
301!
302!^ CALL mp_exchange3d (ng, tile, model, 2, &
303!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
304!^ & NghostPoints, &
305!^ & EWperiodic(ng), NSperiodic(ng), &
306!^ & u(:,:,:,nstp), v(:,:,:,nstp))
307!^
308 CALL mp_exchange3d (ng, tile, model, 2, &
309 & lbi, ubi, lbj, ubj, 1, n(ng), &
310 & nghostpoints, &
311 & ewperiodic(ng), nsperiodic(ng), &
312 & tl_u(:,:,:,nstp), tl_v(:,:,:,nstp))
313# endif
314# endif
315
316# ifdef SOLVE3D
317!
318!-----------------------------------------------------------------------
319! Compute vertically-integrated momentum (tl_ubar, tl_vbar) from
320! initial 3D momentum (tl_u, tl_v).
321!-----------------------------------------------------------------------
322!
323! Here DC(i,1:N) are the grid cell thicknesses, DC(i,0) is the total
324! depth of the water column, and CF(i,0) is the vertical integral.
325!
326 DO j=jstrb,jendb
327 DO i=istrm,iendb
328 dc(i,0)=0.0_r8
329 tl_dc(i,0)=0.0_r8
330 cf(i,0)=0.0_r8
331 tl_cf(i,0)=0.0_r8
332 END DO
333 DO k=1,n(ng)
334 DO i=istrm,iendb
335 dc(i,k)=0.5_r8*(hz(i,j,k)+hz(i-1,j,k))
336 tl_dc(i,k)=0.5_r8*(tl_hz(i,j,k)+tl_hz(i-1,j,k))
337 dc(i,0)=dc(i,0)+dc(i,k)
338 tl_dc(i,0)=tl_dc(i,0)+tl_dc(i,k)
339 cf(i,0)=cf(i,0)+dc(i,k)*u(i,j,k,nstp)
340 tl_cf(i,0)=tl_cf(i,0)+tl_dc(i,k)*u(i,j,k,nstp)+ &
341 & dc(i,k)*tl_u(i,j,k,nstp)- &
342# ifdef TL_IOMS
343 & dc(i,k)*u(i,j,k,nstp)
344# endif
345 END DO
346 END DO
347 DO i=istrm,iendb
348 cff1=1.0_r8/dc(i,0)
349 tl_cff1=-cff1*cff1*tl_dc(i,0)+ &
350# ifdef TL_IOMS
351 & 2.0_r8*cff1
352# endif
353!^ cff2=CF(i,0)*cff1
354!^
355 tl_cff2=tl_cf(i,0)*cff1+cf(i,0)*tl_cff1- &
356# ifdef TL_IOMS
357 & cf(i,0)*cff1
358# endif
359# ifdef MASKING
360!^ cff2=cff2*umask(i,j)
361!^
362 tl_cff2=tl_cff2*umask(i,j)
363# endif
364!^ ubar(i,j,kstp)=cff2
365!^
366 tl_ubar(i,j,kstp)=tl_cff2
367 END DO
368!
369 IF (j.ge.jstrm) THEN
370 DO i=istrb,iendb
371 dc(i,0)=0.0_r8
372 tl_dc(i,0)=0.0_r8
373 cf(i,0)=0.0_r8
374 tl_cf(i,0)=0.0_r8
375 END DO
376 DO k=1,n(ng)
377 DO i=istrb,iendb
378 dc(i,k)=0.5_r8*(hz(i,j,k)+hz(i,j-1,k))
379 tl_dc(i,k)=0.5_r8*(tl_hz(i,j,k)+tl_hz(i,j-1,k))
380 dc(i,0)=dc(i,0)+dc(i,k)
381 tl_dc(i,0)=tl_dc(i,0)+tl_dc(i,k)
382 cf(i,0)=cf(i,0)+dc(i,k)*v(i,j,k,nstp)
383 tl_cf(i,0)=tl_cf(i,0)+tl_dc(i,k)*v(i,j,k,nstp)+ &
384 & dc(i,k)*tl_v(i,j,k,nstp)- &
385# ifdef TL_IOMS
386 & dc(i,k)*v(i,j,k,nstp)
387# endif
388 END DO
389 END DO
390 DO i=istrb,iendb
391 cff1=1.0_r8/dc(i,0)
392 tl_cff1=-cff1*cff1*tl_dc(i,0)+ &
393# ifdef TL_IOMS
394 & 2.0_r8*cff1
395# endif
396!^ cff2=CF(i,0)*cff1
397!^
398 tl_cff2=tl_cf(i,0)*cff1+cf(i,0)*tl_cff1- &
399# ifdef TL_IOMS
400 & cf(i,0)*cff1
401# endif
402# ifdef MASKING
403!^ cff2=cff2*vmask(i,j)
404!^
405 tl_cff2=tl_cff2*vmask(i,j)
406# endif
407!^ vbar(i,j,kstp)=cff2
408!^
409 tl_vbar(i,j,kstp)=tl_cff2
410 END DO
411 END IF
412 END DO
413!
414! Apply boundary conditions.
415!
416 IF (.not.(any(tl_lbc(:,isubar,ng)%radiation).or. &
417 & any(tl_lbc(:,isvbar,ng)%radiation).or. &
418 & any(tl_lbc(:,isubar,ng)%Flather).or. &
419 & any(tl_lbc(:,isvbar,ng)%Flather))) THEN
420!^ CALL u2dbc_tile (ng, tile, &
421!^ & LBi, UBi, LBj, UBj, &
422!^ & IminS, ImaxS, JminS, JmaxS, &
423!^ & krhs, kstp, kstp, &
424!^ & ubar, vbar, zeta)
425!^
426 CALL rp_u2dbc_tile (ng, tile, &
427 & lbi, ubi, lbj, ubj, &
428 & imins, imaxs, jmins, jmaxs, &
429 & krhs, kstp, kstp, &
430 & ubar, vbar, zeta, &
431 & tl_ubar, tl_vbar, tl_zeta)
432!^ CALL v2dbc_tile (ng, tile, &
433!^ & LBi, UBi, LBj, UBj, &
434!^ & IminS, ImaxS, JminS, JmaxS, &
435!^ & krhs, kstp, kstp, &
436!^ & ubar, vbar, zeta)
437!^
438 CALL rp_v2dbc_tile (ng, tile, &
439 & lbi, ubi, lbj, ubj, &
440 & imins, imaxs, jmins, jmaxs, &
441 & krhs, kstp, kstp, &
442 & ubar, vbar, zeta, &
443 & tl_ubar, tl_vbar, tl_zeta)
444 END IF
445!
446 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
447!^ CALL exchange_u2d_tile (ng, tile, &
448!^ & LBi, UBi, LBj, UBj, &
449!^ & ubar(:,:,kstp))
450!^
451 CALL exchange_u2d_tile (ng, tile, &
452 & lbi, ubi, lbj, ubj, &
453 & tl_ubar(:,:,kstp))
454!^ CALL exchange_v2d_tile (ng, tile, &
455!^ & LBi, UBi, LBj, UBj, &
456!^ & vbar(:,:,kstp))
457!^
458 CALL exchange_v2d_tile (ng, tile, &
459 & lbi, ubi, lbj, ubj, &
460 & tl_vbar(:,:,kstp))
461 END IF
462
463# ifdef DISTRIBUTE
464!
465!^ CALL mp_exchange2d (ng, tile, model, 2, &
466!^ & LBi, UBi, LBj, UBj, &
467!^ & NghostPoints, &
468!^ & EWperiodic(ng), NSperiodic(ng), &
469!^ & ubar(:,:,kstp), vbar(:,:,kstp))
470!^
471 CALL mp_exchange2d (ng, tile, model, 2, &
472 & lbi, ubi, lbj, ubj, &
473 & nghostpoints, &
474 & ewperiodic(ng), nsperiodic(ng), &
475 & tl_ubar(:,:,kstp), tl_vbar(:,:,kstp))
476# endif
477
478# else
479!
480!-----------------------------------------------------------------------
481! Initialize other time levels for 2D momentum (shallow-water model).
482!-----------------------------------------------------------------------
483!
484 DO j=jstrb,jendb
485 DO i=istrm,iendb
486!^ cff1=ubar(i,j,kstp)
487!^
488 tl_cff1=tl_ubar(i,j,kstp)
489# ifdef MASKING
490!^ cff1=cff1*umask(i,j)
491!^
492 tl_cff1=tl_cff1*umask(i,j)
493# endif
494!^ ubar(i,j,kstp)=cff1
495!^
496 tl_ubar(i,j,kstp)=tl_cff1
497 END DO
498!
499 IF (j.ge.jstrm) THEN
500 DO i=istrb,iendb
501!^ cff2=vbar(i,j,kstp)
502!^
503 tl_cff2=tl_vbar(i,j,kstp)
504# ifdef MASKING
505!^ cff2=cff2*vmask(i,j)
506!^
507 tl_cff2=tl_cff2*vmask(i,j)
508# endif
509!^ vbar(i,j,kstp)=cff2
510!^
511 tl_vbar(i,j,kstp)=tl_cff2
512 END DO
513 END IF
514 END DO
515!
516! Apply boundary conditions.
517!
518 IF (.not.(any(tl_lbc(:,isubar,ng)%radiation).or. &
519 & any(tl_lbc(:,isvbar,ng)%radiation).or. &
520 & any(tl_lbc(:,isubar,ng)%Flather).or. &
521 & any(tl_lbc(:,isvbar,ng)%Flather))) THEN
522!^ CALL u2dbc_tile (ng, tile, &
523!^ & LBi, UBi, LBj, UBj, &
524!^ & IminS, ImaxS, JminS, JmaxS, &
525!^ & krhs, kstp, kstp, &
526!^ & ubar, vbar, zeta)
527!^
528 CALL rp_u2dbc_tile (ng, tile, &
529 & lbi, ubi, lbj, ubj, &
530 & imins, imaxs, jmins, jmaxs, &
531 & krhs, kstp, kstp, &
532 & ubar, vbar, zeta, &
533 & tl_ubar, tl_vbar, tl_zeta)
534!^ CALL v2dbc_tile (ng, tile, &
535!^ & LBi, UBi, LBj, UBj, &
536!^ & IminS, ImaxS, JminS, JmaxS, &
537!^ & krhs, kstp, kstp, &
538!^ & ubar, vbar, zeta)
539!^
540 CALL rp_v2dbc_tile (ng, tile, &
541 & lbi, ubi, lbj, ubj, &
542 & imins, imaxs, jmins, jmaxs, &
543 & krhs, kstp, kstp, &
544 & ubar, vbar, zeta, &
545 & tl_ubar, tl_vbar, tl_zeta)
546 END IF
547!
548 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
549!^ CALL exchange_u2d_tile (ng, tile, &
550!^ & LBi, UBi, LBj, UBj, &
551!^ & ubar(:,:,kstp))
552!^
553 CALL exchange_u2d_tile (ng, tile, &
554 & lbi, ubi, lbj, ubj, &
555 & tl_ubar(:,:,kstp))
556!^ CALL exchange_v2d_tile (ng, tile, &
557!^ & LBi, UBi, LBj, UBj, &
558!^ & vbar(:,:,kstp))
559!^
560 CALL exchange_v2d_tile (ng, tile, &
561 & lbi, ubi, lbj, ubj, &
562 & tl_vbar(:,:,kstp))
563 END IF
564
565# ifdef DISTRIBUTE
566!
567!^ CALL mp_exchange2d (ng, tile, model, 2, &
568!^ & LBi, UBi, LBj, UBj, &
569!^ & NghostPoints, &
570!^ & EWperiodic(ng), NSperiodic(ng), &
571!^ & ubar(:,:,kstp), vbar(:,:,kstp))
572!^
573 CALL mp_exchange2d (ng, tile, model, 2, &
574 & lbi, ubi, lbj, ubj, &
575 & nghostpoints, &
576 & ewperiodic(ng), nsperiodic(ng), &
577 & tl_ubar(:,:,kstp), tl_vbar(:,:,kstp))
578# endif
579# endif
580
581# ifdef SOLVE3D
582!
583!-----------------------------------------------------------------------
584! Initialize other time levels for tracers.
585!-----------------------------------------------------------------------
586!
587 ic=0
588 DO itrc=1,nt(ng)
589 IF (ltracerclm(itrc,ng).and.lnudgetclm(itrc,ng)) THEN
590 ic=ic+1 ! OBC nudging coefficient index
591 END IF
592 DO k=1,n(ng)
593 DO j=jstrb,jendb
594 DO i=istrb,iendb
595!^ cff1=t(i,j,k,nstp,itrc)
596!^
597 tl_cff1=tl_t(i,j,k,nstp,itrc)
598# ifdef MASKING
599 tl_cff1=tl_cff1*rmask(i,j)
600# endif
601!^ t(i,j,k,nstp,itrc)=cff1
602!^
603 tl_t(i,j,k,nstp,itrc)=tl_cff1
604 END DO
605 END DO
606 END DO
607!
608! Apply boundary conditions.
609!
610!^ CALL t3dbc_tile (ng, tile, itrc, ic, &
611!^ & LBi, UBi, LBj, UBj, N(ng), NT(ng), &
612!^ & IminS, ImaxS, JminS, JmaxS, &
613!^ & nstp, nstp, &
614!^ & t)
615!^
616 CALL rp_t3dbc_tile (ng, tile, itrc, ic, &
617 & lbi, ubi, lbj, ubj, n(ng), nt(ng), &
618 & imins, imaxs, jmins, jmaxs, &
619 & nstp, nstp, &
620 & tl_t)
621!
622 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
623!^ CALL exchange_r3d_tile (ng, tile, &
624!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
625!^ & t(:,:,:,nstp,itrc))
626!^
627 CALL exchange_r3d_tile (ng, tile, &
628 & lbi, ubi, lbj, ubj, 1, n(ng), &
629 & tl_t(:,:,:,nstp,itrc))
630 END IF
631 END DO
632
633# ifdef DISTRIBUTE
634!
635!^ CALL mp_exchange4d (ng, tile, model, 1, &
636!^ & LBi, UBi, LBj, UBj, 1, N(ng), 1, NT(ng), &
637!^ & NghostPoints, &
638!^ & EWperiodic(ng), NSperiodic(ng), &
639!^ & t(:,:,:,nstp,:))
640!^
641 CALL mp_exchange4d (ng, tile, model, 1, &
642 & lbi, ubi, lbj, ubj, 1, n(ng), 1, nt(ng), &
643 & nghostpoints, &
644 & ewperiodic(ng), nsperiodic(ng), &
645 & tl_t(:,:,:,nstp,:))
646# endif
647# endif
648!
649 RETURN
650 END SUBROUTINE rp_ini_fields_tile
651
652!
653!***********************************************************************
654 SUBROUTINE rp_ini_zeta (ng, tile, model)
655!***********************************************************************
656!
657 USE mod_stepping
658!
659! Imported variable declarations.
660!
661 integer, intent(in) :: ng, tile, model
662!
663! Local variable declarations.
664!
665 character (len=*), parameter :: myfile = &
666 & __FILE__//", rp_ini_zeta"
667!
668# include "tile.h"
669!
670# ifdef PROFILE
671 CALL wclock_on (ng, model, 2, __line__, myfile)
672# endif
673 CALL rp_ini_zeta_tile (ng, tile, model, &
674 & lbi, ubi, lbj, ubj, &
675 & imins, imaxs, jmins, jmaxs, &
676 & kstp(ng), krhs(ng), knew(ng), &
677# ifdef MASKING
678 & grid(ng) % rmask, &
679# endif
680# ifdef WET_DRY_NOT_YET
681 & grid(ng) % h, &
682# endif
683# ifdef SOLVE3D
684# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
685 & sedbed(ng) % tl_bed, &
686 & sedbed(ng) % tl_bed_thick0, &
687 & sedbed(ng) % tl_bed_thick, &
688# endif
689 & coupling(ng) % tl_Zt_avg1, &
690# endif
691 & ocean(ng) % zeta, &
692 & ocean(ng) % tl_zeta)
693# ifdef PROFILE
694 CALL wclock_off (ng, model, 2, __line__, myfile)
695# endif
696!
697 RETURN
698 END SUBROUTINE rp_ini_zeta
699!
700!***********************************************************************
701 SUBROUTINE rp_ini_zeta_tile (ng, tile, model, &
702 & LBi, UBi, LBj, UBj, &
703 & IminS, ImaxS, JminS, JmaxS, &
704 & kstp, krhs, knew, &
705# ifdef MASKING
706 & rmask, &
707# endif
708# ifdef WET_DRY_NOT_YET
709 & h, &
710# endif
711# ifdef SOLVE3D
712# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
713 & tl_bed, &
714 & tl_bed_thick0, tl_bed_thick, &
715# endif
716 & tl_Zt_avg1, &
717# endif
718 & zeta, tl_zeta)
719!***********************************************************************
720!
721! Imported variable declarations.
722!
723 integer, intent(in) :: ng, tile, model
724 integer, intent(in) :: LBi, UBi, LBj, UBj
725 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
726 integer, intent(in) :: kstp, krhs, knew
727!
728# ifdef ASSUMED_SHAPE
729# ifdef MASKING
730 real(r8), intent(in) :: rmask(LBi:,LBj:)
731# endif
732# ifdef WET_DRY_NOT_YET
733 real(r8), intent(in) :: h(LBi:,LBj:)
734# endif
735 real(r8), intent(in) :: zeta(LBi:,LBj:,:)
736# ifdef SOLVE3D
737# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
738 real(r8), intent(in) :: tl_bed(LBi:,LBj:,:,:)
739 real(r8), intent(inout) :: tl_bed_thick0(LBi:,LBj:)
740 real(r8), intent(inout) :: tl_bed_thick(LBi:,LBj:,:)
741# endif
742 real(r8), intent(inout) :: tl_Zt_avg1(LBi:,LBj:)
743# endif
744 real(r8), intent(inout) :: tl_zeta(LBi:,LBj:,:)
745
746# else
747
748# ifdef MASKING
749 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
750# endif
751# ifdef WET_DRY_NOT_YET
752 real(r8), intent(in) :: h(LBi:UBi,LBj:UBj)
753# endif
754 real(r8), intent(in) :: zeta(LBi:UBi,LBj:UBj,:)
755# ifdef SOLVE3D
756# if defined SOLVE3D && defined SEDIMENT && defined SED_MORPH
757 real(r8), intent(in) :: tl_bed(LBi:UBi,LBj:UBj,Nbed,MBEDP)
758 real(r8), intent(inout) :: tl_bed_thick0(LBi:UBi,LBj:UBj)
759 real(r8), intent(inout) :: tl_bed_thick(LBi:UBi,LBj:UBj,3)
760# endif
761 real(r8), intent(inout) :: tl_Zt_avg1(LBi:UBi,LBj:UBj)
762# endif
763 real(r8), intent(inout) :: tl_zeta(LBi:UBi,LBj:UBj,:)
764# endif
765!
766! Local variable declarations.
767!
768 integer :: Imin, Imax, Jmin, Jmax
769 integer :: i, j, kbed
770
771 real(r8) :: cff1
772 real(r8) :: tl_cff1
773
774# include "set_bounds.h"
775!
776!-----------------------------------------------------------------------
777! Initialize other time levels for free-surface.
778!-----------------------------------------------------------------------
779!
780 IF (.not.(any(tl_lbc(:,isfsur,ng)%radiation).or. &
781 & any(tl_lbc(:,isfsur,ng)%Chapman_explicit).or. &
782 & any(tl_lbc(:,isfsur,ng)%Chapman_implicit))) THEN
783 imin=istrb
784 imax=iendb
785 jmin=jstrb
786 jmax=jendb
787 ELSE
788 imin=istrt
789 imax=iendt
790 jmin=jstrt
791 jmax=jendt
792 END IF
793 DO j=jmin,jmax
794 DO i=imin,imax
795!^ cff1=zeta(i,j,kstp)
796!^
797 tl_cff1=tl_zeta(i,j,kstp)
798# ifdef MASKING
799!^ cff1=cff1*rmask(i,j)
800!^
801 tl_cff1=tl_cff1*rmask(i,j)
802# endif
803!^ zeta(i,j,kstp)=cff1
804!^
805 tl_zeta(i,j,kstp)=tl_cff1
806 END DO
807 END DO
808!
809! Apply boundary conditions.
810!
811 IF (.not.(any(tl_lbc(:,isfsur,ng)%radiation).or. &
812 & any(tl_lbc(:,isfsur,ng)%Chapman_explicit).or. &
813 & any(tl_lbc(:,isfsur,ng)%Chapman_implicit))) THEN
814!^ CALL zetabc_tile (ng, tile, &
815!^ & LBi, UBi, LBj, UBj, &
816!^ & IminS, ImaxS, JminS, JmaxS, &
817!^ & krhs, kstp, kstp, &
818!^ & zeta)
819!^
820 CALL rp_zetabc_tile (ng, tile, &
821 & lbi, ubi, lbj, ubj, &
822 & imins, imaxs, jmins, jmaxs, &
823 & krhs, kstp, kstp, &
824 & zeta, &
825 & tl_zeta)
826 END IF
827!
828 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
829!^ CALL exchange_r2d_tile (ng, tile, &
830!^ & LBi, UBi, LBj, UBj, &
831!^ & zeta(:,:,kstp))
832!^
833 CALL exchange_r2d_tile (ng, tile, &
834 & lbi, ubi, lbj, ubj, &
835 & tl_zeta(:,:,kstp))
836 END IF
837
838# ifdef DISTRIBUTE
839!
840!^ CALL mp_exchange2d (ng, tile, model, 1, &
841!^ & LBi, UBi, LBj, UBj, &
842!^ & NghostPoints, &
843!^ & EWperiodic(ng), NSperiodic(ng), &
844!^ & zeta(:,:,kstp))
845!^
846 CALL mp_exchange2d (ng, tile, model, 1, &
847 & lbi, ubi, lbj, ubj, &
848 & nghostpoints, &
849 & ewperiodic(ng), nsperiodic(ng), &
850 & tl_zeta(:,:,kstp))
851# endif
852
853# ifdef SOLVE3D
854!
855!-----------------------------------------------------------------------
856! Initialize fast-time averaged free-surface (Zt_avg1) with the inital
857! free-surface
858!-----------------------------------------------------------------------
859!
860 CALL rp_set_zeta_timeavg (ng, tile, model)
861
862# if defined SEDIMENT_NOT_YET && defined SED_MORPH_NOT_YET
863!
864!-----------------------------------------------------------------------
865! Compute initial total thickness for all sediment bed layers.
866!-----------------------------------------------------------------------
867!
868 DO j=jstrt,jendt
869 DO i=istrt,iendt
870!^ bed_thick0(i,j)=0.0_r8
871!^
872 tl_bed_thick0(i,j)=0.0_r8
873 DO kbed=1,nbed
874!^ bed_thick0(i,j)=bed_thick0(i,j)+bed(i,j,kbed,ithck)
875!^
876 tl_bed_thick0(i,j)=tl_bed_thick0(i,j)+tl_bed(i,j,kbed,ithck)
877 END DO
878!^ bed_thick(i,j,1)=bed_thick0(i,j)
879!^ bed_thick(i,j,2)=bed_thick0(i,j)
880!^ bed_thick(i,j,3)=bed_thick0(i,j)
881!^
882 tl_bed_thick(i,j,1)=tl_bed_thick0(i,j)
883 tl_bed_thick(i,j,2)=tl_bed_thick0(i,j)
884 tl_bed_thick(i,j,3)=tl_bed_thick0(i,j)
885 END DO
886 END DO
887!
888 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
889!^ CALL exchange_r2d_tile (ng, tile, &
890!^ & LBi, UBi, LBj, UBj, &
891!^ & bed_thick0)
892!^
893 CALL exchange_r2d_tile (ng, tile, &
894 & lbi, ubi, lbj, ubj, &
895 & tl_bed_thick0)
896!^ CALL exchange_r2d_tile (ng, tile, &
897!^ & LBi, UBi, LBj, UBj, &
898!^ & bed_thick(:,:,1))
899!^
900 CALL exchange_r2d_tile (ng, tile, &
901 & lbi, ubi, lbj, ubj, &
902 & tl_bed_thick(:,:,1))
903!^ CALL exchange_r2d_tile (ng, tile, &
904!^ & LBi, UBi, LBj, UBj, &
905!^ & bed_thick(:,:,2))
906!^
907 CALL exchange_r2d_tile (ng, tile, &
908 & lbi, ubi, lbj, ubj, &
909 & tl_bed_thick(:,:,2))
910!^ CALL exchange_r2d_tile (ng, tile, &
911!^ & LBi, UBi, LBj, UBj, &
912!^ & bed_thick(:,:,3))
913!^
914 CALL exchange_r2d_tile (ng, tile, &
915 & lbi, ubi, lbj, ubj, &
916 & tl_bed_thick(:,:,3))
917 END IF
918
919# ifdef DISTRIBUTE
920!^ CALL mp_exchange2d (ng, tile, model, 4, &
921!^ & LBi, UBi, LBj, UBj, 1, 1, &
922!^ & NghostPoints, &
923!^ & EWperiodic(ng), NSperiodic(ng), &
924!^ & bed_thick0, &
925!^ & bed_thick(:,:,1), &
926!^ & bed_thick(:,:,2), &
927!^ & bed_thick(:,:,3))
928!^
929 CALL mp_exchange2d (ng, tile, model, 4, &
930 & lbi, ubi, lbj, ubj, &
931 & nghostpoints, &
932 & ewperiodic(ng), nsperiodic(ng), &
933 & tl_bed_thick0, &
934 & tl_bed_thick(:,:,1), &
935 & tl_bed_thick(:,:,2), &
936 & tl_bed_thick(:,:,3))
937# endif
938# endif
939# endif
940!
941 RETURN
942 END SUBROUTINE rp_ini_zeta_tile
943
944# ifdef SOLVE3D
945!
946!***********************************************************************
947 SUBROUTINE rp_set_zeta_timeavg (ng, tile, model)
948!***********************************************************************
949!
950 USE mod_stepping
951!
952! Imported variable declarations.
953!
954 integer, intent(in) :: ng, tile, model
955!
956! Local variable declarations.
957!
958 character (len=*), parameter :: myfile = &
959 & __FILE__//", tl_set_zeta_timeavg"
960!
961# include "tile.h"
962!
963# ifdef PROFILE
964 CALL wclock_on (ng, model, 2, __line__, myfile)
965# endif
966 CALL rp_set_zeta_timeavg_tile (ng, tile, model, &
967 & lbi, ubi, lbj, ubj, &
968 & kstp(ng), &
969 & coupling(ng) % tl_Zt_avg1, &
970 & ocean(ng) % tl_zeta)
971# ifdef PROFILE
972 CALL wclock_off (ng, model, 2, __line__, myfile)
973# endif
974!
975 RETURN
976 END SUBROUTINE rp_set_zeta_timeavg
977!
978!***********************************************************************
979 SUBROUTINE rp_set_zeta_timeavg_tile (ng, tile, model, &
980 & LBi, UBi, LBj, UBj, &
981 & kstp, &
982 & tl_Zt_avg1, &
983 & tl_zeta)
984!***********************************************************************
985!
986 USE mod_param
987 USE mod_scalars
988!
990# ifdef DISTRIBUTE
992# endif
993!
994! Imported variable declarations.
995!
996 integer, intent(in) :: ng, tile, model
997 integer, intent(in) :: LBi, UBi, LBj, UBj
998 integer, intent(in) :: kstp
999!
1000# ifdef ASSUMED_SHAPE
1001 real(r8), intent(inout) :: tl_Zt_avg1(LBi:,LBj:)
1002 real(r8), intent(inout) :: tl_zeta(LBi:,LBj:,:)
1003# else
1004 real(r8), intent(inout) :: tl_Zt_avg1(LBi:UBi,LBj:UBj)
1005 real(r8), intent(inout) :: tl_zeta(LBi:UBi,LBj:UBj,:)
1006# endif
1007!
1008! Local variable declarations.
1009!
1010 integer :: i, j
1011!
1012# include "set_bounds.h"
1013!
1014!-----------------------------------------------------------------------
1015! Initialize fast-time averaged free-surface (Zt_avg1) with the inital
1016! free-surface.
1017!-----------------------------------------------------------------------
1018!
1019 DO j=jstrt,jendt
1020 DO i=istrt,iendt
1021!> Zt_avg1(i,j)=zeta(i,j,kstp)
1022!>
1023 tl_zt_avg1(i,j)=tl_zeta(i,j,kstp)
1024 END DO
1025 END DO
1026!
1027 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
1028 CALL exchange_r2d_tile (ng, tile, &
1029 & lbi, ubi, lbj, ubj, &
1030 & tl_zt_avg1)
1031 END IF
1032
1033# ifdef DISTRIBUTE
1034 CALL mp_exchange2d (ng, tile, model, 1, &
1035 & lbi, ubi, lbj, ubj, &
1036 & nghostpoints, &
1037 & ewperiodic(ng), nsperiodic(ng), &
1038 & tl_zt_avg1)
1039# endif
1040!
1041 RETURN
1042 END SUBROUTINE rp_set_zeta_timeavg_tile
1043# endif
1044#endif
1045 END MODULE rp_ini_fields_mod
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)
subroutine exchange_u3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
subroutine exchange_r3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
subroutine exchange_v3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
type(t_coupling), dimension(:), allocatable coupling
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer isvbar
integer isfsur
integer isubar
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
type(t_lbc), dimension(:,:,:), allocatable tl_lbc
Definition mod_param.F:379
integer nghostpoints
Definition mod_param.F:710
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
logical, dimension(:,:), allocatable ltracerclm
logical, dimension(:,:), allocatable lnudgetclm
type(t_sedbed), dimension(:), allocatable sedbed
Definition sedbed_mod.h:157
integer, parameter ithck
integer, dimension(:), allocatable kstp
integer, dimension(:), allocatable knew
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable krhs
integer, dimension(:), allocatable nstp
subroutine mp_exchange4d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, lbt, ubt, nghost, ew_periodic, ns_periodic, a, b, c)
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine, public rp_ini_fields(ng, tile, model)
subroutine, public rp_set_zeta_timeavg(ng, tile, model)
subroutine, public rp_ini_zeta(ng, tile, model)
subroutine rp_ini_zeta_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, kstp, krhs, knew, rmask, h, tl_bed, tl_bed_thick0, tl_bed_thick, tl_zt_avg1, zeta, tl_zeta)
subroutine rp_ini_fields_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, kstp, krhs, knew, nstp, nnew, rmask, umask, vmask, hz, tl_hz, tl_t, u, tl_u, v, tl_v, ubar, tl_ubar, vbar, tl_vbar, zeta, tl_zeta)
subroutine rp_set_zeta_timeavg_tile(ng, tile, model, lbi, ubi, lbj, ubj, kstp, tl_zt_avg1, tl_zeta)
subroutine, public rp_t3dbc_tile(ng, tile, itrc, ic, lbi, ubi, lbj, ubj, ubk, ubt, imins, imaxs, jmins, jmaxs, nstp, nout, tl_t)
Definition rp_t3dbc_im.F:58
subroutine, public rp_u2dbc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, krhs, kstp, kout, ubar, vbar, zeta, tl_ubar, tl_vbar, tl_zeta)
Definition rp_u2dbc_im.F:64
subroutine, public rp_u3dbc_tile(ng, tile, lbi, ubi, lbj, ubj, ubk, imins, imaxs, jmins, jmaxs, nstp, nout, tl_u)
Definition rp_u3dbc_im.F:58
subroutine, public rp_v2dbc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, krhs, kstp, kout, ubar, vbar, zeta, tl_ubar, tl_vbar, tl_zeta)
Definition rp_v2dbc_im.F:64
subroutine, public rp_v3dbc_tile(ng, tile, lbi, ubi, lbj, ubj, ubk, imins, imaxs, jmins, jmaxs, nstp, nout, tl_v)
Definition rp_v3dbc_im.F:59
subroutine, public rp_zetabc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, krhs, kstp, kout, zeta, tl_zeta)
Definition rp_zetabc.F:59
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