ROMS
Loading...
Searching...
No Matches
tl_set_vbc.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 module sets vertical boundary conditons for tangent linear !
13! momentum and tracers. !
14! !
15! BASIC STATE variables needed: stflx, dqdt, t, sss, btflx, u, v, !
16! z_r !
17! !
18! NOTE: stflx and btflx will be over written. !
19! !
20!=======================================================================
21!
22 implicit none
23!
24 PRIVATE
25 PUBLIC :: tl_set_vbc
26!
27 CONTAINS
28
29# ifdef SOLVE3D
30!
31!***********************************************************************
32 SUBROUTINE tl_set_vbc (ng, tile)
33!***********************************************************************
34!
35 USE mod_param
36 USE mod_grid
37 USE mod_forces
38 USE mod_ocean
39 USE mod_stepping
40!
41! Imported variable declarations.
42!
43 integer, intent(in) :: ng, tile
44!
45! Local variable declarations.
46!
47 character (len=*), parameter :: myfile = &
48 & __FILE__
49!
50# include "tile.h"
51!
52# ifdef PROFILE
53 CALL wclock_on (ng, itlm, 6, __line__, myfile)
54# endif
55 CALL tl_set_vbc_tile (ng, tile, &
56 & lbi, ubi, lbj, ubj, &
57 & imins, imaxs, jmins, jmaxs, &
58 & nrhs(ng), &
59 & grid(ng) % Hz, &
60 & grid(ng) % tl_Hz, &
61# if defined UV_LOGDRAG
62 & grid(ng) % ZoBot, &
63# elif defined UV_LDRAG
64 & grid(ng) % rdrag, &
65# elif defined UV_QDRAG
66 & grid(ng) % rdrag2, &
67# endif
68# if !defined BBL_MODEL_NOT_YET || defined ICESHELF
69 & grid(ng) % z_r, &
70 & grid(ng) % z_w, &
71 & grid(ng) % tl_z_r, &
72 & grid(ng) % tl_z_w, &
73# endif
74# if defined ICESHELF
75 & grid(ng) % zice, &
76# endif
77 & ocean(ng) % t, &
78 & ocean(ng) % tl_t, &
79# if !defined BBL_MODEL_NOT_YET || defined ICESHELF
80 & ocean(ng) % u, &
81 & ocean(ng) % v, &
82 & ocean(ng) % tl_u, &
83 & ocean(ng) % tl_v, &
84# endif
85# ifdef QCORRECTION
86 & forces(ng) % dqdt, &
87 & forces(ng) % sst, &
88# endif
89# if defined SCORRECTION || defined SRELAXATION
90 & forces(ng) % sss, &
91# endif
92# if defined ICESHELF
93# ifdef SHORTWAVE
94 & forces(ng) % srflx, &
95# endif
96 & forces(ng) % tl_sustr, &
97 & forces(ng) % tl_svstr, &
98# endif
99# ifndef BBL_MODEL_NOT_YET
100 & forces(ng) % tl_bustr, &
101 & forces(ng) % tl_bvstr, &
102# endif
103 & forces(ng) % stflux, &
104 & forces(ng) % btflux, &
105 & forces(ng) % stflx, &
106 & forces(ng) % btflx, &
107 & forces(ng) % tl_stflx, &
108 & forces(ng) % tl_btflx)
109# ifdef PROFILE
110 CALL wclock_off (ng, itlm, 6, __line__, myfile)
111# endif
112!
113 RETURN
114 END SUBROUTINE tl_set_vbc
115!
116!***********************************************************************
117 SUBROUTINE tl_set_vbc_tile (ng, tile, &
118 & LBi, UBi, LBj, UBj, &
119 & IminS, ImaxS, JminS, JmaxS, &
120 & nrhs, &
121 & Hz, tl_Hz, &
122# if defined UV_LOGDRAG
123 & ZoBot, &
124# elif defined UV_LDRAG
125 & rdrag, &
126# elif defined UV_QDRAG
127 & rdrag2, &
128# endif
129# if !defined BBL_MODEL_NOT_YET || defined ICESHELF
130 & z_r, z_w, &
131 & tl_z_r, tl_z_w, &
132# endif
133# if defined ICESHELF
134 & zice, &
135# endif
136 & t, tl_t, &
137# if !defined BBL_MODEL_NOT_YET || defined ICESHELF
138 & u, v, &
139 & tl_u, tl_v, &
140# endif
141# ifdef QCORRECTION
142 & dqdt, sst, &
143# endif
144# if defined SCORRECTION || defined SRELAXATION
145 & sss, &
146# endif
147# if defined ICESHELF
148# ifdef SHORTWAVE
149 & srflx, &
150# endif
151 & tl_sustr, tl_svstr, &
152# endif
153# ifndef BBL_MODEL_NOT_YET
154 & tl_bustr, tl_bvstr, &
155# endif
156 & stflux, btflux, &
157 & stflx, btflx, &
158 & tl_stflx, tl_btflx)
159!***********************************************************************
160!
161 USE mod_param
162 USE mod_scalars
163!
164 USE bc_2d_mod
165# ifdef DISTRIBUTE
167# endif
168!
169! Imported variable declarations.
170!
171 integer, intent(in) :: ng, tile
172 integer, intent(in) :: LBi, UBi, LBj, UBj
173 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
174 integer, intent(in) :: nrhs
175!
176# ifdef ASSUMED_SHAPE
177 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
178 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
179# if defined UV_LOGDRAG
180 real(r8), intent(in) :: ZoBot(LBi:,LBj:)
181# elif defined UV_LDRAG
182 real(r8), intent(in) :: rdrag(LBi:,LBj:)
183# elif defined UV_QDRAG
184 real(r8), intent(in) :: rdrag2(LBi:,LBj:)
185# endif
186# if !defined BBL_MODEL_NOT_YET || defined ICESHELF
187 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
188 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
189 real(r8), intent(in) :: tl_z_r(LBi:,LBj:,:)
190 real(r8), intent(in) :: tl_z_w(LBi:,LBj:,0:)
191# endif
192# if defined ICESHELF
193 real(r8), intent(in) :: zice(LBi:,LBj:)
194# endif
195 real(r8), intent(in) :: t(LBi:,LBj:,:,:,:)
196 real(r8), intent(in) :: tl_t(LBi:,LBj:,:,:,:)
197# if !defined BBL_MODEL_NOT_YET || defined ICESHELF
198 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
199 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
200 real(r8), intent(in) :: tl_u(LBi:,LBj:,:,:)
201 real(r8), intent(in) :: tl_v(LBi:,LBj:,:,:)
202# endif
203 real(r8), intent(in) :: stflx(LBi:,LBj:,:)
204 real(r8), intent(in) :: btflx(LBi:,LBj:,:)
205# ifdef QCORRECTION
206 real(r8), intent(in) :: dqdt(LBi:,LBj:)
207 real(r8), intent(in) :: sst(LBi:,LBj:)
208# endif
209# if defined SCORRECTION || defined SRELAXATION
210 real(r8), intent(in) :: sss(LBi:,LBj:)
211# endif
212 real(r8), intent(in) :: stflux(LBi:,LBj:,:)
213 real(r8), intent(in) :: btflux(LBi:,LBj:,:)
214# if defined ICESHELF
215# ifdef SHORTWAVE
216 real(r8), intent(inout) :: srflx(LBi:,LBj:)
217# endif
218 real(r8), intent(inout) :: tl_sustr(LBi:,LBj:)
219 real(r8), intent(inout) :: tl_svstr(LBi:,LBj:)
220# endif
221# ifndef BBL_MODEL_NOT_YET
222 real(r8), intent(inout) :: tl_bustr(LBi:,LBj:)
223 real(r8), intent(inout) :: tl_bvstr(LBi:,LBj:)
224# endif
225 real(r8), intent(inout) :: tl_stflx(LBi:,LBj:,:)
226 real(r8), intent(inout) :: tl_btflx(LBi:,LBj:,:)
227# else
228 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
229 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
230# if defined UV_LOGDRAG
231 real(r8), intent(in) :: ZoBot(LBi:UBi,LBj:UBj)
232# elif defined UV_LDRAG
233 real(r8), intent(in) :: rdrag(LBi:UBi,LBj:UBj)
234# elif defined UV_QDRAG
235 real(r8), intent(in) :: rdrag2(LBi:UBi,LBj:UBj)
236# endif
237# if !defined BBL_MODEL_NOT_YET || defined ICESHELF
238 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
239 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
240 real(r8), intent(in) :: tl_z_r(LBi:UBi,LBj:UBj,N(ng))
241 real(r8), intent(in) :: tl_z_w(LBi:UBi,LBj:UBj,0:N(ng))
242# endif
243# if defined ICESHELF
244 real(r8), intent(in) :: zice(LBi:UBi,LBj:UBj)
245# endif
246 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
247 real(r8), intent(in) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
248# if !defined BBL_MODEL_NOT_YET || defined ICESHELF
249 real(r8), intent(in) :: u(LBi:UBi,LBj:UBj,N(ng),2)
250 real(r8), intent(in) :: v(LBi:UBi,LBj:UBj,N(ng),2)
251 real(r8), intent(in) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
252 real(r8), intent(in) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
253# endif
254 real(r8), intent(in) :: stflx(LBi:UBi,LBj:UBj,NT(ng))
255 real(r8), intent(in) :: btflx(LBi:UBi,LBj:UBj,NT(ng))
256# ifdef QCORRECTION
257 real(r8), intent(in) :: dqdt(LBi:UBi,LBj:UBj)
258 real(r8), intent(in) :: sst(LBi:UBi,LBj:UBj)
259# endif
260# if defined SCORRECTION || defined SRELAXATION
261 real(r8), intent(in) :: sss(LBi:UBi,LBj:UBj)
262# endif
263 real(r8), intent(in) :: stflux(LBi:UBi,LBj:UBj,NT(ng))
264 real(r8), intent(in) :: btflux(LBi:UBi,LBj:UBj,NT(ng))
265# if defined ICESHELF
266# ifdef SHORTWAVE
267 real(r8), intent(inout) :: srflx(LBi:UBi,LBj:UBj)
268# endif
269 real(r8), intent(inout) :: tl_sustr(LBi:UBi,LBj:UBj)
270 real(r8), intent(inout) :: tl_svstr(LBi:UBi,LBj:UBj)
271# endif
272# ifndef BBL_MODEL_NOT_YET
273 real(r8), intent(inout) :: tl_bustr(LBi:UBi,LBj:UBj)
274 real(r8), intent(inout) :: tl_bvstr(LBi:UBi,LBj:UBj)
275# endif
276 real(r8), intent(inout) :: tl_stflx(LBi:UBi,LBj:UBj,NT(ng))
277 real(r8), intent(inout) :: tl_btflx(LBi:UBi,LBj:UBj,NT(ng))
278# endif
279!
280! Local variable declarations.
281!
282 integer :: i, j, itrc
283!
284 real(r8) :: EmP, tl_EmP
285# if !defined BBL_MODEL_NOT_YET || defined ICESHELF
286 real(r8) :: cff1, cff2, cff3
287 real(r8) :: tl_cff1, tl_cff2, tl_cff3
288# endif
289
290# if (!defined BBL_MODEL_NOT_YET || defined ICESHELF) && defined UV_LOGDRAG
291 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: wrk
292 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_wrk
293# endif
294
295# include "set_bounds.h"
296!
297!-----------------------------------------------------------------------
298! Tangent linear of loading surface and bottom net heat flux (degC m/s)
299! into state variables "stflx" and "btflx".
300!
301! Notice that the forcing net heat flux stflux(:,:,itemp) is processed
302! elsewhere from data, coupling, bulk flux parameterization,
303! or analytical formulas.
304!
305! During model coupling, we need to make sure that this forcing is
306! unaltered during the coupling interval when ROMS timestep size is
307! smaller. Notice that further manipulations to state variable
308! stflx(:,:,itemp) are allowed below.
309!-----------------------------------------------------------------------
310!
311 DO j=jstrr,jendr
312 DO i=istrr,iendr
313!^ stflx(i,j,itemp)=stflux(i,j,itemp)
314!^
315 tl_stflx(i,j,itemp)=0.0_r8 ! not needed in TLM
316!^ btflx(i,j,itemp)=btflux(i,j,itemp)
317!^
318 tl_btflx(i,j,itemp)=0.0_r8 ! not needed in TLM
319 END DO
320 END DO
321
322# ifdef QCORRECTION
323!
324!-----------------------------------------------------------------------
325! Add in flux correction to surface net heat flux (degC m/s).
326!-----------------------------------------------------------------------
327!
328! Add in net heat flux correction.
329!
330 DO j=jstrr,jendr
331 DO i=istrr,iendr
332!^ stflx(i,j,itemp)=stflx(i,j,itemp)+ &
333!^ & dqdt(i,j)*(t(i,j,N(ng),nrhs,itemp)-sst(i,j))
334!^
335 tl_stflx(i,j,itemp)=tl_stflx(i,j,itemp)+ &
336 & dqdt(i,j)*tl_t(i,j,n(ng),nrhs,itemp)
337 END DO
338 END DO
339# endif
340
341# ifdef LIMIT_STFLX_COOLING
342!
343!-----------------------------------------------------------------------
344! If net heat flux is cooling and SST is at freezing point or below
345! then suppress further cooling. Note: stflx sign convention is that
346! positive means heating the ocean (J Wilkin).
347!-----------------------------------------------------------------------
348!
349! Below the surface heat flux stflx(:,:,itemp) is ZERO if cooling AND
350! the SST is cooler that the threshold. The value is retained if
351! warming.
352!
353! cff3 = 0 if SST warmer than threshold (cff1) - change nothing
354! cff3 = 1 if SST colder than threshold (cff1)
355!
356! 0.5*(cff2-ABS(cff2)) = 0 if flux is warming
357! = stflx(:,:,itemp) if flux is cooling
358!
359 cff1=-2.0_r8 ! nominal SST threshold to cease cooling
360 DO j=jstrr,jendr
361 DO i=istrr,iendr
362 cff2=stflx(i,j,itemp)
363 tl_cff2=tl_stflx(i,j,itemp)
364 cff3=0.5_r8*(1.0_r8+sign(1.0_r8,cff1-t(i,j,n(ng),nrhs,itemp)))
365!^ tl_cff3=0.5_r8*SIGN(1.0_r8, cff1-t(i,j,N(ng),nrhs,itemp))*0.0
366!^ tl_cff3=0.0_r8
367!^
368!^ stflx(i,j,itemp)=cff2-cff3*0.5_r8*(cff2-ABS(cff2))
369!^
370 tl_stflx(i,j,itemp)=(1.0_r8- &
371 & cff3*0.5_r8*(1.0_r8-sign(1.0_r8,cff2)))* &
372 & tl_cff2
373 END DO
374 END DO
375# endif
376
377# ifdef SALINITY
378!
379!-----------------------------------------------------------------------
380! Tangent linear of multiply freshwater fluxes with surface and bottom
381! salinity.
382!
383! If appropriate, apply correction. Notice that input stflux(:,:,isalt)
384! is the net freshwater flux (E-P; m/s) from data, coupling, bulk flux
385! parameterization, or analytical formula. It has not been multiplied
386! by the surface and bottom salinity.
387!-----------------------------------------------------------------------
388!
389 DO j=jstrr,jendr
390 DO i=istrr,iendr
391 emp=stflux(i,j,isalt)
392 tl_emp=0.0_r8
393# if defined SCORRECTION
394!^ stflx(i,j,isalt)=EmP*t(i,j,N(ng),nrhs,isalt)- &
395!^ & Tnudg(isalt,ng)*Hz(i,j,N(ng))* &
396!^ & (t(i,j,N(ng),nrhs,isalt)-sss(i,j))
397!^
398 tl_stflx(i,j,isalt)=emp*tl_t(i,j,n(ng),nrhs,isalt)+ &
399 & tl_emp*t(i,j,n(ng),nrhs,isalt)- &
400 & tnudg(isalt,ng)* &
401 & (tl_hz(i,j,n(ng))* &
402 & (t(i,j,n(ng),nrhs,isalt)-sss(i,j))+ &
403 & hz(i,j,n(ng))* &
404 & tl_t(i,j,n(ng),nrhs,isalt))
405# elif defined SRELAXATION
406!^ stflx(i,j,isalt)=-Tnudg(isalt,ng)*Hz(i,j,N(ng))* &
407!^ & (t(i,j,N(ng),nrhs,isalt)-sss(i,j))
408!^
409 tl_stflx(i,j,isalt)=-tnudg(isalt,ng)* &
410 & (tl_hz(i,j,n(ng))* &
411 & (t(i,j,n(ng),nrhs,isalt)-sss(i,j))+ &
412 & hz(i,j,n(ng))* &
413 & tl_t(i,j,n(ng),nrhs,isalt))
414# else
415!^ stflx(i,j,isalt)=EmP*t(i,j,N(ng),nrhs,isalt)
416!^
417 tl_stflx(i,j,isalt)=emp*tl_t(i,j,n(ng),nrhs,isalt)+ &
418 & tl_emp*t(i,j,n(ng),nrhs,isalt)
419# endif
420!^ btflx(i,j,isalt)=btflx(i,j,isalt)*t(i,j,1,nrhs,isalt)
421!^
422 tl_btflx(i,j,isalt)=btflx(i,j,isalt)* &
423 & tl_t(i,j,1,nrhs,isalt)
424 END DO
425 END DO
426# endif
427
428# if defined BIOLOGY || defined SEDIMENT || defined T_PASSIVE
429!
430!-----------------------------------------------------------------------
431! Load surface and bottom passive tracer fluxes (T m/s).
432!-----------------------------------------------------------------------
433!
434 DO itrc=nat+1,nt(ng)
435 DO j=jstrr,jendr
436 DO i=istrr,iendr
437!^ stflx(i,j,itrc)=stflux(i,j,itrc)
438!^
439 tl_stflx(i,j,itrc)=0.0_r8
440!^ btflx(i,j,itrc)=btflux(i,j,itrc)
441!^
442 tl_btflx(i,j,itrc)=0.0_r8
443 END DO
444 END DO
445 END DO
446# endif
447
448# ifdef ICESHELF
449!
450!-----------------------------------------------------------------------
451! If ice shelf cavities, zero out for now the surface tracer flux
452! over the ice.
453!-----------------------------------------------------------------------
454!
455 DO itrc=1,nt(ng)
456 DO j=jstrr,jendr
457 DO i=istrr,iendr
458 IF (zice(i,j).ne.0.0_r8) THEN
459!^ stflx(i,j,itrc)=0.0_r8
460!^
461 tl_stflx(i,j,itrc)=0.0_r8
462 END IF
463 END DO
464 END DO
465 END DO
466# ifdef SHORTWAVE
467 DO j=jstrr,jendr
468 DO i=istrr,iendr
469 IF (zice(i,j).ne.0.0_r8) THEN
470!^ srflx(i,j)=0.0_r8
471!^
472 srflx(i,j)=0.0_r8
473 END IF
474 END DO
475 END DO
476# endif
477!
478!-----------------------------------------------------------------------
479! If ice shelf cavities, replace surface wind stress with ice shelf
480! cavity stress (m2/s2).
481!-----------------------------------------------------------------------
482
483# if defined UV_LOGDRAG
484!
485! Set logarithmic ice shelf cavity stress.
486!
487 DO j=jstrv-1,jend
488 DO i=istru-1,iend
489 cff1=1.0_r8/log((z_w(i,j,n(ng))-z_r(i,j,n(ng)))/zobot(i,j))
490 tl_cff1=-cff1*cff1*(tl_z_w(i,j,n(ng))-tl_z_r(i,j,n(ng)))/ &
491 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))
492 cff2=vonkar*vonkar*cff1*cff1
493 tl_cff2=vonkar*vonkar*2.0_r8*cff1*tl_cff1
494 cff3=max(cdb_min,cff2)
495 tl_cff3=(0.5_r8-sign(0.5_r8,cdb_min-cff2))*tl_cff2
496 wrk(i,j)=min(cdb_max,cff3)
497 tl_wrk(i,j)=(0.5_r8-sign(0.5_r8,cff3-cdb_max))*tl_cff3
498 END DO
499 END DO
500 DO j=jstr,jend
501 DO i=istru,iend
502 IF (zice(i,j)*zice(i-1,j).ne.0.0_r8) THEN
503 cff1=0.25_r8*(v(i ,j ,n(ng),nrhs)+ &
504 & v(i ,j+1,n(ng),nrhs)+ &
505 & v(i-1,j ,n(ng),nrhs)+ &
506 & v(i-1,j+1,n(ng),nrhs))
507 tl_cff1=0.25_r8*(tl_v(i ,j ,n(ng),nrhs)+ &
508 & tl_v(i ,j+1,n(ng),nrhs)+ &
509 & tl_v(i-1,j ,n(ng),nrhs)+ &
510 & tl_v(i-1,j+1,n(ng),nrhs))
511 cff2=sqrt(u(i,j,n(ng),nrhs)*u(i,j,n(ng),nrhs)+cff1*cff1)
512 IF (cff2.ne.0.0_r8) THEN
513 tl_cff2=(u(i,j,n(ng),nrhs)*tl_u(i,j,n(ng),nrhs)+ &
514 & cff1*tl_cff1)/cff2
515 ELSE
516 tl_cff2=0.0_r8
517 END IF
518!^ sustr(i,j)=-0.5_r8*(wrk(i-1,j)+wrk(i,j))* &
519!^ & u(i,j,N(ng),nrhs)*cff2
520!^
521 tl_sustr(i,j)=-0.5_r8* &
522 & ((tl_wrk(i-1,j)+tl_wrk(i,j))* &
523 & u(i,j,n(ng),nrhs)*cff2+ &
524 & (wrk(i-1,j)+wrk(i,j))* &
525 & (tl_u(i,j,n(ng),nrhs)*cff2+ &
526 & u(i,j,n(ng),nrhs)*tl_cff2))
527 END IF
528 END DO
529 END DO
530 DO j=jstrv,jend
531 DO i=istr,iend
532 IF (zice(i,j)*zice(i,j-1).ne.0.0_r8) THEN
533 cff1=0.25_r8*(u(i ,j ,n(ng),nrhs)+ &
534 & u(i+1,j ,n(ng),nrhs)+ &
535 & u(i ,j-1,n(ng),nrhs)+ &
536 & u(i+1,j-1,n(ng),nrhs))
537 tl_cff1=0.25_r8*(tl_u(i ,j ,n(ng),nrhs)+ &
538 & tl_u(i+1,j ,n(ng),nrhs)+ &
539 & tl_u(i ,j-1,n(ng),nrhs)+ &
540 & tl_u(i+1,j-1,n(ng),nrhs))
541 cff2=sqrt(cff1*cff1+v(i,j,n(ng),nrhs)*v(i,j,n(ng),nrhs))
542 IF (cff2.ne.0.0_r8) THEN
543 tl_cff2=(cff1*tl_cff1+ &
544 & v(i,j,n(ng),nrhs)*tl_v(i,j,n(ng),nrhs))/cff2
545 ELSE
546 tl_cff2=0.0_r8
547 END IF
548!^ svstr(i,j)=-0.5_r8*(wrk(i,j-1)+wrk(i,j))* &
549!^ & v(i,j,N(ng),nrhs)*cff2
550!^
551 tl_svstr(i,j)=-0.5_r8* &
552 & ((tl_wrk(i,j-1)+tl_wrk(i,j))* &
553 & v(i,j,n(ng),nrhs)*cff2+ &
554 & (wrk(i,j-1)+wrk(i,j))* &
555 & (tl_v(i,j,n(ng),nrhs)*cff2+ &
556 & v(i,j,n(ng),nrhs)*tl_cff2))
557 END IF
558 END DO
559 END DO
560# elif defined UV_QDRAG
561!
562! Set quadratic ice shelf cavity stress.
563!
564 DO j=jstr,jend
565 DO i=istru,iend
566 IF (zice(i,j)*zice(i-1,j).ne.0.0_r8) THEN
567 cff1=0.25_r8*(v(i ,j ,n(ng),nrhs)+ &
568 & v(i ,j+1,n(ng),nrhs)+ &
569 & v(i-1,j ,n(ng),nrhs)+ &
570 & v(i-1,j+1,n(ng),nrhs))
571 tl_cff1=0.25_r8*(tl_v(i ,j ,n(ng),nrhs)+ &
572 & tl_v(i ,j+1,n(ng),nrhs)+ &
573 & tl_v(i-1,j ,n(ng),nrhs)+ &
574 & tl_v(i-1,j+1,n(ng),nrhs))
575 & cff2=sqrt(u(i,j,n(ng),nrhs)*u(i,j,n(ng),nrhs)+cff1*cff1)
576 IF (cff2.ne.0.0_r8) THEN
577 tl_cff2=(u(i,j,n(ng),nrhs)*tl_u(i,j,n(ng),nrhs)+ &
578 & cff1*tl_cff1)/cff2
579 ELSE
580 tl_cff2=0.0_r8
581 END IF
582!^ sustr(i,j)=-0.5_r8*(rdrag2(i-1,j)+rdrag2(i,j))* &
583!^ & u(i,j,N(ng),nrhs)*cff2
584!^
585 tl_sustr(i,j)=-0.5_r8*(rdrag2(i-1,j)+rdrag2(i,j))* &
586 & (tl_u(i,j,n(ng),nrhs)*cff2+ &
587 & u(i,j,n(ng),nrhs)*tl_cff2)
588 END IF
589 END DO
590 END DO
591 DO j=jstrv,jend
592 DO i=istr,iend
593 IF (zice(i,j)*zice(i,j-1).ne.0.0_r8) THEN
594 cff1=0.25_r8*(u(i ,j ,n(ng),nrhs)+ &
595 & u(i+1,j ,n(ng),nrhs)+ &
596 & u(i ,j-1,n(ng),nrhs)+ &
597 & u(i+1,j-1,n(ng),nrhs))
598 tl_cff1=0.25_r8*(tl_u(i ,j ,n(ng),nrhs)+ &
599 & tl_u(i+1,j ,n(ng),nrhs)+ &
600 & tl_u(i ,j-1,n(ng),nrhs)+ &
601 & tl_u(i+1,j-1,n(ng),nrhs))
602 cff2=sqrt(cff1*cff1+v(i,j,n(ng),nrhs)*v(i,j,n(ng),nrhs))
603 IF (cff2.ne.0.0_r8) THEN
604 tl_cff2=(cff1*tl_cff1+ &
605 & v(i,j,n(ng),nrhs)*tl_v(i,j,n(ng),nrhs))/cff2
606 ELSE
607 tl_cff2=0.0_r8
608 END IF
609!^ svstr(i,j)=-0.5_r8*(rdrag2(i,j-1)+rdrag2(i,j))* &
610!^ & v(i,j,N(ng),nrhs)*cff2
611!^
612 tl_svstr(i,j)=-0.5_r8*(rdrag2(i,j-1)+rdrag2(i,j))* &
613 & (tl_v(i,j,n(ng),nrhs)*cff2+ &
614 & v(i,j,n(ng),nrhs)*tl_cff2)
615 END IF
616 END DO
617 END DO
618# elif defined UV_LDRAG
619!
620! Set linear ice shelf cavity stress.
621!
622 DO j=jstr,jend
623 DO i=istru,iend
624 IF (zice(i,j)*zice(i-1,j).ne.0.0_r8) THEN
625!^ sustr(i,j)=-0.5_r8*(rdrag(i-1,j)+rdrag(i,j))* &
626!^ & u(i,j,N(ng),nrhs)
627!^
628 tl_sustr(i,j)=-0.5_r8*(rdrag(i-1,j)+rdrag(i,j))* &
629 & tl_u(i,j,n(ng),nrhs)
630 END IF
631 END DO
632 END DO
633 DO j=jstrv,jend
634 DO i=istr,iend
635 IF (zice(i,j)*zice(i,j-1).ne.0.0_r8) THEN
636!^ svstr(i,j)=-0.5_r8*(rdrag(i,j-1)+rdrag(i,j))* &
637!^ & v(i,j,N(ng),nrhs)
638!^
639 tl_svstr(i,j)=-0.5_r8*(rdrag(i,j-1)+rdrag(i,j))* &
640 & tl_v(i,j,n(ng),nrhs)
641 END IF
642 END DO
643 END DO
644# else
645 DO j=jstr,jend
646 DO i=istru,iend
647 IF (zice(i,j)*zice(i-1,j).ne.0.0_r8) THEN
648!^ sustr(i,j)=0.0_r8
649!^
650 tl_sustr(i,j)=0.0_r8
651 END IF
652 END DO
653 END DO
654 DO j=jstrv,jend
655 DO i=istr,iend
656 IF (zice(i,j)*zice(i,j-1).ne.0.0_r8) THEN
657!^ svstr(i,j)=0.0_r8
658!^
659 tl_svstr(i,j)=0.0_r8
660 END IF
661 END DO
662 END DO
663# endif
664!
665! Apply periodic or gradient boundary conditions for output
666! purposes only.
667!
668!^ CALL bc_u2d_tile (ng, tile, &
669!^ & LBi, UBi, LBj, UBj, &
670!^ & sustr)
671!^
672 CALL bc_u2d_tile (ng, tile, &
673 & lbi, ubi, lbj, ubj, &
674 & tl_sustr)
675!^ CALL bc_v2d_tile (ng, tile, &
676!^ & LBi, UBi, LBj, UBj, &
677!^ & svstr)
678!^
679 CALL bc_v2d_tile (ng, tile, &
680 & lbi, ubi, lbj, ubj, &
681 & tl_svstr)
682
683# ifdef DISTRIBUTE
684!^ CALL mp_exchange2d (ng, tile, iNLM, 2, &
685!^ & LBi, UBi, LBj, UBj, &
686!^ & NghostPoints, &
687!^ & EWperiodic(ng), NSperiodic(ng), &
688!^ & sustr, svstr)
689!^
690 CALL mp_exchange2d (ng, tile, itlm, 2, &
691 & lbi, ubi, lbj, ubj, &
692 & nghostpoints, &
693 & ewperiodic(ng), nsperiodic(ng), &
694 & tl_sustr, tl_svstr)
695# endif
696# endif
697# ifndef BBL_MODEL_NOT_YET
698!
699!-----------------------------------------------------------------------
700! Set kinematic bottom momentum flux (m2/s2).
701!-----------------------------------------------------------------------
702
703# if defined UV_LOGDRAG
704!
705! Set logarithmic bottom stress.
706!
707 DO j=jstrv-1,jend
708 DO i=istru-1,iend
709 cff1=1.0_r8/log((z_r(i,j,1)-z_w(i,j,0))/zobot(i,j))
710 tl_cff1=-cff1*cff1*(tl_z_r(i,j,1)-tl_z_w(i,j,0))/ &
711 & (z_r(i,j,1)-z_w(i,j,0))
712 cff2=vonkar*vonkar*cff1*cff1
713 tl_cff2=vonkar*vonkar*2.0_r8*cff1*tl_cff1
714 cff3=max(cdb_min,cff2)
715 tl_cff3=(0.5_r8-sign(0.5_r8,cdb_min-cff2))*tl_cff2
716 wrk(i,j)=min(cdb_max,cff3)
717 tl_wrk(i,j)=(0.5_r8-sign(0.5_r8,cff3-cdb_max))*tl_cff3
718 END DO
719 END DO
720 DO j=jstr,jend
721 DO i=istru,iend
722 cff1=0.25_r8*(v(i ,j ,1,nrhs)+ &
723 & v(i ,j+1,1,nrhs)+ &
724 & v(i-1,j ,1,nrhs)+ &
725 & v(i-1,j+1,1,nrhs))
726 tl_cff1=0.25_r8*(tl_v(i ,j ,1,nrhs)+ &
727 & tl_v(i ,j+1,1,nrhs)+ &
728 & tl_v(i-1,j ,1,nrhs)+ &
729 & tl_v(i-1,j+1,1,nrhs))
730 cff2=sqrt(u(i,j,1,nrhs)*u(i,j,1,nrhs)+cff1*cff1)
731 IF (cff2.ne.0.0_r8) THEN
732 tl_cff2=(u(i,j,1,nrhs)*tl_u(i,j,1,nrhs)+cff1*tl_cff1)/cff2
733 ELSE
734 tl_cff2=0.0_r8
735 END IF
736!^ bustr(i,j)=0.5_r8*(wrk(i-1,j)+wrk(i,j))* &
737!^ & u(i,j,1,nrhs)*cff2
738!^
739 tl_bustr(i,j)=0.5_r8* &
740 & ((tl_wrk(i-1,j)+tl_wrk(i,j))* &
741 & u(i,j,1,nrhs)*cff2+ &
742 & (wrk(i-1,j)+wrk(i,j))* &
743 & (tl_u(i,j,1,nrhs)*cff2+ &
744 & u(i,j,1,nrhs)*tl_cff2))
745 END DO
746 END DO
747 DO j=jstrv,jend
748 DO i=istr,iend
749 cff1=0.25_r8*(u(i ,j ,1,nrhs)+ &
750 & u(i+1,j ,1,nrhs)+ &
751 & u(i ,j-1,1,nrhs)+ &
752 & u(i+1,j-1,1,nrhs))
753 tl_cff1=0.25_r8*(tl_u(i ,j ,1,nrhs)+ &
754 & tl_u(i+1,j ,1,nrhs)+ &
755 & tl_u(i ,j-1,1,nrhs)+ &
756 & tl_u(i+1,j-1,1,nrhs))
757 cff2=sqrt(cff1*cff1+v(i,j,1,nrhs)*v(i,j,1,nrhs))
758 IF (cff2.ne.0.0_r8) THEN
759 tl_cff2=(cff1*tl_cff1+v(i,j,1,nrhs)*tl_v(i,j,1,nrhs))/cff2
760 ELSE
761 tl_cff2=0.0_r8
762 END IF
763!^ bvstr(i,j)=0.5_r8*(wrk(i,j-1)+wrk(i,j))* &
764!^ & v(i,j,1,nrhs)*cff2
765!^
766 tl_bvstr(i,j)=0.5_r8* &
767 & ((tl_wrk(i,j-1)+tl_wrk(i,j))* &
768 & v(i,j,1,nrhs)*cff2+ &
769 & (wrk(i,j-1)+wrk(i,j))* &
770 & (tl_v(i,j,1,nrhs)*cff2+ &
771 & v(i,j,1,nrhs)*tl_cff2))
772 END DO
773 END DO
774# elif defined UV_QDRAG
775!
776! Set quadratic bottom stress.
777!
778 DO j=jstr,jend
779 DO i=istru,iend
780 cff1=0.25_r8*(v(i ,j ,1,nrhs)+ &
781 & v(i ,j+1,1,nrhs)+ &
782 & v(i-1,j ,1,nrhs)+ &
783 & v(i-1,j+1,1,nrhs))
784 tl_cff1=0.25_r8*(tl_v(i ,j ,1,nrhs)+ &
785 & tl_v(i ,j+1,1,nrhs)+ &
786 & tl_v(i-1,j ,1,nrhs)+ &
787 & tl_v(i-1,j+1,1,nrhs))
788 cff2=sqrt(u(i,j,1,nrhs)*u(i,j,1,nrhs)+cff1*cff1)
789 IF (cff2.ne.0.0_r8) THEN
790 tl_cff2=(u(i,j,1,nrhs)*tl_u(i,j,1,nrhs)+cff1*tl_cff1)/cff2
791 ELSE
792 tl_cff2=0.0_r8
793 END IF
794!^ bustr(i,j)=0.5_r8*(rdrag2(i-1,j)+rdrag2(i,j))* &
795!^ & u(i,j,1,nrhs)*cff2
796!^
797 tl_bustr(i,j)=0.5_r8*(rdrag2(i-1,j)+rdrag2(i,j))* &
798 & (tl_u(i,j,1,nrhs)*cff2+ &
799 & u(i,j,1,nrhs)*tl_cff2)
800 END DO
801 END DO
802 DO j=jstrv,jend
803 DO i=istr,iend
804 cff1=0.25_r8*(u(i ,j ,1,nrhs)+ &
805 & u(i+1,j ,1,nrhs)+ &
806 & u(i ,j-1,1,nrhs)+ &
807 & u(i+1,j-1,1,nrhs))
808 tl_cff1=0.25_r8*(tl_u(i ,j ,1,nrhs)+ &
809 & tl_u(i+1,j ,1,nrhs)+ &
810 & tl_u(i ,j-1,1,nrhs)+ &
811 & tl_u(i+1,j-1,1,nrhs))
812 cff2=sqrt(cff1*cff1+v(i,j,1,nrhs)*v(i,j,1,nrhs))
813 IF (cff2.ne.0.0_r8) THEN
814 tl_cff2=(cff1*tl_cff1+v(i,j,1,nrhs)*tl_v(i,j,1,nrhs))/cff2
815 ELSE
816 tl_cff2=0.0_r8
817 END IF
818!^ bvstr(i,j)=0.5_r8*(rdrag2(i,j-1)+rdrag2(i,j))* &
819!^ & v(i,j,1,nrhs)*cff2
820!^
821 tl_bvstr(i,j)=0.5_r8*(rdrag2(i,j-1)+rdrag2(i,j))* &
822 & (tl_v(i,j,1,nrhs)*cff2+ &
823 & v(i,j,1,nrhs)*tl_cff2)
824 END DO
825 END DO
826# elif defined UV_LDRAG
827!
828! Set linear bottom stress.
829!
830 DO j=jstr,jend
831 DO i=istru,iend
832!^ bustr(i,j)=0.5_r8*(rdrag(i-1,j)+rdrag(i,j))* &
833!^ & u(i,j,1,nrhs)
834!^
835 tl_bustr(i,j)=0.5_r8*(rdrag(i-1,j)+rdrag(i,j))* &
836 & tl_u(i,j,1,nrhs)
837 END DO
838 END DO
839 DO j=jstrv,jend
840 DO i=istr,iend
841!^ bvstr(i,j)=0.5_r8*(rdrag(i,j-1)+rdrag(i,j))* &
842!^ & v(i,j,1,nrhs)
843!^
844 tl_bvstr(i,j)=0.5_r8*(rdrag(i,j-1)+rdrag(i,j))* &
845 & tl_v(i,j,1,nrhs)
846 END DO
847 END DO
848# endif
849!
850! Apply boundary conditions.
851!
852!^ CALL bc_u2d_tile (ng, tile, &
853!^ & LBi, UBi, LBj, UBj, &
854!^ & bustr)
855!^
856 CALL bc_u2d_tile (ng, tile, &
857 & lbi, ubi, lbj, ubj, &
858 & tl_bustr)
859!^ CALL bc_v2d_tile (ng, tile, &
860!^ & LBi, UBi, LBj, UBj, &
861!^ & bvstr)
862!^
863 CALL bc_v2d_tile (ng, tile, &
864 & lbi, ubi, lbj, ubj, &
865 & tl_bvstr)
866
867# ifdef DISTRIBUTE
868!^ CALL mp_exchange2d (ng, tile, iNLM, 2, &
869!^ & LBi, UBi, LBj, UBj, &
870!^ & NghostPoints, &
871!^ & EWperiodic(ng), NSperiodic(ng), &
872!^ & bustr, bvstr)
873!^
874 CALL mp_exchange2d (ng, tile, itlm, 2, &
875 & lbi, ubi, lbj, ubj, &
876 & nghostpoints, &
877 & ewperiodic(ng), nsperiodic(ng), &
878 & tl_bustr, tl_bvstr)
879# endif
880# endif
881!
882 RETURN
883 END SUBROUTINE tl_set_vbc_tile
884
885# else
886
887!
888!***********************************************************************
889 SUBROUTINE tl_set_vbc (ng, tile)
890!***********************************************************************
891!
892 USE mod_param
893 USE mod_forces
894 USE mod_grid
895 USE mod_ocean
896 USE mod_stepping
897!
898! Imported variable declarations.
899!
900 integer, intent(in) :: ng, tile
901!
902! Local variable declarations.
903!
904 character (len=*), parameter :: MyFile = &
905 & __FILE__
906!
907# include "tile.h"
908!
909# ifdef PROFILE
910 CALL wclock_on (ng, itlm, 6, __line__, myfile)
911# endif
912 CALL tl_set_vbc_tile (ng, tile, &
913 & lbi, ubi, lbj, ubj, &
914 & imins, imaxs, jmins, jmaxs, &
915 & krhs(ng), kstp(ng), knew(ng), &
916# if defined UV_LDRAG
917 & grid(ng) % rdrag, &
918# elif defined UV_QDRAG
919 & grid(ng) % rdrag2, &
920# endif
921 & ocean(ng) % ubar, &
922 & ocean(ng) % vbar, &
923 & ocean(ng) % tl_ubar, &
924 & ocean(ng) % tl_vbar, &
925 & forces(ng) % tl_bustr, &
926 & forces(ng) % tl_bvstr)
927# ifdef PROFILE
928 CALL wclock_off (ng, itlm, 6, __line__, myfile)
929# endif
930!
931 RETURN
932 END SUBROUTINE tl_set_vbc
933!
934!***********************************************************************
935 SUBROUTINE tl_set_vbc_tile (ng, tile, &
936 & LBi, UBi, LBj, UBj, &
937 & IminS, ImaxS, JminS, JmaxS, &
938 & krhs, kstp, knew, &
939# if defined UV_LDRAG
940 & rdrag, &
941# elif defined UV_QDRAG
942 & rdrag2, &
943# endif
944 & ubar, vbar, &
945 & tl_ubar, tl_vbar, &
946 & tl_bustr, tl_bvstr)
947!***********************************************************************
948!
949 USE mod_param
950 USE mod_scalars
951!
952 USE bc_2d_mod
953# ifdef DISTRIBUTE
955# endif
956!
957! Imported variable declarations.
958!
959 integer, intent(in) :: ng, tile
960 integer, intent(in) :: LBi, UBi, LBj, UBj
961 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
962 integer, intent(in) :: krhs, kstp, knew
963!
964# ifdef ASSUMED_SHAPE
965# if defined UV_LDRAG
966 real(r8), intent(in) :: rdrag(LBi:,LBj:)
967# elif defined UV_QDRAG
968 real(r8), intent(in) :: rdrag2(LBi:,LBj:)
969# endif
970 real(r8), intent(in) :: ubar(LBi:,LBj:,:)
971 real(r8), intent(in) :: vbar(LBi:,LBj:,:)
972 real(r8), intent(in) :: tl_ubar(LBi:,LBj:,:)
973 real(r8), intent(in) :: tl_vbar(LBi:,LBj:,:)
974 real(r8), intent(inout) :: tl_bustr(LBi:,LBj:)
975 real(r8), intent(inout) :: tl_bvstr(LBi:,LBj:)
976# else
977# if defined UV_LDRAG
978 real(r8), intent(in) :: rdrag(LBi:UBi,LBj:UBj)
979# elif defined UV_QDRAG
980 real(r8), intent(in) :: rdrag2(LBi:UBi,LBj:UBj)
981# endif
982 real(r8), intent(in) :: ubar(LBi:UBi,LBj:UBj,:)
983 real(r8), intent(in) :: vbar(LBi:UBi,LBj:UBj,:)
984 real(r8), intent(in) :: tl_ubar(LBi:UBi,LBj:UBj,:)
985 real(r8), intent(in) :: tl_vbar(LBi:UBi,LBj:UBj,:)
986 real(r8), intent(inout) :: tl_bustr(LBi:UBi,LBj:UBj)
987 real(r8), intent(inout) :: tl_bvstr(LBi:UBi,LBj:UBj)
988# endif
989!
990! Local variable declarations.
991!
992 integer :: i, j
993!
994 real(r8) :: cff1, cff2
995 real(r8) :: tl_cff1, tl_cff2
996
997# include "set_bounds.h"
998!
999!-----------------------------------------------------------------------
1000! Set kinematic barotropic bottom momentum stress (m2/s2).
1001!-----------------------------------------------------------------------
1002
1003# if defined UV_LDRAG
1004!
1005! Set linear bottom stress.
1006!
1007 DO j=jstr,jend
1008 DO i=istru,iend
1009!^ bustr(i,j)=0.5_r8*(rdrag(i-1,j)+rdrag(i,j))* &
1010!^ & ubar(i,j,krhs)
1011!^
1012 tl_bustr(i,j)=0.5_r8*(rdrag(i-1,j)+rdrag(i,j))* &
1013 & tl_ubar(i,j,krhs)
1014 END DO
1015 END DO
1016 DO j=jstrv,jend
1017 DO i=istr,iend
1018!^ bvstr(i,j)=0.5_r8*(rdrag(i,j-1)+rdrag(i,j))* &
1019!^ & vbar(i,j,krhs)
1020!^
1021 tl_bvstr(i,j)=0.5_r8*(rdrag(i,j-1)+rdrag(i,j))* &
1022 & tl_vbar(i,j,krhs)
1023 END DO
1024 END DO
1025# elif defined UV_QDRAG
1026!
1027! Set quadratic bottom stress.
1028!
1029 DO j=jstr,jend
1030 DO i=istru,iend
1031 cff1=0.25_r8*(vbar(i ,j ,krhs)+ &
1032 & vbar(i ,j+1,krhs)+ &
1033 & vbar(i-1,j ,krhs)+ &
1034 & vbar(i-1,j+1,krhs))
1035 tl_cff1=0.25_r8*(tl_vbar(i ,j ,krhs)+ &
1036 & tl_vbar(i ,j+1,krhs)+ &
1037 & tl_vbar(i-1,j ,krhs)+ &
1038 & tl_vbar(i-1,j+1,krhs))
1039 cff2=sqrt(ubar(i,j,krhs)*ubar(i,j,krhs)+cff1*cff1)
1040 IF (cff2.ne.0.0_r8) THEN
1041 tl_cff2=(ubar(i,j,krhs)*tl_ubar(i,j,krhs)+cff1*tl_cff1)/cff2
1042 ELSE
1043 tl_cff2=0.0_r8
1044 END IF
1045!^ bustr(i,j)=0.5_r8*(rdrag2(i-1,j)+rdrag2(i,j))* &
1046!^ & ubar(i,j,krhs)*cff2
1047!^
1048 tl_bustr(i,j)=0.5_r8*(rdrag2(i-1,j)+rdrag2(i,j))* &
1049 & (tl_ubar(i,j,krhs)*cff2+ &
1050 & ubar(i,j,krhs)*tl_cff2)
1051 END DO
1052 END DO
1053 DO j=jstrv,jend
1054 DO i=istr,iend
1055 cff1=0.25_r8*(ubar(i ,j ,krhs)+ &
1056 & ubar(i+1,j ,krhs)+ &
1057 & ubar(i ,j-1,krhs)+ &
1058 & ubar(i+1,j-1,krhs))
1059 tl_cff1=0.25_r8*(tl_ubar(i ,j ,krhs)+ &
1060 & tl_ubar(i+1,j ,krhs)+ &
1061 & tl_ubar(i ,j-1,krhs)+ &
1062 & tl_ubar(i+1,j-1,krhs))
1063 cff2=sqrt(cff1*cff1+vbar(i,j,krhs)*vbar(i,j,krhs))
1064 IF (cff2.ne.0.0_r8) THEN
1065 tl_cff2=(cff1*tl_cff1+vbar(i,j,krhs)*tl_vbar(i,j,krhs))/cff2
1066 ELSE
1067 tl_cff2=0.0_r8
1068 END IF
1069!^ bvstr(i,j)=0.5_r8*(rdrag2(i,j-1)+rdrag2(i,j))* &
1070!^ & vbar(i,j,krhs)*cff2
1071!^
1072 tl_bvstr(i,j)=0.5_r8*(rdrag2(i,j-1)+rdrag2(i,j))* &
1073 & (tl_vbar(i,j,krhs)*cff2+ &
1074 & vbar(i,j,krhs)*tl_cff2)
1075 END DO
1076 END DO
1077# endif
1078!
1079! Apply boundary conditions.
1080!
1081!^ CALL bc_u2d_tile (ng, tile, &
1082!^ & LBi, UBi, LBj, UBj, &
1083!^ & bustr)
1084!^
1085 CALL bc_u2d_tile (ng, tile, &
1086 & lbi, ubi, lbj, ubj, &
1087 & tl_bustr)
1088!^ CALL bc_v2d_tile (ng, tile, &
1089!^ & LBi, UBi, LBj, UBj, &
1090!^ & bvstr)
1091!^
1092 CALL bc_v2d_tile (ng, tile, &
1093 & lbi, ubi, lbj, ubj, &
1094 & tl_bvstr)
1095
1096# ifdef DISTRIBUTE
1097!^ CALL mp_exchange2d (ng, tile, iNLM, 2, &
1098!^ & LBi, UBi, LBj, UBj, &
1099!^ & NghostPoints, &
1100!^ & EWperiodic(ng), NSperiodic(ng), &
1101!^ & bustr, bvstr)
1102!^
1103 CALL mp_exchange2d (ng, tile, itlm, 2, &
1104 & lbi, ubi, lbj, ubj, &
1105 & nghostpoints, &
1106 & ewperiodic(ng), nsperiodic(ng), &
1107 & tl_bustr, tl_bvstr)
1108# endif
1109!
1110 RETURN
1111 END SUBROUTINE tl_set_vbc_tile
1112# endif
1113#endif
1114 END MODULE tl_set_vbc_mod
subroutine bc_v2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition bc_2d.F:345
subroutine bc_u2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition bc_2d.F:167
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer nat
Definition mod_param.F:499
integer nghostpoints
Definition mod_param.F:710
integer, parameter itlm
Definition mod_param.F:663
real(dp) cdb_min
real(dp) vonkar
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
real(dp), dimension(:,:), allocatable tnudg
integer isalt
integer itemp
real(dp) cdb_max
integer, dimension(:), allocatable kstp
integer, dimension(:), allocatable knew
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable krhs
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine tl_set_vbc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, hz, tl_hz, zobot, z_r, z_w, tl_z_r, tl_z_w, zice, t, tl_t, u, v, tl_u, tl_v, dqdt, sst, sss, tl_sustr, tl_svstr, stflux, btflux, stflx, btflx, tl_stflx, tl_btflx)
Definition tl_set_vbc.F:159
subroutine, public tl_set_vbc(ng, tile)
Definition tl_set_vbc.F:33
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