ROMS
Loading...
Searching...
No Matches
ana_hmixcoef.h
Go to the documentation of this file.
1!!
2 SUBROUTINE ana_sponge (ng, tile, model)
3!
4!! git $Id$
5!!================================================= Hernan G. Arango ===
6!! Copyright (c) 2002-2025 The ROMS Group !
7!! Licensed under a MIT/X style license !
8!! See License_ROMS.md !
9!=======================================================================
10! !
11! This routine rescales horizontal mixing coefficients according !
12! to the grid size. Also, if applicable, increases horizontal !
13! in sponge areas. !
14! !
15! WARNING: All biharmonic coefficients are assumed to have the !
16! square root taken and have m^2 s^-1/2 units. This !
17! will allow multiplying the biharmonic coefficient !
18! to harmonic operator. !
19! !
20!=======================================================================
21!
22 USE mod_param
23!
24! Imported variable declarations.
25!
26 integer, intent(in) :: ng, tile, model
27!
28! Local variable declarations.
29!
30 character (len=*), parameter :: MyFile = &
31 & __FILE__
32!
33#include "tile.h"
34!
35 CALL ana_sponge_tile (ng, tile, model, &
36 & lbi, ubi, lbj, ubj, &
37 & imins, imaxs, jmins, jmaxs)
38!
39! Set analytical header file name used.
40!
41#ifdef DISTRIBUTE
42 IF (lanafile) THEN
43#else
44 IF (lanafile.and.(tile.eq.0)) THEN
45#endif
46 ananame( 8)=myfile
47 END IF
48!
49 RETURN
50 END SUBROUTINE ana_sponge
51!
52!***********************************************************************
53 SUBROUTINE ana_sponge_tile (ng, tile, model, &
54 & LBi, UBi, LBj, UBj, &
55 & IminS, ImaxS, JminS, JmaxS)
56!***********************************************************************
57!
58 USE mod_param
59 USE mod_grid
60 USE mod_mixing
61 USE mod_scalars
62!
64#ifdef DISTRIBUTE
66# ifdef SOLVE3D
68# endif
69#endif
70!
71! Imported variable declarations.
72!
73 integer, intent(in) :: ng, tile, model
74 integer, intent(in) :: LBi, UBi, LBj, UBj
75 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
76!
77! Local variable declarations.
78!
79 integer :: Iwrk, i, itrc, j
80!
81 real(r8) :: cff, cff1, cff2, fac
82#ifdef WC13
83 real(r8) :: cff_t, cff_s, cff1_t, cff2_t, cff1_s, cff2_s
84#endif
85
86#include "set_bounds.h"
87!
88!-----------------------------------------------------------------------
89! Increase horizontal mixing in the sponge areas.
90!-----------------------------------------------------------------------
91
92#if defined ADRIA02
93!
94! Adriatic Sea southern sponge areas.
95!
96 fac=4.0_r8
97# if defined UV_VIS2
98 DO i=istrt,iendt
99 DO j=jstrt,min(6,jendt)
100 cff=visc2(ng)+real(6-j,r8)*(fac*visc2(ng)-visc2(ng))/6.0_r8
101 mixing(ng) % visc2_r(i,j)=cff
102 mixing(ng) % visc2_p(i,j)=cff
103 END DO
104 DO j=max(jstrt,7),jendt
105 mixing(ng) % visc2_r(i,j)=0.0_r8
106 mixing(ng) % visc2_p(i,j)=0.0_r8
107 END DO
108 END DO
109# endif
110
111# if defined TS_DIF2
112 DO itrc=1,nat
113 DO i=istrt,iendt
114 DO j=jstrt,min(6,jendt)
115 cff=tnu2(itrc,ng)+ &
116 & real(6-j,r8)*(fac*tnu2(itemp,ng)-tnu2(itemp,ng))/6.0_r8
117 mixing(ng) % diff2(i,j,itrc)=cff
118 END DO
119 DO j=max(jstrt,7),jendt
120 mixing(ng) % diff2(i,j,itrc)=0.0_r8
121 END DO
122 END DO
123 END DO
124# endif
125
126#elif defined WC13
127!
128! US West Coast sponge areas.
129!
130 iwrk=int(user(1)) ! same for sponge and nudging layers
131
132# if defined UV_VIS2
133!
134! Momentum sponge regions: sponge viscosities as in Marchesiello
135! et al 2003.
136!
137 cff1=visc2(ng)
138 cff2=100.0_r8
139!
140! Southern edge.
141!
142 DO j=jstrt,min(iwrk,jendt)
143 cff=cff1+real(iwrk-j,r8)*(cff2-cff1)/real(iwrk,r8)
144 DO i=istrt,iendt
145 mixing(ng)%visc2_r(i,j)=max(min(cff,cff2),cff1)
146 mixing(ng)%visc2_p(i,j)=max(min(cff,cff2),cff1)
147 END DO
148 END DO
149!
150! Northern edge.
151!
152 DO j=max(jstrt,mm(ng)+1-iwrk),jendt
153 cff=cff2-real(mm(ng)+1-j,r8)*(cff2-cff1)/real(iwrk,r8)
154 DO i=istrt,iendt
155 mixing(ng) % visc2_r(i,j)=max(min(cff,cff2),cff1)
156 mixing(ng) % visc2_p(i,j)=max(min(cff,cff2),cff1)
157 END DO
158 END DO
159!
160! Western edge.
161!
162 DO i=istrt,min(iwrk,iendt)
163 DO j=max(jstrt,i),min(mm(ng)+1-i,jendt)
164 cff=cff1+real(iwrk-i,r8)*(cff2-cff1)/real(iwrk,r8)
165 mixing(ng) % visc2_r(i,j)=max(min(cff,cff2),cff1)
166 mixing(ng) % visc2_p(i,j)=max(min(cff,cff2),cff1)
167 END DO
168 END DO
169# endif
170
171# if defined TS_DIF2
172!
173! Tracer sponge regions: sponge diffusivities as in Marchesiello
174! et al 2003.
175!
176 cff1_t=tnu2(itemp,ng)
177# ifdef SALINITY
178 cff1_s=tnu2(isalt,ng)
179# endif
180 cff2_t=50.0_r8
181 cff2_s=50.0_r8
182!
183! Southern edge.
184!
185 DO j=jstrt,min(iwrk,jendt)
186 cff_t=cff1_t+real(iwrk-j,r8)*(cff2_t-cff1_t)/real(iwrk,r8)
187# ifdef SALINITY
188 cff_s=cff1_s+real(iwrk-j,r8)*(cff2_s-cff1_s)/real(iwrk,r8)
189# endif
190 DO i=istrt,iendt
191 mixing(ng) % diff2(i,j,itemp)=max(min(cff_t,cff2_t),cff1_t)
192# ifdef SALINITY
193 mixing(ng) % diff2(i,j,isalt)=max(min(cff_s,cff2_s),cff1_s)
194# endif
195 END DO
196 END DO
197!
198! Northern edge.
199!
200 DO j=max(jstrt,mm(ng)+1-iwrk),jendt
201 cff_t=cff2_t-real(mm(ng)+1-j,r8)*(cff2_t-cff1_t)/real(iwrk,r8)
202# ifdef SALINITY
203 cff_s=cff2_s-real(mm(ng)+1-j,r8)*(cff2_s-cff1_s)/real(iwrk,r8)
204# endif
205 DO i=istrt,iendt
206 mixing(ng) % diff2(i,j,itemp)=max(min(cff_t,cff2_t),cff1_t)
207# ifdef SALINITY
208 mixing(ng) % diff2(i,j,isalt)=max(min(cff_s,cff2_s),cff1_s)
209# endif
210 END DO
211 END DO
212!
213! Western edge.
214!
215 DO i=istrt,min(iwrk,iendt)
216 DO j=max(jstrt,i),min(mm(ng)+1-i,jendt)
217 cff_t=cff1_t+real(iwrk-i,r8)*(cff2_t-cff1_t)/real(iwrk,r8)
218# ifdef SALINITY
219 cff_s=cff1_s+real(iwrk-i,r8)*(cff2_s-cff1_s)/real(iwrk,r8)
220# endif
221 mixing(ng) % diff2(i,j,itemp)=max(min(cff_t,cff2_t),cff1_t)
222# ifdef SALINITY
223 mixing(ng) % diff2(i,j,isalt)=max(min(cff_s,cff2_s),cff1_s)
224# endif
225 END DO
226 END DO
227# endif
228#endif
229!
230!-----------------------------------------------------------------------
231! Exchange boundary data.
232!-----------------------------------------------------------------------
233!
234!! WARNING: This section is generic for all applications. Please do not
235!! change the code below.
236!!
237#ifdef UV_VIS2
238 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
239 CALL exchange_r2d_tile (ng, tile, &
240 & lbi, ubi, lbj, ubj, &
241 & mixing(ng) % visc2_r)
242 CALL exchange_p2d_tile (ng, tile, &
243 & lbi, ubi, lbj, ubj, &
244 & mixing(ng) % visc2_p)
245 END IF
246#endif
247
248#ifdef UV_VIS4
249 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
250 CALL exchange_r2d_tile (ng, tile, &
251 & lbi, ubi, lbj, ubj, &
252 & mixing(ng) % visc4_r)
253 CALL exchange_p2d_tile (ng, tile, &
254 & lbi, ubi, lbj, ubj, &
255 & mixing(ng) % visc4_p)
256 END IF
257#endif
258
259#ifdef SOLVE3D
260# ifdef TS_DIF2
261 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
262 DO itrc=1,nt(ng)
263 CALL exchange_r2d_tile (ng, tile, &
264 & lbi, ubi, lbj, ubj, &
265 & mixing(ng) % diff2(:,:,itrc))
266 END DO
267 END IF
268# endif
269
270# ifdef TS_DIF4
271 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
272 DO itrc=1,nt(ng)
273 CALL exchange_r2d_tile (ng, tile, &
274 & lbi, ubi, lbj, ubj, &
275 & mixing(ng) % diff4(:,:,itrc))
276 END DO
277 END IF
278# endif
279#endif
280
281#ifdef DISTRIBUTE
282!
283# ifdef UV_VIS2
284 CALL mp_exchange2d (ng, tile, model, 2, &
285 & lbi, ubi, lbj, ubj, &
286 & nghostpoints, &
287 & ewperiodic(ng), nsperiodic(ng), &
288 & mixing(ng) % visc2_r, &
289 & mixing(ng) % visc2_p)
290# endif
291
292# ifdef UV_VIS4
293 CALL mp_exchange2d (ng, tile, model, 2, &
294 & lbi, ubi, lbj, ubj, &
295 & nghostpoints, &
296 & ewperiodic(ng), nsperiodic(ng), &
297 & mixing(ng) % visc4_r, &
298 & mixing(ng) % visc4_p)
299# endif
300
301# ifdef SOLVE3D
302# ifdef TS_DIF2
303 CALL mp_exchange3d (ng, tile, model, 1, &
304 & lbi, ubi, lbj, ubj, 1, nt(ng), &
305 & nghostpoints, &
306 & ewperiodic(ng), nsperiodic(ng), &
307 & mixing(ng) % diff2)
308# endif
309
310# ifdef TS_DIF4
311 CALL mp_exchange3d (ng, tile, model, 1, &
312 & lbi, ubi, lbj, ubj, 1, nt(ng), &
313 & nghostpoints, &
314 & ewperiodic(ng), nsperiodic(ng), &
315 & mixing(ng) % diff4)
316# endif
317# endif
318#endif
319!
320 RETURN
321 END SUBROUTINE ana_sponge_tile
subroutine ana_sponge(ng, tile, model)
Definition ana_hmixcoef.h:3
subroutine ana_sponge_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs)
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_p2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition exchange_2d.F:66
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
integer nat
Definition mod_param.F:499
integer nghostpoints
Definition mod_param.F:710
integer, dimension(:), allocatable nt
Definition mod_param.F:489
integer, dimension(:), allocatable mm
Definition mod_param.F:456
real(r8), dimension(:,:), allocatable tnu2
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
real(r8), dimension(:), allocatable user
real(r8), dimension(:), allocatable visc2
integer isalt
integer itemp
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)