ROMS
Loading...
Searching...
No Matches
ana_smflux.h
Go to the documentation of this file.
1!!
2 SUBROUTINE ana_smflux (ng, tile, model)
3!
4!! git $Id$
5!!======================================================================
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 sets kinematic surface momentum flux (wind stress) !
12! "sustr" and "svstr" (m2/s2) using an analytical expression. !
13! !
14!=======================================================================
15!
16 USE mod_param
17 USE mod_forces
18 USE mod_grid
19 USE mod_ncparam
20!
21! Imported variable declarations.
22!
23 integer, intent(in) :: ng, tile, model
24!
25! Local variable declarations.
26!
27 character (len=*), parameter :: MyFile = &
28 & __FILE__
29!
30#include "tile.h"
31!
32 CALL ana_smflux_tile (ng, tile, model, &
33 & lbi, ubi, lbj, ubj, &
34 & imins, imaxs, jmins, jmaxs, &
35 & grid(ng) % angler, &
36#ifdef SPHERICAL
37 & grid(ng) % lonr, &
38 & grid(ng) % latr, &
39#else
40 & grid(ng) % xr, &
41 & grid(ng) % yr, &
42#endif
43#ifdef TL_IOMS
44 & forces(ng) % tl_sustr, &
45 & forces(ng) % tl_svstr, &
46#endif
47 & forces(ng) % sustr, &
48 & forces(ng) % svstr)
49!
50! Set analytical header file name used.
51!
52#ifdef DISTRIBUTE
53 IF (lanafile) THEN
54#else
55 IF (lanafile.and.(tile.eq.0)) THEN
56#endif
57 ananame(24)=myfile
58 END IF
59!
60 RETURN
61 END SUBROUTINE ana_smflux
62!
63!***********************************************************************
64 SUBROUTINE ana_smflux_tile (ng, tile, model, &
65 & LBi, UBi, LBj, UBj, &
66 & IminS, ImaxS, JminS, JmaxS, &
67 & angler, &
68#ifdef SPHERICAL
69 & lonr, latr, &
70#else
71 & xr, yr, &
72#endif
73#ifdef TL_IOMS
74 & tl_sustr, tl_svstr, &
75#endif
76 & sustr, svstr)
77!***********************************************************************
78!
79 USE mod_param
80 USE mod_scalars
81!
83#ifdef DISTRIBUTE
85#endif
86!
87! Imported variable declarations.
88!
89 integer, intent(in) :: ng, tile, model
90 integer, intent(in) :: LBi, UBi, LBj, UBj
91 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
92!
93#ifdef ASSUMED_SHAPE
94 real(r8), intent(in) :: angler(LBi:,LBj:)
95# ifdef SPHERICAL
96 real(r8), intent(in) :: lonr(LBi:,LBj:)
97 real(r8), intent(in) :: latr(LBi:,LBj:)
98# else
99 real(r8), intent(in) :: xr(LBi:,LBj:)
100 real(r8), intent(in) :: yr(LBi:,LBj:)
101# endif
102 real(r8), intent(out) :: sustr(LBi:,LBj:)
103 real(r8), intent(out) :: svstr(LBi:,LBj:)
104# ifdef TL_IOMS
105 real(r8), intent(out) :: tl_sustr(LBi:,LBj:)
106 real(r8), intent(out) :: tl_svstr(LBi:,LBj:)
107# endif
108#else
109 real(r8), intent(in) :: angler(LBi:UBi,LBj:UBj)
110# ifdef SPHERICAL
111 real(r8), intent(in) :: lonr(LBi:UBi,LBj:UBj)
112 real(r8), intent(in) :: latr(LBi:UBi,LBj:UBj)
113# else
114 real(r8), intent(in) :: xr(LBi:UBi,LBj:UBj)
115 real(r8), intent(in) :: yr(LBi:UBi,LBj:UBj)
116# endif
117 real(r8), intent(out) :: sustr(LBi:UBi,LBj:UBj)
118 real(r8), intent(out) :: svstr(LBi:UBi,LBj:UBj)
119# ifdef TL_IOMS
120 real(r8), intent(out) :: tl_sustr(LBi:UBi,LBj:UBj)
121 real(r8), intent(out) :: tl_svstr(LBi:UBi,LBj:UBj)
122# endif
123#endif
124!
125! Local variable declarations.
126!
127 integer :: i, j
128!
129 real(r8) :: Ewind, Nwind, cff, val1, val2, windamp, winddir
130#if defined LAKE_SIGNELL
131 real(r8) :: cff1, mxst, ramp_u, ramp_time, ramp_d
132#endif
133
134#include "set_bounds.h"
135!
136!-----------------------------------------------------------------------
137! Set kinematic surface momentum flux (wind stress) component in the
138! XI-direction (m2/s2) at horizontal U-points.
139!-----------------------------------------------------------------------
140!
141#ifdef BASIN
142 val1=5.0e-05_r8*(1.0_r8+tanh((time(ng)-6.0_r8*86400.0_r8)/ &
143 & (3.0_r8*86400.0_r8)))
144 val2=2.0_r8*pi/el(ng)
145 DO j=jstrt,jendt
146 DO i=istrp,iendt
147 sustr(i,j)=-val1*cos(val2*yr(i,j))
148# ifdef TL_IOMS
149 tl_sustr(i,j)=-val1*cos(val2*yr(i,j))
150# endif
151 END DO
152 END DO
153#elif defined BL_TEST
154 ewind=0.0_r8/rho0
155 nwind=0.3_r8/rho0
156 DO j=jstrt,jendt
157 DO i=istrt,iendt
158 sustr(i,j)=ewind
159# ifdef TL_IOMS
160 tl_sustr(i,j)=ewind
161# endif
162 END DO
163 END DO
164#elif defined CANYON
165 DO j=jstrt,jendt
166 DO i=istrp,iendt
167 sustr(i,j)=5.0e-05_r8*sin(2.0_r8*pi*tdays(ng)/10.0_r8)* &
168 & (1.0_r8-tanh((yr(i,j)-0.5_r8*el(ng))/10000.0_r8))
169# ifdef TL_IOMS
170 tl_sustr(i,j)=5.0e-05_r8*sin(2.0_r8*pi*tdays(ng)/10.0_r8)* &
171 & (1.0_r8-tanh((yr(i,j)-0.5_r8*el(ng))/10000.0_r8))
172# endif
173 END DO
174 END DO
175#elif defined CHANNEL_NECK
176!! IF ((tdays(ng)-dstart).le.4.0_r8) THEN
177!! windamp=-0.01_r8*SIN(pi*(tdays(ng)-dstart)/8.0_r8)/rho0
178!! ELSE
179 windamp=-0.01_r8/rho0
180!! END IF
181 DO j=jstrt,jendt
182 DO i=istrp,iendt
183 sustr(i,j)=windamp
184# ifdef TL_IOMS
185 tl_sustr(i,j)=windamp
186# endif
187 END DO
188 END DO
189#elif defined MIXED_LAYER
190 DO j=jstrt,jendt
191 DO i=istrp,iendt
192 sustr(i,j)=0.0001_r8 ! m2/s2
193# ifdef TL_IOMS
194 tl_sustr(i,j)=0.0001_r8 ! m2/s2
195# endif
196 END DO
197 END DO
198#elif defined DOUBLE_GYRE
199!! windamp=user(1)/rho0
200 windamp=-0.05_r8/rho0
201 val1=2.0_r8*pi/el(ng)
202 DO j=jstrt,jendt
203 DO i=istrp,iendt
204 sustr(i,j)=windamp*cos(val1*yr(i,j))
205# ifdef TL_IOMS
206 tl_sustr(i,j)=windamp*cos(val1*yr(i,j))
207# endif
208 END DO
209 END DO
210#elif defined FLT_TEST
211 DO j=jstrt,jendt
212 DO i=istrp,iendt
213 sustr(i,j)=1.0e-03_r8
214# ifdef TL_IOMS
215 tl_sustr(i,j)=1.0e-03_r8
216# endif
217 END DO
218 END DO
219#elif defined LAKE_SIGNELL
220 mxst=0.2500_r8 ! N/m2
221 ramp_u=15.0_r8 ! start ramp UP at RAMP_UP hours
222 ramp_time=10.0_r8 ! ramp from 0 to 1 over RAMP_TIME hours
223 ramp_d=50.0_r8 ! start ramp DOWN at RAMP_DOWN hours
224 DO j=jstrt,jendt
225 DO i=istrp,iendt
226 cff1=min((0.5_r8*(tanh((time(ng)/3600.0_r8-ramp_u)/ &
227 & (ramp_time/5.0_r8))+1.0_r8)), &
228 & (1.0_r8-(0.5_r8*(tanh((time(ng)/3600.0_r8-ramp_d)/ &
229 & (ramp_time/5.0_r8))+1.0_r8))))
230 sustr(i,j)=mxst/rho0*cff1
231# ifdef TL_IOMS
232 tl_sustr(i,j)=mxst/rho0*cff1
233# endif
234 END DO
235 END DO
236#elif defined LMD_TEST
237 IF (time(ng).le.57600.0_r8) THEN
238 windamp=-0.6_r8*sin(pi*time(ng)/57600.0_r8)* &
239 & sin(2.0_r8*pi*time(ng)/57600.0_r8)/rho0
240 ELSE
241 windamp=0.0_r8
242 END IF
243 DO j=jstrt,jendt
244 DO i=istrp,iendt
245 sustr(i,j)=windamp
246# ifdef TL_IOMS
247 tl_sustr(i,j)=windamp
248# endif
249 END DO
250 END DO
251#elif defined NJ_BIGHT
252!! windamp=0.086824313_r8
253!! winddir=0.5714286_r8
254!! if ((tdays(ng)-dstart).le.0.5_r8) then
255!! Ewind=windamp*winddir*SIN(pi*(tdays(ng)-dstart))/rho0
256!! Nwind=windamp*SIN(pi*(tdays(ng)-dstart))/rho0
257!! else
258!! Ewind=windamp*winddir/rho0
259!! Nwind=windamp/rho0
260!! endif
261 IF ((tdays(ng)-dstart).le.3.0_r8) THEN
262 winddir=60.0_r8
263 windamp=0.1_r8
264 ELSE IF (((tdays(ng)-dstart).gt.3.0_r8).and. &
265 & ((tdays(ng)-dstart).le.4.0_r8)) THEN
266 winddir= 60.0_r8*((tdays(ng)-dstart)-2.0_r8)- &
267 & 120.0_r8*((tdays(ng)-dstart)-2.0_r8)
268 windamp=0.0_r8
269 ELSE
270 winddir=-120.0_r8
271 windamp=0.0_r8
272 END IF
273 ewind=windamp*cos(pi*winddir/180.0_r8)/rho0
274 nwind=windamp*sin(pi*winddir/180.0_r8)/rho0
275 DO j=jstrt,jendt
276 DO i=istrp,iendt
277 val1=0.5_r8*(angler(i-1,j)+angler(i,j))
278 sustr(i,j)=ewind*cos(val1)+nwind*sin(val1)
279# ifdef TL_IOMS
280 tl_sustr(i,j)=ewind*cos(val1)+nwind*sin(val1)
281# endif
282 END DO
283 END DO
284#elif defined SED_TOY
285 DO j=jstrt,jendt
286 DO i=istrp,iendt
287 cff=0.0001_r8
288 IF (time(ng).gt.3000.0_r8) THEN
289 cff=0.0_r8
290 END IF
291 sustr(i,j)=cff
292# ifdef TL_IOMS
293 tl_sustr(i,j)=cff
294# endif
295 END DO
296 END DO
297#elif defined SHOREFACE
298 DO j=jstrt,jendt
299 DO i=istrp,iendt
300 sustr(i,j)=0.0_r8
301# ifdef TL_IOMS
302 tl_sustr(i,j)=0.0_r8
303# endif
304 END DO
305 END DO
306#elif defined UPWELLING
307 IF (nsperiodic(ng)) THEN
308 DO j=jstrt,jendt
309 DO i=istrp,iendt
310 sustr(i,j)=0.0_r8
311# ifdef TL_IOMS
312 tl_sustr(i,j)=0.0_r8
313# endif
314 END DO
315 END DO
316 ELSE IF (ewperiodic(ng)) THEN
317 IF ((tdays(ng)-dstart).le.2.0_r8) THEN
318 windamp=-0.1_r8*sin(pi*(tdays(ng)-dstart)/4.0_r8)/rho0
319 ELSE
320 windamp=-0.1_r8/rho0
321 END IF
322 DO j=jstrt,jendt
323 DO i=istrp,iendt
324 sustr(i,j)=windamp
325# ifdef TL_IOMS
326 tl_sustr(i,j)=windamp
327# endif
328 END DO
329 END DO
330 END IF
331#elif defined WINDBASIN
332 IF ((tdays(ng)-dstart).le.2.0_r8) THEN
333 windamp=-0.1_r8*sin(pi*(tdays(ng)-dstart)/4.0_r8)/rho0
334 ELSE
335 windamp=-0.1_r8/rho0
336 END IF
337 DO j=jstrt,jendt
338 DO i=istrp,iendt
339 sustr(i,j)=windamp
340# ifdef TL_IOMS
341 tl_sustr(i,j)=windamp
342# endif
343 END DO
344 END DO
345#else
346 DO j=jstrt,jendt
347 DO i=istrp,iendt
348 sustr(i,j)=0.0_r8
349# ifdef TL_IOMS
350 tl_sustr(i,j)=0.0_r8
351# endif
352 END DO
353 END DO
354#endif
355!
356!-----------------------------------------------------------------------
357! Set kinematic surface momentum flux (wind stress) component in the
358! ETA-direction (m2/s2) at horizontal V-points.
359!-----------------------------------------------------------------------
360!
361#if defined BL_TEST
362 DO j=jstrt,jendt
363 DO i=istrt,iendt
364 svstr(i,j)=nwind
365# ifdef TL_IOMS
366 tl_svstr(i,j)=nwind
367# endif
368 END DO
369 END DO
370#elif defined LMD_TEST
371 IF (time(ng).le.57600.0_r8) THEN
372 windamp=-0.6_r8*sin(pi*time(ng)/57600.0_r8)* &
373 & cos(2.0_r8*pi*time(ng)/57600.0_r8)/rho0
374 ELSE
375 windamp=0.0_r8
376 END IF
377 DO j=jstrp,jendt
378 DO i=istrt,iendt
379 svstr(i,j)=windamp
380# ifdef TL_IOMS
381 tl_svstr(i,j)=windamp
382# endif
383 END DO
384 END DO
385#elif defined NJ_BIGHT
386 DO j=jstrp,jendt
387 DO i=istrt,iendt
388 val1=0.5_r8*(angler(i,j)+angler(i,j-1))
389 svstr(i,j)=-ewind*sin(val1)+nwind*cos(val1)
390# ifdef TL_IOMS
391 tl_svstr(i,j)=-ewind*sin(val1)+nwind*cos(val1)
392# endif
393 END DO
394 END DO
395#elif defined SED_TOY
396 DO j=jstrp,jendt
397 DO i=istrt,iendt
398 svstr(i,j)=0.0_r8
399# ifdef TL_IOMS
400 tl_svstr(i,j)=0.0_r8
401# endif
402 END DO
403 END DO
404#elif defined SHOREFACE
405 DO j=jstrp,jendt
406 DO i=istrt,iendt
407 svstr(i,j)=0.0_r8
408# ifdef TL_IOMS
409 tl_svstr(i,j)=0.0_r8
410# endif
411 END DO
412 END DO
413#elif defined UPWELLING
414 IF (nsperiodic(ng)) THEN
415 IF ((tdays(ng)-dstart).le.2.0_r8) THEN
416 windamp=-0.1_r8*sin(pi*(tdays(ng)-dstart)/4.0_r8)/rho0
417 ELSE
418 windamp=-0.1_r8/rho0
419 END IF
420 DO j=jstrp,jendt
421 DO i=istrt,iendt
422 svstr(i,j)=windamp
423# ifdef TL_IOMS
424 tl_svstr(i,j)=windamp
425# endif
426 END DO
427 END DO
428 ELSE IF (ewperiodic(ng)) THEN
429 DO j=jstrp,jendt
430 DO i=istrt,iendt
431 svstr(i,j)=0.0_r8
432# ifdef TL_IOMS
433 tl_svstr(i,j)=0.0_r8
434# endif
435 END DO
436 END DO
437 END IF
438#else
439 DO j=jstrp,jendt
440 DO i=istrt,iendt
441 svstr(i,j)=0.0_r8
442# ifdef TL_IOMS
443 tl_svstr(i,j)=0.0_r8
444# endif
445 END DO
446 END DO
447#endif
448!
449! Exchange boundary data.
450!
451 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
452 CALL exchange_u2d_tile (ng, tile, &
453 & lbi, ubi, lbj, ubj, &
454 & sustr)
455#ifdef TL_IOMS
456 CALL exchange_u2d_tile (ng, tile, &
457 & lbi, ubi, lbj, ubj, &
458 & tl_sustr)
459#endif
460 CALL exchange_v2d_tile (ng, tile, &
461 & lbi, ubi, lbj, ubj, &
462 & svstr)
463#ifdef TL_IOMS
464 CALL exchange_v2d_tile (ng, tile, &
465 & lbi, ubi, lbj, ubj, &
466 & tl_svstr)
467#endif
468 END IF
469
470#ifdef DISTRIBUTE
471 CALL mp_exchange2d (ng, tile, model, 2, &
472 & lbi, ubi, lbj, ubj, &
473 & nghostpoints, &
474 & ewperiodic(ng), nsperiodic(ng), &
475 & sustr, svstr)
476# ifdef TL_IOMS
477 CALL mp_exchange2d (ng, tile, model, 2, &
478 & lbi, ubi, lbj, ubj, &
479 & nghostpoints, &
480 & ewperiodic(ng), nsperiodic(ng), &
481 & tl_sustr, tl_svstr)
482# endif
483#endif
484!
485 RETURN
486 END SUBROUTINE ana_smflux_tile
subroutine ana_smflux(ng, tile, model)
Definition ana_smflux.h:3
subroutine ana_smflux_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, angler, lonr, latr, xr, yr, tl_sustr, tl_svstr, sustr, svstr)
Definition ana_smflux.h:77
subroutine exchange_u2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_v2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
logical lanafile
character(len=256), dimension(39) ananame
integer nghostpoints
Definition mod_param.F:710
real(r8), dimension(:), allocatable el
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
real(dp), dimension(:), allocatable tdays
real(dp) dstart
real(dp), dimension(:), allocatable time
real(dp) rho0
real(dp), parameter pi
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)