/[MITgcm]/MITgcm/pkg/obcs/obcs_sponge.F
ViewVC logotype

Contents of /MITgcm/pkg/obcs/obcs_sponge.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph

Revision 1.14 - (show annotations) (download)
Sat Sep 27 00:13:40 2014 UTC (9 years, 7 months ago) by dimitri
Branch: MAIN
CVS Tags: checkpoint66g, checkpoint66f, checkpoint66e, checkpoint66d, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint66o, checkpoint66n, checkpoint66m, checkpoint66l, checkpoint66k, checkpoint66j, checkpoint66i, checkpoint66h, checkpoint65z, checkpoint65x, checkpoint65y, checkpoint65r, checkpoint65s, checkpoint65p, checkpoint65q, checkpoint65v, checkpoint65w, checkpoint65t, checkpoint65u, checkpoint65j, checkpoint65k, checkpoint65h, checkpoint65i, checkpoint65n, checkpoint65o, checkpoint65l, checkpoint65m, checkpoint65f, checkpoint65g, checkpoint65e, HEAD
Changes since 1.13: +527 -607 lines 
replaced the just-added OBCS_OPTIONS.h CPP options with run-time variables
see pkg/obcs/OBCS_PARAMS.h for details)
1 C $Header: /u/gcmpack/MITgcm/pkg/obcs/obcs_sponge.F,v 1.13 2014/09/25 00:37:05 dimitri Exp $
2 C $Name: $
3
4 #include "OBCS_OPTIONS.h"
5
6 C-- File obcs_sponge.F:
7 C-- Contents:
8 C-- o OBCS_SPONGE_U
9 C-- o OBCS_SPONGE_V
10 C-- o OBCS_SPONGE_T
11 C-- o OBCS_SPONGE_S
12
13 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
14
15 CStartOfInterface
16 SUBROUTINE OBCS_SPONGE_U(
17 U gU_arr,
18 I iMin,iMax,jMin,jMax, k, bi, bj,
19 I myTime, myIter, myThid )
20 C *==========================================================*
21 C | S/R OBCS_SPONGE_U
22 C | o Contains problem specific forcing for zonal velocity.
23 C *==========================================================*
24 C | Adds a relaxation term to gU near Open-Boundaries
25 C *==========================================================*
26 IMPLICIT NONE
27
28 C == Global data ==
29 #include "SIZE.h"
30 #include "EEPARAMS.h"
31 #include "PARAMS.h"
32 #include "GRID.h"
33 #include "DYNVARS.h"
34 #include "OBCS_PARAMS.h"
35 #include "OBCS_GRID.h"
36 #include "OBCS_FIELDS.h"
37
38 C == Routine arguments ==
39 C gU_arr :: the tendency array
40 C iMin,iMax :: Working range of x-index for applying forcing.
41 C jMin,jMax :: Working range of y-index for applying forcing.
42 C k :: Current vertical level index
43 C bi,bj :: Current tile indices
44 _RL gU_arr(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
45 INTEGER iMin, iMax, jMin, jMax
46 INTEGER k, bi, bj
47 _RL myTime
48 INTEGER myIter
49 INTEGER myThid
50 CEndOfInterface
51
52 #if defined(ALLOW_OBCS) && defined(ALLOW_OBCS_SPONGE)
53 C == Local variables ==
54 C Loop counters
55 INTEGER i, j, isl, jsl
56 _RL urelax, lambda_obcs_u
57
58 IF ( useOBCSsponge .AND. spongeThickness.NE.0 ) THEN
59
60 C Northern Open Boundary
61 #ifdef ALLOW_OBCS_NORTH
62 IF ( OBCSsponge_N .AND. OBCSsponge_UatNS ) THEN
63 DO i=iMin,iMax
64 IF ( OB_Jn(i,bi,bj).NE.OB_indexNone ) THEN
65 DO jsl= 1,spongeThickness
66 j=OB_Jn(i,bi,bj)-jsl
67 IF ((j.ge.jmin).and.(j.le.jmax)) THEN
68 IF (useLinearSponge) THEN
69 urelax = OBNu(i,k,bi,bj)
70 ELSE
71 urelax=(
72 & float(spongeThickness-jsl)*OBNu(i,k,bi,bj)
73 & + float(jsl)*uVel(i,j,k,bi,bj) )
74 & / float(spongeThickness)
75 ENDIF
76 lambda_obcs_u = (
77 & float(spongeThickness-jsl)*Vrelaxobcsbound
78 & + float(jsl)*Vrelaxobcsinner)
79 & / float(spongeThickness)
80 IF (lambda_obcs_u.ne.0.) THEN
81 lambda_obcs_u = 1. _d 0 / lambda_obcs_u
82 ELSE
83 lambda_obcs_u = 0. _d 0
84 ENDIF
85 gU_arr(i,j) = gU_arr(i,j)
86 & - _maskW(i,j,k,bi,bj) * lambda_obcs_u
87 & * ( uVel(i,j,k,bi,bj) - urelax )
88 ENDIF
89 ENDDO
90 ENDIF
91 ENDDO
92 ENDIF
93 #endif /* ALLOW_OBCS_NORTH */
94
95 C Southern Open Boundary
96 #ifdef ALLOW_OBCS_SOUTH
97 IF ( OBCSsponge_S .AND. OBCSsponge_UatNS ) THEN
98 DO i=iMin,iMax
99 IF ( OB_Js(i,bi,bj).NE.OB_indexNone ) THEN
100 DO jsl= 1,spongeThickness
101 j=OB_Js(i,bi,bj)+jsl
102 IF ((j.ge.jmin).and.(j.le.jmax)) THEN
103 IF (useLinearSponge) THEN
104 urelax= OBSu(i,k,bi,bj)
105 ELSE
106 urelax=(
107 & float(spongeThickness-jsl)*OBSu(i,k,bi,bj)
108 & + float(jsl)*uVel(i,j,k,bi,bj) )
109 & / float(spongeThickness)
110 ENDIF
111 lambda_obcs_u = (
112 & float(spongeThickness-jsl)*Vrelaxobcsbound
113 & + float(jsl)*Vrelaxobcsinner)
114 & / float(spongeThickness)
115 IF (lambda_obcs_u.ne.0.) THEN
116 lambda_obcs_u = 1. _d 0 / lambda_obcs_u
117 ELSE
118 lambda_obcs_u = 0. _d 0
119 ENDIF
120 gU_arr(i,j) = gU_arr(i,j)
121 & - _maskW(i,j,k,bi,bj) * lambda_obcs_u
122 & * ( uVel(i,j,k,bi,bj) - urelax )
123 ENDIF
124 ENDDO
125 ENDIF
126 ENDDO
127 ENDIF
128 #endif /* ALLOW_OBCS_SOUTH */
129
130 C Eastern Open Boundary
131 #ifdef ALLOW_OBCS_EAST
132 IF ( OBCSsponge_E .AND. OBCSsponge_UatEW ) THEN
133 DO j=jMin,jMax
134 IF ( OB_Ie(j,bi,bj).NE.OB_indexNone ) THEN
135 DO isl= 1,spongeThickness
136 i=OB_Ie(j,bi,bj)-isl
137 IF ((i.ge.imin).and.(i.le.imax)) THEN
138 IF (useLinearSponge) THEN
139 urelax=OBEu(j,k,bi,bj)
140 ELSE
141 urelax=(
142 & float(spongeThickness-isl)*OBEu(j,k,bi,bj)
143 & + float(isl)*uVel(i,j,k,bi,bj) )
144 & / float(spongeThickness)
145 ENDIF
146 lambda_obcs_u = (
147 & float(spongeThickness-isl)*Urelaxobcsbound
148 & + float(isl)*Urelaxobcsinner)
149 & / float(spongeThickness)
150 IF (lambda_obcs_u.ne.0.) THEN
151 lambda_obcs_u = 1. _d 0 / lambda_obcs_u
152 ELSE
153 lambda_obcs_u = 0. _d 0
154 ENDIF
155 gU_arr(i,j) = gU_arr(i,j)
156 & - _maskW(i,j,k,bi,bj) * lambda_obcs_u
157 & * ( uVel(i,j,k,bi,bj) - urelax )
158 ENDIF
159 ENDDO
160 ENDIF
161 ENDDO
162 ENDIF
163 #endif /* ALLOW_OBCS_EAST */
164
165 C Western Open Boundary
166 #ifdef ALLOW_OBCS_WEST
167 IF ( OBCSsponge_W .AND. OBCSsponge_UatEW ) THEN
168 DO j=jMin,jMax
169 IF ( OB_Iw(j,bi,bj).NE.OB_indexNone ) THEN
170 DO isl= 1,spongeThickness
171 i=OB_Iw(j,bi,bj)+isl+1
172 IF ((i.ge.imin).and.(i.le.imax)) THEN
173 IF (useLinearSponge) THEN
174 urelax= OBWu(j,k,bi,bj)
175 ELSE
176 urelax=(
177 & float(spongeThickness-isl)*OBWu(j,k,bi,bj)
178 & + float(isl)*uVel(i,j,k,bi,bj) )
179 & / float(spongeThickness)
180 ENDIF
181 lambda_obcs_u= (
182 & float(spongeThickness-isl)*Urelaxobcsbound
183 & + float(isl)*Urelaxobcsinner)
184 & / float(spongeThickness)
185 IF (lambda_obcs_u.ne.0.) THEN
186 lambda_obcs_u = 1. _d 0 / lambda_obcs_u
187 ELSE
188 lambda_obcs_u = 0. _d 0
189 ENDIF
190 gU_arr(i,j) = gU_arr(i,j)
191 & - _maskW(i,j,k,bi,bj) * lambda_obcs_u
192 & * ( uVel(i,j,k,bi,bj) - urelax )
193 ENDIF
194 ENDDO
195 ENDIF
196 ENDDO
197 ENDIF
198 #endif /* ALLOW_OBCS_WEST */
199
200 ENDIF
201
202 #endif /* ALLOW_OBCS & ALLOW_OBCS_SPONGE */
203
204 RETURN
205 END
206
207 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
208
209 CStartOfInterface
210 SUBROUTINE OBCS_SPONGE_V(
211 U gV_arr,
212 I iMin,iMax,jMin,jMax, k, bi, bj,
213 I myTime, myIter, myThid )
214 C *==========================================================*
215 C | S/R OBCS_SPONGE_V
216 C | o Contains problem specific forcing for merid velocity.
217 C *==========================================================*
218 C | Adds a relaxation term to gV near Open-Boundaries
219 C *==========================================================*
220 IMPLICIT NONE
221
222 C == Global data ==
223 #include "SIZE.h"
224 #include "EEPARAMS.h"
225 #include "PARAMS.h"
226 #include "GRID.h"
227 #include "DYNVARS.h"
228 #include "OBCS_PARAMS.h"
229 #include "OBCS_GRID.h"
230 #include "OBCS_FIELDS.h"
231
232 C == Routine arguments ==
233 C gV_arr :: the tendency array
234 C iMin,iMax :: Working range of x-index for applying forcing.
235 C jMin,jMax :: Working range of y-index for applying forcing.
236 C k :: Current vertical level index
237 C bi,bj :: Current tile indices
238 _RL gV_arr(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
239 INTEGER iMin, iMax, jMin, jMax
240 INTEGER k, bi, bj
241 _RL myTime
242 INTEGER myIter
243 INTEGER myThid
244 CEndOfInterface
245
246 #if defined(ALLOW_OBCS) && defined(ALLOW_OBCS_SPONGE)
247 C == Local variables ==
248 C Loop counters
249 INTEGER i, j, isl, jsl
250 _RL vrelax,lambda_obcs_v
251
252 IF ( useOBCSsponge .AND. spongeThickness.NE.0 ) THEN
253
254 C Northern Open Boundary
255 #ifdef ALLOW_OBCS_NORTH
256 IF ( OBCSsponge_N .AND. OBCSsponge_VatNS ) THEN
257 DO i=iMin,iMax
258 IF ( OB_Jn(i,bi,bj).NE.OB_indexNone ) THEN
259 DO jsl= 1,spongeThickness
260 j=OB_Jn(i,bi,bj)-jsl
261 IF ((j.ge.jmin).and.(j.le.jmax)) THEN
262 IF (useLinearSponge) THEN
263 vrelax= OBNv(i,k,bi,bj)
264 ELSE
265 vrelax=(
266 & float(spongeThickness-jsl)*OBNv(i,k,bi,bj)
267 & + float(jsl)*vVel(i,j,k,bi,bj) )
268 & / float(spongeThickness)
269 ENDIF
270 lambda_obcs_v = (
271 & float(spongeThickness-jsl)*Vrelaxobcsbound
272 & + float(jsl)*Vrelaxobcsinner)
273 & / float(spongeThickness)
274 IF (lambda_obcs_v.ne.0.) THEN
275 lambda_obcs_v = 1. _d 0 / lambda_obcs_v
276 ELSE
277 lambda_obcs_v = 0. _d 0
278 ENDIF
279 gV_arr(i,j) = gV_arr(i,j)
280 & - _maskS(i,j,k,bi,bj) * lambda_obcs_v
281 & * ( vVel(i,j,k,bi,bj) - vrelax )
282 ENDIF
283 ENDDO
284 ENDIF
285 ENDDO
286 ENDIF
287 #endif /* ALLOW_OBCS_NORTH */
288
289 C Southern Open Boundary
290 #ifdef ALLOW_OBCS_SOUTH
291 IF ( OBCSsponge_S .AND. OBCSsponge_VatNS ) THEN
292 DO i=iMin,iMax
293 IF ( OB_Js(i,bi,bj).NE.OB_indexNone ) THEN
294 DO jsl= 1,spongeThickness
295 j=OB_Js(i,bi,bj)+jsl+1
296 IF ((j.ge.jmin).and.(j.le.jmax)) THEN
297 IF (useLinearSponge) THEN
298 vrelax= OBSv(i,k,bi,bj)
299 ELSE
300 vrelax=(
301 & float(spongeThickness-jsl)*OBSv(i,k,bi,bj)
302 & + float(jsl)*vVel(i,j,k,bi,bj) )
303 & / float(spongeThickness)
304 ENDIF
305 lambda_obcs_v = (
306 & float(spongeThickness-jsl)*Vrelaxobcsbound
307 & + float(jsl)*Vrelaxobcsinner)
308 & / float(spongeThickness)
309 IF (lambda_obcs_v.ne.0.) THEN
310 lambda_obcs_v = 1. _d 0 / lambda_obcs_v
311 ELSE
312 lambda_obcs_v = 0. _d 0
313 ENDIF
314 gV_arr(i,j) = gV_arr(i,j)
315 & - _maskS(i,j,k,bi,bj) * lambda_obcs_v
316 & * ( vVel(i,j,k,bi,bj) - vrelax )
317 ENDIF
318 ENDDO
319 ENDIF
320 ENDDO
321 ENDIF
322 #endif /* ALLOW_OBCS_SOUTH */
323
324 C Eastern Open Boundary
325 #ifdef ALLOW_OBCS_EAST
326 IF ( OBCSsponge_E .AND. OBCSsponge_VatEW ) THEN
327 DO j=jMin,jMax
328 IF ( OB_Ie(j,bi,bj).NE.OB_indexNone ) THEN
329 DO isl= 1,spongeThickness
330 i=OB_Ie(j,bi,bj)-isl
331 IF ((i.ge.imin).and.(i.le.imax)) THEN
332 IF (useLinearSponge) THEN
333 vrelax= OBEv(j,k,bi,bj)
334 ELSE
335 vrelax=(
336 & float(spongeThickness-isl)*OBEv(j,k,bi,bj)
337 & + float(isl)*vVel(i,j,k,bi,bj) )
338 & / float(spongeThickness)
339 ENDIF
340 lambda_obcs_v = (
341 & float(spongeThickness-isl)*Urelaxobcsbound
342 & + float(isl)*Urelaxobcsinner)
343 & / float(spongeThickness)
344 If (lambda_obcs_v.ne.0.) THEN
345 lambda_obcs_v = 1. _d 0 / lambda_obcs_v
346 ELSE
347 lambda_obcs_v = 0. _d 0
348 ENDIF
349 gV_arr(i,j) = gV_arr(i,j)
350 & - _maskS(i,j,k,bi,bj) * lambda_obcs_v
351 & * ( vVel(i,j,k,bi,bj) - vrelax )
352 ENDIF
353 ENDDO
354 ENDIF
355 ENDDO
356 ENDIF
357 #endif /* ALLOW_OBCS_EAST */
358
359 C Western Open Boundary
360 #ifdef ALLOW_OBCS_WEST
361 IF ( OBCSsponge_W .AND. OBCSsponge_VatEW ) THEN
362 DO j=jMin,jMax
363 IF ( OB_Iw(j,bi,bj).NE.OB_indexNone ) THEN
364 DO isl= 1,spongeThickness
365 i=OB_Iw(j,bi,bj)+isl
366 IF ((i.ge.imin).and.(i.le.imax)) THEN
367 IF (useLinearSponge) THEN
368 vrelax = OBWv(j,k,bi,bj)
369 ELSE
370 vrelax=(
371 & float(spongeThickness-isl)*OBWv(j,k,bi,bj)
372 & + float(isl)*vVel(i,j,k,bi,bj) )
373 & / float(spongeThickness)
374 ENDIF
375 lambda_obcs_v = (
376 & float(spongeThickness-isl)*Urelaxobcsbound
377 & + float(isl)*Urelaxobcsinner)
378 & / float(spongeThickness)
379 IF (lambda_obcs_v.ne.0.) THEN
380 lambda_obcs_v = 1. _d 0 / lambda_obcs_v
381 ELSE
382 lambda_obcs_v = 0. _d 0
383 ENDIF
384 gV_arr(i,j) = gV_arr(i,j)
385 & - _maskS(i,j,k,bi,bj) * lambda_obcs_v
386 & * ( vVel(i,j,k,bi,bj) - vrelax )
387 ENDIF
388 ENDDO
389 ENDIF
390 ENDDO
391 ENDIF
392 #endif /* ALLOW_OBCS_WEST */
393
394 ENDIF
395
396 #endif /* ALLOW_OBCS & ALLOW_OBCS_SPONGE */
397
398 RETURN
399 END
400
401 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
402
403 CStartOfInterface
404 SUBROUTINE OBCS_SPONGE_T(
405 U gT_arr,
406 I iMin,iMax,jMin,jMax, k, bi, bj,
407 I myTime, myIter, myThid )
408 C *==========================================================*
409 C | S/R OBCS_SPONGE_T
410 C | o Contains problem specific forcing for temperature.
411 C *==========================================================*
412 C | Adds a relaxation term to gT near Open-Boundaries
413 C *==========================================================*
414 IMPLICIT NONE
415
416 C == Global data ==
417 #include "SIZE.h"
418 #include "EEPARAMS.h"
419 #include "PARAMS.h"
420 #include "GRID.h"
421 #include "DYNVARS.h"
422 #include "OBCS_PARAMS.h"
423 #include "OBCS_GRID.h"
424 #include "OBCS_FIELDS.h"
425 #ifdef ALLOW_AUTODIFF_TAMC
426 # include "tamc.h"
427 # include "tamc_keys.h"
428 #endif
429
430 C == Routine arguments ==
431 C gT_arr :: the tendency array
432 C iMin,iMax :: Working range of x-index for applying forcing.
433 C jMin,jMax :: Working range of y-index for applying forcing.
434 C k :: Current vertical level index
435 C bi,bj :: Current tile indices
436 _RL gT_arr(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
437 INTEGER iMin, iMax, jMin, jMax
438 INTEGER k, bi, bj
439 _RL myTime
440 INTEGER myIter
441 INTEGER myThid
442 CEndOfInterface
443
444 #if defined(ALLOW_OBCS) && defined(ALLOW_OBCS_SPONGE)
445 C == Local variables ==
446 C Loop counters
447 INTEGER i, j, isl, jsl
448 _RL trelax, lambda_obcs_t
449
450 IF ( useOBCSsponge .AND. spongeThickness.NE.0 ) THEN
451
452 #ifdef ALLOW_AUTODIFF_TAMC
453 act1 = bi - myBxLo(myThid)
454 max1 = myBxHi(myThid) - myBxLo(myThid) + 1
455 act2 = bj - myByLo(myThid)
456 max2 = myByHi(myThid) - myByLo(myThid) + 1
457 act3 = myThid - 1
458 max3 = nTx*nTy
459 act4 = ikey_dynamics - 1
460 ikey = (act1 + 1) + act2*max1
461 & + act3*max1*max2
462 & + act4*max1*max2*max3
463 kkey = (ikey-1)*Nr + k
464 #endif /* ALLOW_AUTODIFF_TAMC */
465
466 C Northern Open Boundary
467 #ifdef ALLOW_OBCS_NORTH
468 IF ( OBCSsponge_N .AND. OBCSsponge_Theta ) THEN
469 #ifdef ALLOW_AUTODIFF_TAMC
470 CADJ STORE OBNt(:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
471 #endif
472 DO i=iMin,iMax
473 IF ( OB_Jn(i,bi,bj).NE.OB_indexNone ) THEN
474 DO jsl= 1,spongeThickness
475 j=OB_Jn(i,bi,bj)-jsl
476 IF ((j.ge.jmin).and.(j.le.jmax)) THEN
477 IF (OBNt(i,k,bi,bj).ne. 0.d0) then
478 IF (useLinearSponge) THEN
479 trelax = OBNt(i,k,bi,bj)
480 ELSE
481 trelax=(
482 & float(spongeThickness-jsl)*OBNt(i,k,bi,bj)
483 & + float(jsl)*theta(i,j,k,bi,bj) )
484 & / float(spongeThickness)
485 ENDIF
486 lambda_obcs_t = (
487 & float(spongeThickness-jsl)*Vrelaxobcsbound
488 & + float(jsl)*Vrelaxobcsinner)
489 & / float(spongeThickness)
490 IF (lambda_obcs_t.ne.0.) THEN
491 lambda_obcs_t = 1. _d 0 / lambda_obcs_t
492 ELSE
493 lambda_obcs_t = 0. _d 0
494 ENDIF
495 gT_arr(i,j) = gT_arr(i,j)
496 & - maskC(i,j,k,bi,bj) * lambda_obcs_t
497 & * ( theta(i,j,k,bi,bj) - trelax )
498 ENDIF
499 ENDIF
500 ENDDO
501 ENDIF
502 ENDDO
503 ENDIF
504 #endif /* ALLOW_OBCS_NORTH */
505
506 C Southern Open Boundary
507 #ifdef ALLOW_OBCS_SOUTH
508 IF ( OBCSsponge_S .AND. OBCSsponge_Theta ) THEN
509 #ifdef ALLOW_AUTODIFF_TAMC
510 CADJ STORE OBSt(:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
511 #endif
512 DO i=iMin,iMax
513 IF ( OB_Js(i,bi,bj).NE.OB_indexNone ) THEN
514 DO jsl= 1,spongeThickness
515 j=OB_Js(i,bi,bj)+jsl
516 IF ((j.ge.jmin).and.(j.le.jmax)) THEN
517 IF (OBSt(i,k,bi,bj).ne. 0.d0) then
518 IF (useLinearSponge) THEN
519 trelax= OBSt(i,k,bi,bj)
520 ELSE
521 trelax=(
522 & float(spongeThickness-jsl)*OBSt(i,k,bi,bj)
523 & + float(jsl)*theta(i,j,k,bi,bj) )
524 & / float(spongeThickness)
525 ENDIF
526 lambda_obcs_t = (
527 & float(spongeThickness-jsl)*Vrelaxobcsbound
528 & + float(jsl)*Vrelaxobcsinner)
529 & / float(spongeThickness)
530 IF (lambda_obcs_t.ne.0.) THEN
531 lambda_obcs_t = 1. _d 0 / lambda_obcs_t
532 ELSE
533 lambda_obcs_t = 0. _d 0
534 ENDIF
535 gT_arr(i,j) = gT_arr(i,j)
536 & - maskC(i,j,k,bi,bj) * lambda_obcs_t
537 & * ( theta(i,j,k,bi,bj) - trelax )
538 ENDIF
539 ENDIF
540 ENDDO
541 ENDIF
542 ENDDO
543 ENDIF
544 #endif /* ALLOW_OBCS_SOUTH */
545
546 C Eastern Open Boundary
547 #ifdef ALLOW_OBCS_EAST
548 IF ( OBCSsponge_E .AND. OBCSsponge_Theta ) THEN
549 #ifdef ALLOW_AUTODIFF_TAMC
550 CADJ STORE OBEt(:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
551 #endif
552 DO j=jMin,jMax
553 IF ( OB_Ie(j,bi,bj).NE.OB_indexNone ) THEN
554 DO isl= 1,spongeThickness
555 i=OB_Ie(j,bi,bj)-isl
556 IF ((i.ge.imin).and.(i.le.imax)) THEN
557 IF (OBEt(j,k,bi,bj).ne. 0.d0) then
558 IF (useLinearSponge) THEN
559 trelax = OBEt(j,k,bi,bj)
560 ELSE
561 trelax=(
562 & float(spongeThickness-isl)*OBEt(j,k,bi,bj)
563 & + float(isl)*theta(i,j,k,bi,bj) )
564 & / float(spongeThickness)
565 ENDIF
566 lambda_obcs_t = (
567 & float(spongeThickness-isl)*Urelaxobcsbound
568 & + float(isl)*Urelaxobcsinner)
569 & / float(spongeThickness)
570 IF (lambda_obcs_t.ne.0.) THEN
571 lambda_obcs_t = 1. _d 0 / lambda_obcs_t
572 ELSE
573 lambda_obcs_t = 0. _d 0
574 ENDIF
575 gT_arr(i,j) = gT_arr(i,j)
576 & - maskC(i,j,k,bi,bj) * lambda_obcs_t
577 & * ( theta(i,j,k,bi,bj) - trelax )
578 ENDIF
579 ENDIF
580 ENDDO
581 ENDIF
582 ENDDO
583 ENDIF
584 #endif /* ALLOW_OBCS_EAST */
585
586 C Western Open Boundary
587 #ifdef ALLOW_OBCS_WEST
588 IF ( OBCSsponge_W .AND. OBCSsponge_Theta ) THEN
589 #ifdef ALLOW_AUTODIFF_TAMC
590 CADJ STORE OBWt(:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
591 #endif
592 DO j=jMin,jMax
593 IF ( OB_Iw(j,bi,bj).NE.OB_indexNone ) THEN
594 DO isl= 1,spongeThickness
595 i=OB_Iw(j,bi,bj)+isl
596 IF ((i.ge.imin).and.(i.le.imax)) THEN
597 IF (OBWt(j,k,bi,bj).ne. 0.d0) THEN
598 IF (useLinearSponge) THEN
599 trelax= OBWt(j,k,bi,bj)
600 ELSE
601 trelax=(
602 & float(spongeThickness-isl)*OBWt(j,k,bi,bj)
603 & + float(isl)*theta(i,j,k,bi,bj) )
604 & / float(spongeThickness)
605 ENDIF
606 lambda_obcs_t= (
607 & float(spongeThickness-isl)*Urelaxobcsbound
608 & + float(isl)*Urelaxobcsinner)
609 & / float(spongeThickness)
610 IF (lambda_obcs_t .ne. 0.) THEN
611 lambda_obcs_t = 1. _d 0 / lambda_obcs_t
612 ELSE
613 lambda_obcs_t = 0. _d 0
614 ENDIF
615 gT_arr(i,j) = gT_arr(i,j)
616 & - maskC(i,j,k,bi,bj) * lambda_obcs_t
617 & * ( theta(i,j,k,bi,bj) - trelax )
618 ENDIF
619 ENDIF
620 ENDDO
621 ENDIF
622 ENDDO
623 ENDIF
624 #endif /* ALLOW_OBCS_WEST */
625
626 ENDIF
627
628 #endif /* ALLOW_OBCS & ALLOW_OBCS_SPONGE */
629
630 RETURN
631 END
632
633 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
634
635 CStartOfInterface
636 SUBROUTINE OBCS_SPONGE_S(
637 U gS_arr,
638 I iMin,iMax,jMin,jMax, k, bi, bj,
639 I myTime, myIter, myThid )
640 C *==========================================================*
641 C | S/R OBCS_SPONGE_S
642 C | o Contains problem specific forcing for salinity.
643 C *==========================================================*
644 C | Adds a relaxation term to gS near Open-Boundaries
645 C *==========================================================*
646 IMPLICIT NONE
647
648 C == Global data ==
649 #include "SIZE.h"
650 #include "EEPARAMS.h"
651 #include "PARAMS.h"
652 #include "GRID.h"
653 #include "DYNVARS.h"
654 #include "OBCS_PARAMS.h"
655 #include "OBCS_GRID.h"
656 #include "OBCS_FIELDS.h"
657 #ifdef ALLOW_AUTODIFF_TAMC
658 # include "tamc.h"
659 # include "tamc_keys.h"
660 #endif
661
662 C == Routine arguments ==
663 C gS_arr :: the tendency array
664 C iMin,iMax :: Working range of x-index for applying forcing.
665 C jMin,jMax :: Working range of y-index for applying forcing.
666 C k :: Current vertical level index
667 C bi,bj :: Current tile indices
668 _RL gS_arr(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
669 INTEGER iMin, iMax, jMin, jMax
670 INTEGER k, bi, bj
671 _RL myTime
672 INTEGER myIter
673 INTEGER myThid
674 CEndOfInterface
675
676 #if defined(ALLOW_OBCS) && defined(ALLOW_OBCS_SPONGE)
677 C == Local variables ==
678 C Loop counters
679 INTEGER i, j, isl, jsl
680 _RL srelax, lambda_obcs_s
681
682 IF ( useOBCSsponge .AND. spongeThickness.NE.0 ) THEN
683
684 #ifdef ALLOW_AUTODIFF_TAMC
685 act1 = bi - myBxLo(myThid)
686 max1 = myBxHi(myThid) - myBxLo(myThid) + 1
687 act2 = bj - myByLo(myThid)
688 max2 = myByHi(myThid) - myByLo(myThid) + 1
689 act3 = myThid - 1
690 max3 = nTx*nTy
691 act4 = ikey_dynamics - 1
692 ikey = (act1 + 1) + act2*max1
693 & + act3*max1*max2
694 & + act4*max1*max2*max3
695 kkey = (ikey-1)*Nr + k
696 #endif /* ALLOW_AUTODIFF_TAMC */
697
698 C Northern Open Boundary
699 #ifdef ALLOW_OBCS_NORTH
700 IF ( OBCSsponge_N .AND. OBCSsponge_Salt ) THEN
701 #ifdef ALLOW_AUTODIFF_TAMC
702 CADJ STORE OBNs(:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
703 #endif
704 DO i=iMin,iMax
705 IF ( OB_Jn(i,bi,bj).NE.OB_indexNone ) THEN
706 DO jsl= 1,spongeThickness
707 j=OB_Jn(i,bi,bj)-jsl
708 IF ((j.ge.jmin).and.(j.le.jmax)) THEN
709 IF (OBNs(i,k,bi,bj).ne. 0.d0) then
710 IF (useLinearSponge) THEN
711 srelax= OBNs(i,k,bi,bj)
712 ELSE
713 srelax=(
714 & float(spongeThickness-jsl)*OBNs(i,k,bi,bj)
715 & + float(jsl)*salt(i,j,k,bi,bj) )
716 & / float(spongeThickness)
717 ENDIF
718 lambda_obcs_s = (
719 & float(spongeThickness-jsl)*Vrelaxobcsbound
720 & + float(jsl)*Vrelaxobcsinner)
721 & / float(spongeThickness)
722 IF (lambda_obcs_s.ne.0.) THEN
723 lambda_obcs_s = 1. _d 0 / lambda_obcs_s
724 ELSE
725 lambda_obcs_s = 0. _d 0
726 ENDIF
727 gS_arr(i,j) = gS_arr(i,j)
728 & - maskC(i,j,k,bi,bj) * lambda_obcs_s
729 & * ( salt(i,j,k,bi,bj) - srelax )
730 ENDIF
731 ENDIF
732 ENDDO
733 ENDIF
734 ENDDO
735 ENDIF
736 #endif /* ALLOW_OBCS_NORTH */
737
738 C Southern Open Boundary
739 #ifdef ALLOW_OBCS_SOUTH
740 IF ( OBCSsponge_S .AND. OBCSsponge_Salt ) THEN
741 #ifdef ALLOW_AUTODIFF_TAMC
742 CADJ STORE OBSs(:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
743 #endif
744 DO i=iMin,iMax
745 IF ( OB_Js(i,bi,bj).NE.OB_indexNone ) THEN
746 DO jsl= 1,spongeThickness
747 j=OB_Js(i,bi,bj)+jsl
748 IF ((j.ge.jmin).and.(j.le.jmax)) THEN
749 IF (OBSs(i,k,bi,bj).ne. 0.d0) THEN
750 IF (useLinearSponge) THEN
751 srelax= OBSs(i,k,bi,bj)
752 ELSE
753 srelax=(
754 & float(spongeThickness-jsl)*OBSs(i,k,bi,bj)
755 & + float(jsl)*salt(i,j,k,bi,bj) )
756 & / float(spongeThickness)
757 ENDIF
758 lambda_obcs_s = (
759 & float(spongeThickness)*Vrelaxobcsbound
760 & + float(jsl)*Vrelaxobcsinner)
761 & / float(spongeThickness)
762 IF (lambda_obcs_s.ne.0.) THEN
763 lambda_obcs_s = 1. _d 0 / lambda_obcs_s
764 ELSE
765 lambda_obcs_s = 0. _d 0
766 ENDIF
767 gS_arr(i,j) = gS_arr(i,j)
768 & - maskC(i,j,k,bi,bj) * lambda_obcs_s
769 & * ( salt(i,j,k,bi,bj) - srelax )
770 ENDIF
771 ENDIF
772 ENDDO
773 ENDIF
774 ENDDO
775 ENDIF
776 #endif /* ALLOW_OBCS_SOUTH */
777
778 C Eastern Open Boundary
779 #ifdef ALLOW_OBCS_EAST
780 IF ( OBCSsponge_E .AND. OBCSsponge_Salt ) THEN
781 #ifdef ALLOW_AUTODIFF_TAMC
782 CADJ STORE OBEs(:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
783 #endif
784 DO j=jMin,jMax
785 IF ( OB_Ie(j,bi,bj).NE.OB_indexNone ) THEN
786 DO isl= 1,spongeThickness
787 i=OB_Ie(j,bi,bj)-isl
788 IF ((i.ge.imin).and.(i.le.imax)) THEN
789 IF (OBEs(j,k,bi,bj).ne. 0.d0) THEN
790 IF (useLinearSponge) THEN
791 srelax= OBEs(j,k,bi,bj)
792 ELSE
793 srelax=(
794 & float(spongeThickness-isl)*OBEs(j,k,bi,bj)
795 & + float(isl)*salt(i,j,k,bi,bj) )
796 & / float(spongeThickness)
797 ENDIF
798 lambda_obcs_s = (
799 & float(spongeThickness-isl)*Urelaxobcsbound
800 & + float(isl)*Urelaxobcsinner)
801 & / float(spongeThickness)
802 IF (lambda_obcs_s.ne.0.) THEN
803 lambda_obcs_s = 1. _d 0 / lambda_obcs_s
804 ELSE
805 lambda_obcs_s = 0. _d 0
806 ENDIF
807 gS_arr(i,j) = gS_arr(i,j)
808 & - maskC(i,j,k,bi,bj) * lambda_obcs_s
809 & * ( salt(i,j,k,bi,bj) - srelax )
810 ENDIF
811 ENDIF
812 ENDDO
813 ENDIF
814 ENDDO
815 ENDIF
816 #endif /* ALLOW_OBCS_EAST */
817
818 C Western Open Boundary
819 #ifdef ALLOW_OBCS_WEST
820 IF ( OBCSsponge_W .AND. OBCSsponge_Salt ) THEN
821 #ifdef ALLOW_AUTODIFF_TAMC
822 CADJ STORE OBWs(:,k,bi,bj) = comlev1_bibj_k, key=kkey, byte=isbyte
823 #endif
824 DO j=jMin,jMax
825 IF ( OB_Iw(j,bi,bj).NE.OB_indexNone ) THEN
826 DO isl= 1,spongeThickness
827 i=OB_Iw(j,bi,bj)+isl
828 IF ((i.ge.imin).and.(i.le.imax)) THEN
829 IF (OBWs(j,k,bi,bj).ne. 0.d0) then
830 IF (useLinearSponge) THEN
831 srelax= OBWs(j,k,bi,bj)
832 ELSE
833 srelax=(
834 & float(spongeThickness-isl)*OBWs(j,k,bi,bj)
835 & + float(isl)*salt(i,j,k,bi,bj) )
836 & / float(spongeThickness)
837 ENDIF
838 lambda_obcs_s= (
839 & float(spongeThickness-isl)*Urelaxobcsbound
840 & + float(isl)*Urelaxobcsinner)
841 & / float(spongeThickness)
842 IF (lambda_obcs_s.ne.0.) THEN
843 lambda_obcs_s = 1. _d 0 / lambda_obcs_s
844 ELSE
845 lambda_obcs_s = 0. _d 0
846 ENDIF
847 gS_arr(i,j) = gS_arr(i,j)
848 & - maskC(i,j,k,bi,bj) * lambda_obcs_s
849 & * ( salt(i,j,k,bi,bj) - srelax )
850 ENDIF
851 ENDIF
852 ENDDO
853 ENDIF
854 ENDDO
855 ENDIF
856 #endif /* ALLOW_OBCS_WEST */
857
858 ENDIF
859
860 #endif /* ALLOW_OBCS & ALLOW_OBCS_SPONGE */
861
862 RETURN
863 END
  ViewVC Help
Powered by ViewVC 1.1.22