ROMS
Loading...
Searching...
No Matches
tl_set_depth.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#ifdef TANGENT
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 computes the time evolving depths of the model grid !
13! and its associated vertical transformation metric (thickness). !
14! !
15! Currently, two vertical coordinate transformations are available !
16! with various possible vertical stretching, C(s), functions, (see !
17! routine "set_scoord.F" for details). !
18! !
19! BASIC STATE variables needed: NONE !
20! Independent Variables: tl_Hz, tl_z_r, tl_z_w !
21! !
22!=======================================================================
23!
24 implicit none
25!
26 PRIVATE
27# ifdef SOLVE3D
29# ifdef ADJUST_BOUNDARY
30 PUBLIC :: tl_set_depth_bry
31# endif
32# endif
33 PUBLIC :: tl_bath, tl_bath_tile
34!
35 CONTAINS
36
37# ifdef SOLVE3D
38!
39!***********************************************************************
40 SUBROUTINE tl_set_depth (ng, tile, model)
41!***********************************************************************
42!
43 USE mod_param
44 USE mod_coupling
45 USE mod_grid
46 USE mod_ocean
47 USE mod_stepping
48!
49! Imported variable declarations.
50!
51 integer, intent(in) :: ng, tile, model
52!
53! Local variable declarations.
54!
55 character (len=*), parameter :: myfile = &
56 & __FILE__
57!
58# include "tile.h"
59!
60# ifdef PROFILE
61 CALL wclock_on (ng, model, 12, __line__, myfile)
62# endif
63 CALL tl_set_depth_tile (ng, tile, model, &
64 & lbi, ubi, lbj, ubj, &
65 & imins, imaxs, jmins, jmaxs, &
66 & nstp(ng), nnew(ng), &
67 & grid(ng) % h, &
68 & grid(ng) % tl_h, &
69# ifdef ICESHELF
70 & grid(ng) % zice, &
71# endif
72 & coupling(ng) % Zt_avg1, &
73 & coupling(ng) % tl_Zt_avg1, &
74 & grid(ng) % tl_Hz, &
75 & grid(ng) % tl_z_r, &
76 & grid(ng) % tl_z_w)
77# ifdef PROFILE
78 CALL wclock_off (ng, model, 12, __line__, myfile)
79# endif
80!
81 RETURN
82 END SUBROUTINE tl_set_depth
83!
84!***********************************************************************
85 SUBROUTINE tl_set_depth_tile (ng, tile, model, &
86 & LBi, UBi, LBj, UBj, &
87 & IminS, ImaxS, JminS, JmaxS, &
88 & nstp, nnew, &
89 & h, tl_h, &
90# ifdef ICESHELF
91 & zice, &
92# endif
93 & Zt_avg1, tl_Zt_avg1, &
94 & tl_Hz, tl_z_r, tl_z_w)
95!***********************************************************************
96!
97 USE mod_param
98 USE mod_scalars
99!
102# ifdef DISTRIBUTE
104# endif
105!
106! Imported variable declarations.
107!
108 integer, intent(in) :: ng, tile, model
109 integer, intent(in) :: lbi, ubi, lbj, ubj
110 integer, intent(in) :: imins, imaxs, jmins, jmaxs
111 integer, intent(in) :: nstp, nnew
112!
113# ifdef ASSUMED_SHAPE
114 real(r8), intent(in) :: h(lbi:,lbj:)
115# ifdef ICESHELF
116 real(r8), intent(in) :: zice(lbi:,lbj:)
117# endif
118 real(r8), intent(in) :: zt_avg1(lbi:,lbj:)
119 real(r8), intent(in) :: tl_zt_avg1(lbi:,lbj:)
120 real(r8), intent(inout) :: tl_h(lbi:,lbj:)
121 real(r8), intent(out) :: tl_hz(lbi:,lbj:,:)
122 real(r8), intent(out) :: tl_z_r(lbi:,lbj:,:)
123 real(r8), intent(out) :: tl_z_w(lbi:,lbj:,0:)
124# else
125 real(r8), intent(in) :: h(lbi:ubi,lbj:ubj)
126# ifdef ICESHELF
127 real(r8), intent(in) :: zice(lbi:ubi,lbj:ubj)
128# endif
129 real(r8), intent(in) :: zt_avg1(lbi:ubi,lbj:ubj)
130 real(r8), intent(in) :: tl_zt_avg1(lbi:ubi,lbj:ubj)
131 real(r8), intent(inout) :: tl_h(lbi:ubi,lbj:ubj)
132 real(r8), intent(out) :: tl_hz(lbi:ubi,lbj:ubj,n(ng))
133 real(r8), intent(out) :: tl_z_r(lbi:ubi,lbj:ubj,n(ng))
134 real(r8), intent(out) :: tl_z_w(lbi:ubi,lbj:ubj,0:n(ng))
135# endif
136!
137! Local variable declarations.
138!
139 integer :: i, j, k
140
141 real(r8) :: cff, cff_r, cff1_r, cff2_r, cff_w, cff1_w, cff2_w
142 real(r8) :: hinv, hwater, z_r0, z_w0
143 real(r8) :: tl_cff2_r, tl_cff2_w
144 real(r8) :: tl_hinv, tl_hwater, tl_z_r0, tl_z_w0
145
146# ifdef WET_DRY
147 real(r8), parameter :: eps = 1.0e-14_r8
148# endif
149
150# include "set_bounds.h"
151!
152!-----------------------------------------------------------------------
153! Original formulation: Compute vertical depths (meters, negative) at
154! RHO- and W-points, and vertical grid
155! thicknesses. Various stretching functions are possible.
156!
157! z_w(x,y,s,t) = Zo_w + zeta(x,y,t) * [1.0 + Zo_w / h(x,y)]
158!
159! Zo_w = hc * [s(k) - C(k)] + C(k) * h(x,y)
160!
161!-----------------------------------------------------------------------
162!
163 IF (vtransform(ng).eq.1) THEN
164 DO j=jstrt,jendt
165 DO i=istrt,iendt
166# if defined WET_DRY
167 IF (h(i,j).eq.0.0_r8) THEN
168!^ h(i,j)=eps
169!^
170 tl_h(i,j)=0.0_r8
171 END IF
172# endif
173!^ z_w(i,j,0)=-h(i,j)
174!^
175 tl_z_w(i,j,0)=-tl_h(i,j)
176 END DO
177 DO k=1,n(ng)
178 cff_r=hc(ng)*(scalars(ng)%sc_r(k)-scalars(ng)%Cs_r(k))
179 cff_w=hc(ng)*(scalars(ng)%sc_w(k)-scalars(ng)%Cs_w(k))
180 cff1_r=scalars(ng)%Cs_r(k)
181 cff1_w=scalars(ng)%Cs_w(k)
182 DO i=istrt,iendt
183 hwater=h(i,j)
184# ifdef ICESHELF
185 hwater=hwater-abs(zice(i,j))
186# endif
187 tl_hwater=tl_h(i,j)
188 hinv=1.0_r8/hwater
189 tl_hinv=-hinv*hinv*tl_hwater
190
191 z_w0=cff_w+cff1_w*hwater
192 tl_z_w0=cff1_w*tl_hwater
193
194!^ z_w(i,j,k)=z_w0+Zt_avg1(i,j)*(1.0_r8+z_w0*hinv)
195!^
196 tl_z_w(i,j,k)=tl_z_w0+ &
197 & tl_zt_avg1(i,j)*(1.0_r8+z_w0*hinv)+ &
198 & zt_avg1(i,j)*(tl_z_w0*hinv+z_w0*tl_hinv)
199
200 z_r0=cff_r+cff1_r*hwater
201 tl_z_r0=cff1_r*tl_hwater
202
203!^ z_r(i,j,k)=z_r0+Zt_avg1(i,j)*(1.0_r8+z_r0*hinv)
204!^
205 tl_z_r(i,j,k)=tl_z_r0+ &
206 & tl_zt_avg1(i,j)*(1.0_r8+z_r0*hinv)+ &
207 & zt_avg1(i,j)*(tl_z_r0*hinv+z_r0*tl_hinv)
208# ifdef ICESHELF
209!^ z_w(i,j,k)=z_w(i,j,k)-ABS(zice(i,j))
210!^ z_r(i,j,k)=z_r(i,j,k)-ABS(zice(i,j))
211# endif
212!^ Hz(i,j,k)=z_w(i,j,k)-z_w(i,j,k-1)
213!^
214 tl_hz(i,j,k)=tl_z_w(i,j,k)-tl_z_w(i,j,k-1)
215 END DO
216 END DO
217 END DO
218!
219!-----------------------------------------------------------------------
220! New formulation: Compute vertical depths (meters, negative) at
221! RHO- and W-points, and vertical grid thicknesses.
222! Various stretching functions are possible.
223!
224! z_w(x,y,s,t) = zeta(x,y,t) + [zeta(x,y,t)+ h(x,y)] * Zo_w
225!
226! Zo_w = [hc * s(k) + C(k) * h(x,y)] / [hc + h(x,y)]
227!
228!-----------------------------------------------------------------------
229!
230 ELSE IF (vtransform(ng).eq.2) THEN
231 DO j=jstrt,jendt
232 DO i=istrt,iendt
233# if defined WET_DRY
234 IF (h(i,j).eq.0.0_r8) THEN
235!^ h(i,j)=eps
236!^
237 tl_h(i,j)=0.0_r8
238 END IF
239# endif
240!^ z_w(i,j,0)=-h(i,j)
241!^
242 tl_z_w(i,j,0)=-tl_h(i,j)
243 END DO
244 DO k=1,n(ng)
245 cff_r=hc(ng)*scalars(ng)%sc_r(k)
246 cff_w=hc(ng)*scalars(ng)%sc_w(k)
247 cff1_r=scalars(ng)%Cs_r(k)
248 cff1_w=scalars(ng)%Cs_w(k)
249 DO i=istrt,iendt
250 hwater=h(i,j)
251# ifdef ICESHELF
252 hwater=hwater-abs(zice(i,j))
253# endif
254 tl_hwater=tl_h(i,j)
255 hinv=1.0_r8/(hc(ng)+hwater)
256 tl_hinv=-hinv*hinv*tl_hwater
257
258 cff2_w=(cff_w+cff1_w*hwater)*hinv
259 tl_cff2_w=cff1_w*tl_hwater*hinv+ &
260 & (cff_w+cff1_w*hwater)*tl_hinv
261
262!^ z_w(i,j,k)=Zt_avg1(i,j)+ &
263!^ & (Zt_avg1(i,j)+hwater)*cff2_w
264!^
265 tl_z_w(i,j,k)=tl_zt_avg1(i,j)+ &
266 & (tl_zt_avg1(i,j)+tl_hwater)*cff2_w+ &
267 & (zt_avg1(i,j)+hwater)*tl_cff2_w
268
269 cff2_r=(cff_r+cff1_r*hwater)*hinv
270 tl_cff2_r=cff1_r*tl_hwater*hinv+ &
271 & (cff_r+cff1_r*hwater)*tl_hinv
272
273!^ z_r(i,j,k)=Zt_avg1(i,j)+ &
274!^ & (Zt_avg1(i,j)+hwater)*cff2_r
275!^
276 tl_z_r(i,j,k)=tl_zt_avg1(i,j)+ &
277 & (tl_zt_avg1(i,j)+tl_hwater)*cff2_r+ &
278 & (zt_avg1(i,j)+hwater)*tl_cff2_r
279
280# ifdef ICESHELF
281!^ z_w(i,j,k)=z_w(i,j,k)-ABS(zice(i,j))
282!^ z_r(i,j,k)=z_r(i,j,k)-ABS(zice(i,j))
283# endif
284!^ Hz(i,j,k)=z_w(i,j,k)-z_w(i,j,k-1)
285!^
286 tl_hz(i,j,k)=tl_z_w(i,j,k)-tl_z_w(i,j,k-1)
287 END DO
288 END DO
289 END DO
290 END IF
291!
292!-----------------------------------------------------------------------
293! Exchange boundary information.
294!-----------------------------------------------------------------------
295!
296 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
297!^ CALL exchange_r2d_tile (ng, tile, &
298!^ & LBi, UBi, LBj, UBj, &
299!^ & h)
300!^
301 CALL exchange_r2d_tile (ng, tile, &
302 & lbi, ubi, lbj, ubj, &
303 & tl_h)
304!^ CALL exchange_w3d_tile (ng, tile, &
305!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
306!^ & z_w)
307!^
308 CALL exchange_w3d_tile (ng, tile, &
309 & lbi, ubi, lbj, ubj, 0, n(ng), &
310 & tl_z_w)
311!^ CALL exchange_r3d_tile (ng, tile, &
312!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
313!^ & z_r)
314!^
315 CALL exchange_r3d_tile (ng, tile, &
316 & lbi, ubi, lbj, ubj, 1, n(ng), &
317 & tl_z_r)
318!^ CALL exchange_r3d_tile (ng, tile, &
319!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
320!^ & Hz)
321!^
322 CALL exchange_r3d_tile (ng, tile, &
323 & lbi, ubi, lbj, ubj, 1, n(ng), &
324 & tl_hz)
325 END IF
326
327# ifdef DISTRIBUTE
328!^ CALL mp_exchange2d (ng, tile, model, 1, &
329!^ & LBi, UBi, LBj, UBj, &
330!^ & NghostPoints, &
331!^ & EWperiodic(ng), NSperiodic(ng), &
332!^ & h)
333!^
334 CALL mp_exchange2d (ng, tile, model, 1, &
335 & lbi, ubi, lbj, ubj, &
336 & nghostpoints, &
337 & ewperiodic(ng), nsperiodic(ng), &
338 & tl_h)
339!^ CALL mp_exchange3d (ng, tile, model, 1, &
340!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
341!^ & NghostPoints, &
342!^ & EWperiodic(ng), NSperiodic(ng), &
343!^ & z_w)
344!^
345 CALL mp_exchange3d (ng, tile, model, 1, &
346 & lbi, ubi, lbj, ubj, 0, n(ng), &
347 & nghostpoints, &
348 & ewperiodic(ng), nsperiodic(ng), &
349 & tl_z_w)
350!^ CALL mp_exchange3d (ng, tile, model, 2, &
351!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
352!^ & NghostPoints, &
353!^ & EWperiodic(ng), NSperiodic(ng), &
354!^ & z_r, Hz)
355!^
356 CALL mp_exchange3d (ng, tile, model, 2, &
357 & lbi, ubi, lbj, ubj, 1, n(ng), &
358 & nghostpoints, &
359 & ewperiodic(ng), nsperiodic(ng), &
360 & tl_z_r, tl_hz)
361# endif
362!
363 RETURN
364 END SUBROUTINE tl_set_depth_tile
365
366# ifdef ADJUST_BOUNDARY
367!
368!***********************************************************************
369 SUBROUTINE tl_set_depth_bry (ng, tile, model)
370!***********************************************************************
371!
372 USE mod_param
373 USE mod_grid
374!
375! Imported variable declarations.
376!
377 integer, intent(in) :: ng, tile, model
378!
379! Local variable declarations.
380!
381 character (len=*), parameter :: myfile = &
382 & __FILE__//", tl_set_depth_bry"
383!
384# include "tile.h"
385!
386# ifdef PROFILE
387 CALL wclock_on (ng, model, 12, __line__, myfile)
388# endif
389 CALL tl_set_depth_bry_tile (ng, tile, model, &
390 & lbi, ubi, lbj, ubj, lbij, ubij, &
391 & imins, imaxs, jmins, jmaxs, &
392 & grid(ng) % h, &
393 & grid(ng) % tl_h, &
394# ifdef ICESHELF
395 & grid(ng) % zice, &
396# endif
397 & grid(ng) % tl_Hz_bry)
398# ifdef PROFILE
399 CALL wclock_off (ng, model, 12, __line__, myfile)
400# endif
401!
402 RETURN
403 END SUBROUTINE tl_set_depth_bry
404!
405!***********************************************************************
406 SUBROUTINE tl_set_depth_bry_tile (ng, tile, model, &
407 & LBi, UBi, LBj, UBj, LBij, UBij, &
408 & IminS, ImaxS, JminS, JmaxS, &
409 & h, tl_h, &
410# ifdef ICESHELF
411 & zice, &
412# endif
413 & tl_Hz_bry)
414!***********************************************************************
415!
416 USE mod_param
417 USE mod_boundary
418 USE mod_ncparam
419 USE mod_scalars
420!
421# ifdef DISTRIBUTE
423# endif
424!
425! Imported variable declarations.
426!
427 integer, intent(in) :: ng, tile, model
428 integer, intent(in) :: LBi, UBi, LBj, UBj, LBij, UBij
429 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
430!
431# ifdef ASSUMED_SHAPE
432 real(r8), intent(in) :: h(LBi:,LBj:)
433 real(r8), intent(in) :: tl_h(LBi:,LBj:)
434# ifdef ICESHELF
435 real(r8), intent(in) :: zice(LBi:,LBj:)
436# endif
437 real(r8), intent(out) :: tl_Hz_bry(LBij:,:,:)
438# else
439 real(r8), intent(in) :: h(LBi:UBi,LBj:UBj)
440 real(r8), intent(in) :: tl_h(LBi:UBi,LBj:UBj)
441# ifdef ICESHELF
442 real(r8), intent(in) :: zice(LBi:UBi,LBj:UBj)
443# endif
444 real(r8), intent(out) :: tl_Hz_bry(LBij:UBij,N(ng),4)
445# endif
446!
447! Local variable declarations.
448!
449 integer :: i, ibry, j, k
450
451 real(r8) :: cff_w, cff1_w, cff2_w
452 real(r8) :: hinv, hwater, z_w0
453 real(r8) :: tl_cff2_w, tl_hinv, tl_hwater, tl_z_w0
454
455 real(r8), dimension(0:N(ng)) :: tl_Zw
456
457# include "set_bounds.h"
458!
459!-----------------------------------------------------------------------
460! Original formulation: Compute vertical depths (meters, negative) at
461! RHO- and W-points, and vertical grid
462! thicknesses. Various stretching functions are possible.
463!
464! z_w(x,y,s,t) = Zo_w + zeta(x,y,t) * [1.0 + Zo_w / h(x,y)]
465!
466! Zo_w = hc * [s(k) - C(k)] + C(k) * h(x,y)
467!
468!-----------------------------------------------------------------------
469!
470 IF (vtransform(ng).eq.1) THEN
471
472 IF (tl_lbc(iwest,isfsur,ng)%acquire.and. &
473 & domain(ng)%Western_Edge(tile)) THEN
474 i=bounds(ng)%edge(iwest,r2dvar)
475 DO j=jstrt,jendt
476 hwater=h(i,j)
477# ifdef ICESHELF
478 hwater=hwater-abs(zice(i,j))
479# endif
480 tl_hwater=tl_h(i,j)
481 hinv=1.0_r8/hwater
482 tl_hinv=-hinv*hinv*tl_hwater
483!^ Zw(0)=-h(i,j)
484!^
485 tl_zw(0)=-tl_h(i,j)
486 DO k=1,n(ng)
487 cff_w=hc(ng)*(scalars(ng)%sc_w(k)-scalars(ng)%Cs_w(k))
488 cff1_w=scalars(ng)%Cs_w(k)
489 z_w0=cff_w+cff1_w*hwater
490 tl_z_w0=cff1_w*tl_hwater
491!^ Zw(k)=z_w0+BOUNDARY(ng)%zeta_west(j)*(1.0_r8+z_w0*hinv)
492!^
493 tl_zw(k)=tl_z_w0+ &
494 & boundary(ng)%tl_zeta_west(j)* &
495 & (1.0_r8+z_w0*hinv)+ &
496 & boundary(ng)%zeta_west(j)* &
497 & (tl_z_w0*hinv+z_w0*tl_hinv)
498# ifdef ICESHELF
499!^ Zw(k)=Zw(k)-ABS(zice(i,j))
500# endif
501!^ Hz_bry(j,k,iwest)=Zw(k)-Zw(k-1)
502!^
503 tl_hz_bry(j,k,iwest)=tl_zw(k)-tl_zw(k-1)
504 END DO
505 END DO
506 END IF
507
508 IF (tl_lbc(ieast,isfsur,ng)%acquire.and. &
509 & domain(ng)%Eastern_Edge(tile)) THEN
510 i=bounds(ng)%edge(ieast,r2dvar)
511 DO j=jstrt,jendt
512 hwater=h(i,j)
513# ifdef ICESHELF
514 hwater=hwater-abs(zice(i,j))
515# endif
516 tl_hwater=tl_h(i,j)
517 hinv=1.0_r8/hwater
518 tl_hinv=-hinv*hinv*tl_hwater
519!^ Zw(0)=-h(i,j)
520!^
521 tl_zw(0)=-tl_h(i,j)
522 DO k=1,n(ng)
523 cff_w=hc(ng)*(scalars(ng)%sc_w(k)-scalars(ng)%Cs_w(k))
524 cff1_w=scalars(ng)%Cs_w(k)
525 z_w0=cff_w+cff1_w*hwater
526 tl_z_w0=cff1_w*tl_hwater
527!^ Zw(k)=z_w0+BOUNDARY(ng)%zeta_east(j)*(1.0_r8+z_w0*hinv)
528!^
529 tl_zw(k)=tl_z_w0+ &
530 & boundary(ng)%tl_zeta_east(j)* &
531 & (1.0_r8+z_w0*hinv)+ &
532 & boundary(ng)%zeta_east(j)* &
533 & (tl_z_w0*hinv+z_w0*tl_hinv)
534# ifdef ICESHELF
535!^ Zw(k)=Zw(k)-ABS(zice(i,j))
536# endif
537!^ Hz_bry(j,k,ieast)=Zw(k)-Zw(k-1)
538!^
539 tl_hz_bry(j,k,ieast)=tl_zw(k)-tl_zw(k-1)
540 END DO
541 END DO
542 END IF
543
544 IF (tl_lbc(isouth,isfsur,ng)%acquire.and. &
545 & domain(ng)%Southern_Edge(tile)) THEN
546 j=bounds(ng)%edge(isouth,r2dvar)
547 DO i=istrt,iendt
548 hwater=h(i,j)
549# ifdef ICESHELF
550 hwater=hwater-abs(zice(i,j))
551# endif
552 tl_hwater=tl_h(i,j)
553 hinv=1.0_r8/hwater
554 tl_hinv=-hinv*hinv*tl_hwater
555!^ Zw(0)=-h(i,j)
556!^
557 tl_zw(0)=-tl_h(i,j)
558 DO k=1,n(ng)
559 cff_w=hc(ng)*(scalars(ng)%sc_w(k)-scalars(ng)%Cs_w(k))
560 cff1_w=scalars(ng)%Cs_w(k)
561 z_w0=cff_w+cff1_w*hwater
562 tl_z_w0=cff1_w*tl_hwater
563!^ Zw(k)=z_w0+BOUNDARY(ng)%zeta_south(i)*(1.0_r8+z_w0*hinv)
564!^
565 tl_zw(k)=tl_z_w0+ &
566 & boundary(ng)%tl_zeta_south(i)* &
567 & (1.0_r8+z_w0*hinv)+ &
568 & boundary(ng)%zeta_south(i)* &
569 & (tl_z_w0*hinv+z_w0*tl_hinv)
570# ifdef ICESHELF
571!^ Zw(k)=Zw(k)-ABS(zice(i,j))
572# endif
573!^ Hz_bry(i,k,isouth)=Zw(k)-Zw(k-1)
574!^
575 tl_hz_bry(i,k,isouth)=tl_zw(k)-tl_zw(k-1)
576 END DO
577 END DO
578 END IF
579
580 IF (tl_lbc(inorth,isfsur,ng)%acquire.and. &
581 & domain(ng)%Northern_Edge(tile)) THEN
582 j=bounds(ng)%edge(inorth,r2dvar)
583 DO i=istrt,iendt
584 hwater=h(i,j)
585# ifdef ICESHELF
586 hwater=hwater-abs(zice(i,j))
587# endif
588 tl_hwater=tl_h(i,j)
589 hinv=1.0_r8/hwater
590 tl_hinv=-hinv*hinv*tl_hwater
591!^ Zw(0)=-h(i,j)
592!^
593 tl_zw(0)=-tl_h(i,j)
594 DO k=1,n(ng)
595 cff_w=hc(ng)*(scalars(ng)%sc_w(k)-scalars(ng)%Cs_w(k))
596 cff1_w=scalars(ng)%Cs_w(k)
597 z_w0=cff_w+cff1_w*hwater
598 tl_z_w0=cff1_w*tl_hwater
599!^ Zw(k)=z_w0+BOUNDARY(ng)%zeta_north(i)*(1.0_r8+z_w0*hinv)
600!^
601 tl_zw(k)=tl_z_w0+ &
602 & boundary(ng)%tl_zeta_north(i)* &
603 & (1.0_r8+z_w0*hinv)+ &
604 & boundary(ng)%zeta_north(i)* &
605 & (tl_z_w0*hinv+z_w0*tl_hinv)
606# ifdef ICESHELF
607!^ Zw(k)=Zw(k)-ABS(zice(i,j))
608# endif
609!^ Hz_bry(i,k,inorth)=Zw(k)-Zw(k-1)
610!^
611 tl_hz_bry(i,k,inorth)=tl_zw(k)-tl_zw(k-1)
612 END DO
613 END DO
614 END IF
615!
616!-----------------------------------------------------------------------
617! New formulation: Compute vertical depths (meters, negative) at
618! RHO- and W-points, and vertical grid thicknesses.
619! Various stretching functions are possible.
620!
621! z_w(x,y,s,t) = zeta(x,y,t) + [zeta(x,y,t)+ h(x,y)] * Zo_w
622!
623! Zo_w = [hc * s(k) + C(k) * h(x,y)] / [hc + h(x,y)]
624!
625!-----------------------------------------------------------------------
626!
627 ELSE IF (vtransform(ng).eq.2) THEN
628
629 IF (tl_lbc(iwest,isfsur,ng)%acquire.and. &
630 & domain(ng)%Western_Edge(tile)) THEN
631 i=bounds(ng)%edge(iwest,r2dvar)
632 DO j=jstrt,jendt
633 hwater=h(i,j)
634# ifdef ICESHELF
635 hwater=hwater-abs(zice(i,j))
636# endif
637 tl_hwater=tl_h(i,j)
638 hinv=1.0_r8/(hc(ng)+hwater)
639 tl_hinv=-hinv*hinv*tl_hwater
640!^ Zw(0)=-h(i,j)
641!^
642 tl_zw(0)=-tl_h(i,j)
643 DO k=1,n(ng)
644 cff_w=hc(ng)*scalars(ng)%sc_w(k)
645 cff1_w=scalars(ng)%Cs_w(k)
646 cff2_w=(cff_w+cff1_w*hwater)*hinv
647 tl_cff2_w=cff1_w*tl_hwater*hinv+ &
648 & (cff_w+cff1_w*hwater)*tl_hinv
649!^ Zw(k)=BOUNDARY(ng)%zeta_west(j)+ &
650!^ & (BOUNDARY(ng)%zeta_west(j)+hwater)*cff2_w
651!^
652 tl_zw(k)=boundary(ng)%tl_zeta_west(j)+ &
653 & (boundary(ng)%tl_zeta_west(j)+ &
654 & tl_hwater)*cff2_w+ &
655 & (boundary(ng)%zeta_west(j)+ &
656 & hwater)*tl_cff2_w
657# ifdef ICESHELF
658!^ Zw(k)=Zw(k)-ABS(zice(i,j))
659# endif
660!^ Hz_bry(j,k,iwest)=Zw(k)-Zw(k-1)
661!^
662 tl_hz_bry(j,k,iwest)=tl_zw(k)-tl_zw(k-1)
663 END DO
664 END DO
665 END IF
666
667 IF (tl_lbc(ieast,isfsur,ng)%acquire.and. &
668 & domain(ng)%Eastern_Edge(tile)) THEN
669 i=bounds(ng)%edge(ieast,r2dvar)
670 DO j=jstrt,jendt
671 hwater=h(i,j)
672# ifdef ICESHELF
673 hwater=hwater-abs(zice(i,j))
674# endif
675 tl_hwater=tl_h(i,j)
676 hinv=1.0_r8/(hc(ng)+hwater)
677 tl_hinv=-hinv*hinv*tl_hwater
678!^ Zw(0)=-h(i,j)
679!^
680 tl_zw(0)=-tl_h(i,j)
681 DO k=1,n(ng)
682 cff_w=hc(ng)*scalars(ng)%sc_w(k)
683 cff1_w=scalars(ng)%Cs_w(k)
684 cff2_w=(cff_w+cff1_w*hwater)*hinv
685 tl_cff2_w=cff1_w*tl_hwater*hinv+ &
686 & (cff_w+cff1_w*hwater)*tl_hinv
687!^ Zw(k)=BOUNDARY(ng)%zeta_east(j)+ &
688!^ & (BOUNDARY(ng)%zeta_east(j)+hwater)*cff2_w
689!^
690 tl_zw(k)=boundary(ng)%tl_zeta_east(j)+ &
691 & (boundary(ng)%tl_zeta_east(j)+ &
692 & tl_hwater)*cff2_w+ &
693 & (boundary(ng)%zeta_east(j)+ &
694 & hwater)*tl_cff2_w
695# ifdef ICESHELF
696!^ Zw(k)=Zw(k)-ABS(zice(i,j))
697# endif
698!^ Hz_bry(j,k,ieast)=Zw(k)-Zw(k-1)
699!^
700 tl_hz_bry(j,k,ieast)=tl_zw(k)-tl_zw(k-1)
701 END DO
702 END DO
703 END IF
704
705 IF (tl_lbc(isouth,isfsur,ng)%acquire.and. &
706 & domain(ng)%Southern_Edge(tile)) THEN
707 j=bounds(ng)%edge(isouth,r2dvar)
708 DO i=istrt,iendt
709 hwater=h(i,j)
710# ifdef ICESHELF
711 hwater=hwater-abs(zice(i,j))
712# endif
713 tl_hwater=tl_h(i,j)
714 hinv=1.0_r8/(hc(ng)+hwater)
715 tl_hinv=-hinv*hinv*tl_hwater
716!^ Zw(0)=-h(i,j)
717!^
718 tl_zw(0)=-tl_h(i,j)
719 DO k=1,n(ng)
720 cff_w=hc(ng)*scalars(ng)%sc_w(k)
721 cff1_w=scalars(ng)%Cs_w(k)
722 cff2_w=(cff_w+cff1_w*hwater)*hinv
723 tl_cff2_w=cff1_w*tl_hwater*hinv+ &
724 & (cff_w+cff1_w*hwater)*tl_hinv
725!^ Zw(k)=BOUNDARY(ng)%zeta_south(i)+ &
726!^ & (BOUNDARY(ng)%zeta_south(i)+hwater)*cff2_w
727!^
728 tl_zw(k)=boundary(ng)%tl_zeta_south(i)+ &
729 & (boundary(ng)%tl_zeta_south(i)+ &
730 & tl_hwater)*cff2_w+ &
731 & (boundary(ng)%zeta_south(i)+ &
732 & hwater)*tl_cff2_w
733# ifdef ICESHELF
734!^ Zw(k)=Zw(k)-ABS(zice(i,j))
735# endif
736!^ Hz_bry(i,k,isouth)=Zw(k)-Zw(k-1)
737!^
738 tl_hz_bry(i,k,isouth)=tl_zw(k)-tl_zw(k-1)
739 END DO
740 END DO
741 END IF
742
743 IF (tl_lbc(inorth,isfsur,ng)%acquire.and. &
744 & domain(ng)%Northern_Edge(tile)) THEN
745 j=bounds(ng)%edge(inorth,r2dvar)
746 DO i=istrt,iendt
747 hwater=h(i,j)
748# ifdef ICESHELF
749 hwater=hwater-abs(zice(i,j))
750# endif
751 tl_hwater=tl_h(i,j)
752 hinv=1.0_r8/(hc(ng)+hwater)
753 tl_hinv=-hinv*hinv*tl_hwater
754!^ Zw(0)=-h(i,j)
755!^
756 tl_zw(0)=-tl_h(i,j)
757 DO k=1,n(ng)
758 cff_w=hc(ng)*scalars(ng)%sc_w(k)
759 cff1_w=scalars(ng)%Cs_w(k)
760 cff2_w=(cff_w+cff1_w*hwater)*hinv
761 tl_cff2_w=cff1_w*tl_hwater*hinv+ &
762 & (cff_w+cff1_w*hwater)*tl_hinv
763!^ Zw(k)=BOUNDARY(ng)%zeta_north(i)+ &
764!^ & (BOUNDARY(ng)%zeta_north(i)+hwater)*cff2_w
765!^
766 tl_zw(k)=boundary(ng)%tl_zeta_north(i)+ &
767 & (boundary(ng)%tl_zeta_north(i)+ &
768 & tl_hwater)*cff2_w+ &
769 & (boundary(ng)%zeta_north(i)+ &
770 & hwater)*tl_cff2_w
771# ifdef ICESHELF
772!^ Zw(k)=Zw(k)-ABS(zice(i,j))
773# endif
774!^ Hz_bry(i,k,inorth)=Zw(k)-Zw(k-1)
775!^
776 tl_hz_bry(i,k,inorth)=tl_zw(k)-tl_zw(k-1)
777 END DO
778 END DO
779 END IF
780 END IF
781
782# ifdef DISTRIBUTE
783!
784!-----------------------------------------------------------------------
785! Exchange boundary information.
786!-----------------------------------------------------------------------
787!
788 DO ibry=1,4
789!^ CALL mp_exchange3d_bry (ng, tile, model, 1, ibry, &
790!^ & LBij, UBij, 1, N(ng), &
791!^ & NghostPoints, &
792!^ & EWperiodic(ng), NSperiodic(ng), &
793!^ & Hz_bry(:,:,ibry))
794!^
795 CALL mp_exchange3d_bry (ng, tile, model, 1, ibry, &
796 & lbij, ubij, 1, n(ng), &
797 & nghostpoints, &
798 & ewperiodic(ng), nsperiodic(ng), &
799 & tl_hz_bry(:,:,ibry))
800 END DO
801# endif
802!
803 RETURN
804 END SUBROUTINE tl_set_depth_bry_tile
805# endif
806# endif
807!
808!***********************************************************************
809 SUBROUTINE tl_bath (ng, tile)
810!***********************************************************************
811!
812 USE mod_param
813 USE mod_grid
814!
815! Imported variable declarations.
816!
817 integer, intent(in) :: ng, tile
818!
819! Local variable declarations.
820!
821 character (len=*), parameter :: myfile = &
822 & __FILE__//", tl_bath"
823!
824# include "tile.h"
825!
826# ifdef PROFILE
827 CALL wclock_on (ng, itlm, 12, __line__, myfile)
828# endif
829 CALL tl_bath_tile (ng, tile, &
830 & lbi, ubi, lbj, ubj, &
831 & imins, imaxs, jmins, jmaxs, &
832 & grid(ng) % tl_h)
833# ifdef PROFILE
834 CALL wclock_off (ng, itlm, 12, __line__, myfile)
835# endif
836!
837 RETURN
838 END SUBROUTINE tl_bath
839!
840!***********************************************************************
841 SUBROUTINE tl_bath_tile (ng, tile, &
842 & LBi, UBi, LBj, UBj, &
843 & IminS, ImaxS, JminS, JmaxS, &
844 & tl_h)
845!***********************************************************************
846!
847 USE mod_param
848 USE mod_scalars
849!
851# ifdef DISTRIBUTE
853# endif
854!
855! Imported variable declarations.
856!
857 integer, intent(in) :: ng, tile
858 integer, intent(in) :: lbi, ubi, lbj, ubj
859 integer, intent(in) :: imins, imaxs, jmins, jmaxs
860!
861# ifdef ASSUMED_SHAPE
862 real(r8), intent(out) :: tl_h(lbi:,lbj:)
863# else
864 real(r8), intent(out) :: tl_h(lbi:ubi,lbj:ubj)
865# endif
866!
867! Local variable declarations.
868!
869 integer :: i, j
870
871# include "set_bounds.h"
872!
873!-----------------------------------------------------------------------
874! Initialize tangent linear bathymetry to zero.
875!-----------------------------------------------------------------------
876!
877 DO j=jstrt,jendt
878 DO i=istrt,iendt
879 tl_h(i,j)=0.0_r8
880 END DO
881 END DO
882
883 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
884 CALL exchange_r2d_tile (ng, tile, &
885 & lbi, ubi, lbj, ubj, &
886 & tl_h)
887 END IF
888
889# ifdef DISTRIBUTE
890 CALL mp_exchange2d (ng, tile, itlm, 1, &
891 & lbi, ubi, lbj, ubj, &
892 & nghostpoints, &
893 & ewperiodic(ng), nsperiodic(ng), &
894 & tl_h)
895# endif
896!
897 RETURN
898 END SUBROUTINE tl_bath_tile
899#endif
900 END MODULE tl_set_depth_mod
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_r3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
subroutine exchange_w3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
type(t_boundary), dimension(:), allocatable boundary
type(t_coupling), dimension(:), allocatable coupling
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer isfsur
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
type(t_lbc), dimension(:,:,:), allocatable tl_lbc
Definition mod_param.F:379
integer nghostpoints
Definition mod_param.F:710
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, parameter itlm
Definition mod_param.F:663
integer, parameter r2dvar
Definition mod_param.F:717
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
real(dp), dimension(:), allocatable hc
integer, parameter isouth
type(t_scalars), dimension(:), allocatable scalars
Definition mod_scalars.F:65
integer, parameter ieast
integer, parameter inorth
integer, dimension(:), allocatable vtransform
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
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 mp_exchange3d_bry(ng, tile, model, nvar, boundary, lbij, ubij, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine, public tl_set_depth_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nstp, nnew, h, tl_h, zice, zt_avg1, tl_zt_avg1, tl_hz, tl_z_r, tl_z_w)
subroutine tl_set_depth_bry_tile(ng, tile, model, lbi, ubi, lbj, ubj, lbij, ubij, imins, imaxs, jmins, jmaxs, h, tl_h, zice, tl_hz_bry)
subroutine, public tl_bath_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, tl_h)
subroutine, public tl_bath(ng, tile)
subroutine, public tl_set_depth_bry(ng, tile, model)
subroutine, public tl_set_depth(ng, tile, model)
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