Ticket #33: ana_smflux.h

File ana_smflux.h, 13.5 KB (added by m.hadfield, 18 years ago)
Line 
1 SUBROUTINE ana_smflux (ng, tile, model)
2!
3!! svn $Id: ana_smflux.h 34 2007-04-27 04:40:21Z arango $
4!!======================================================================
5!! Copyright (c) 2002-2007 The ROMS/TOMS Group !
6!! Licensed under a MIT/X style license !
7!! See License_ROMS.txt !
8!! !
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#include "tile.h"
26!
27 CALL ana_smflux_tile (ng, model, Istr, Iend, Jstr, Jend, &
28 & LBi, UBi, LBj, UBj, &
29 & GRID(ng) % angler, &
30#ifdef SPHERICAL
31 & GRID(ng) % lonr, &
32 & GRID(ng) % latr, &
33#else
34 & GRID(ng) % xr, &
35 & GRID(ng) % yr, &
36#endif
37#ifdef TL_IOMS
38 & FORCES(ng) % tl_sustr, &
39 & FORCES(ng) % tl_svstr, &
40#endif
41 & FORCES(ng) % sustr, &
42 & FORCES(ng) % svstr)
43!
44! Set analytical header file name used.
45!
46 IF (Lanafile) THEN
47 ANANAME(24)='ROMS/Functionals/ana_smflux.h'
48 END IF
49
50 RETURN
51 END SUBROUTINE ana_smflux
52!
53!***********************************************************************
54 SUBROUTINE ana_smflux_tile (ng, model, Istr, Iend, Jstr, Jend, &
55 & LBi, UBi, LBj, UBj, &
56 & angler, &
57#ifdef SPHERICAL
58 & lonr, latr, &
59#else
60 & xr, yr, &
61#endif
62#ifdef TL_IOMS
63 & tl_sustr, tl_svstr, &
64#endif
65 & sustr, svstr)
66!***********************************************************************
67!
68 USE mod_param
69 USE mod_scalars
70!
71#if defined EW_PERIODIC || defined NS_PERIODIC
72 USE exchange_2d_mod
73#endif
74#ifdef DISTRIBUTE
75 USE mp_exchange_mod, ONLY : mp_exchange2d
76#endif
77!
78! Imported variable declarations.
79!
80 integer, intent(in) :: ng, model, Iend, Istr, Jend, Jstr
81 integer, intent(in) :: LBi, UBi, LBj, UBj
82!
83#ifdef ASSUMED_SHAPE
84 real(r8), intent(in) :: angler(LBi:,LBj:)
85# ifdef SPHERICAL
86 real(r8), intent(in) :: lonr(LBi:,LBj:)
87 real(r8), intent(in) :: latr(LBi:,LBj:)
88# else
89 real(r8), intent(in) :: xr(LBi:,LBj:)
90 real(r8), intent(in) :: yr(LBi:,LBj:)
91# endif
92 real(r8), intent(out) :: sustr(LBi:,LBj:)
93 real(r8), intent(out) :: svstr(LBi:,LBj:)
94# ifdef TL_IOMS
95 real(r8), intent(out) :: tl_sustr(LBi:,LBj:)
96 real(r8), intent(out) :: tl_svstr(LBi:,LBj:)
97# endif
98#else
99 real(r8), intent(in) :: angler(LBi:UBi,LBj:UBj)
100# ifdef SPHERICAL
101 real(r8), intent(in) :: lonr(LBi:UBi,LBj:UBj)
102 real(r8), intent(in) :: latr(LBi:UBi,LBj:UBj)
103# else
104 real(r8), intent(in) :: xr(LBi:UBi,LBj:UBj)
105 real(r8), intent(in) :: yr(LBi:UBi,LBj:UBj)
106# endif
107 real(r8), intent(out) :: sustr(LBi:UBi,LBj:UBj)
108 real(r8), intent(out) :: svstr(LBi:UBi,LBj:UBj)
109# ifdef TL_IOMS
110 real(r8), intent(out) :: tl_sustr(LBi:UBi,LBj:UBj)
111 real(r8), intent(out) :: tl_svstr(LBi:UBi,LBj:UBj)
112# endif
113#endif
114!
115! Local variable declarations.
116!
117#ifdef DISTRIBUTE
118# ifdef EW_PERIODIC
119 logical :: EWperiodic=.TRUE.
120# else
121 logical :: EWperiodic=.FALSE.
122# endif
123# ifdef NS_PERIODIC
124 logical :: NSperiodic=.TRUE.
125# else
126 logical :: NSperiodic=.FALSE.
127# endif
128#endif
129 integer :: IstrR, IendR, JstrR, JendR, IstrU, JstrV
130 integer :: i, j
131 real(r8) :: Ewind, Nwind, cff, val1, val2, windamp, winddir
132#if defined LAKE_SIGNELL
133 real(r8) :: cff1, mxst, ramp_u, ramp_time, ramp_d
134#endif
135
136#include "set_bounds.h"
137!
138!-----------------------------------------------------------------------
139! Set kinematic surface momentum flux (wind stress) component in the
140! XI-direction (m2/s2) at horizontal U-points.
141!-----------------------------------------------------------------------
142!
143#ifdef BASIN
144 val1=5.0E-05_r8*(1.0_r8+TANH((time(ng)-6.0_r8*86400.0_r8)/ &
145 & (3.0_r8*86400.0_r8)))
146 val2=2.0_r8*pi/el(ng)
147 DO j=JstrR,JendR
148 DO i=Istr,IendR
149 sustr(i,j)=-val1*COS(val2*yr(i,j))
150# ifdef TL_IOMS
151 tl_sustr(i,j)=-val1*COS(val2*yr(i,j))
152# endif
153 END DO
154 END DO
155#elif defined BL_TEST
156 Ewind=0.0_r8/rho0
157 Nwind=0.3_r8/rho0
158 DO j=JstrR,JendR
159 DO i=IstrR,IendR
160 sustr(i,j)=Ewind
161# ifdef TL_IOMS
162 tl_sustr(i,j)=Ewind
163# endif
164 END DO
165 END DO
166#elif defined CANYON
167 DO j=JstrR,JendR
168 DO i=Istr,IendR
169 sustr(i,j)=5.0E-05_r8*SIN(2.0_r8*pi*tdays(ng)/10.0_r8)* &
170 & (1.0_r8-TANH((yr(i,j)-0.5_r8*el(ng))/10000.0_r8))
171# ifdef TL_IOMS
172 tl_sustr(i,j)=5.0E-05_r8*SIN(2.0_r8*pi*tdays(ng)/10.0_r8)* &
173 & (1.0_r8-TANH((yr(i,j)-0.5_r8*el(ng))/10000.0_r8))
174# endif
175 END DO
176 END DO
177#elif defined CHANNEL_NECK
178!! IF ((tdays(ng)-dstart).le.4.0_r8) THEN
179!! windamp=-0.01_r8*SIN(pi*(tdays(ng)-dstart)/8.0_r8)/rho0
180!! ELSE
181 windamp=-0.01_r8/rho0
182!! END IF
183 DO j=JstrR,JendR
184 DO i=Istr,IendR
185 sustr(i,j)=windamp
186# ifdef TL_IOMS
187 tl_sustr(i,j)=windamp
188# endif
189 END DO
190 END DO
191#elif defined MIXED_LAYER
192 DO j=JstrR,JendR
193 DO i=Istr,IendR
194 sustr(i,j)=0.0001_r8 ! m2/s2
195# ifdef TL_IOMS
196 tl_sustr(i,j)=0.0001_r8 ! m2/s2
197# endif
198 END DO
199 END DO
200#elif defined DOUBLE_GYRE
201!! windamp=user(1)/rho0
202 windamp=-0.05_r8/rho0
203 val1=2.0_r8*pi/el(ng)
204 DO j=JstrR,JendR
205 DO i=Istr,IendR
206 sustr(i,j)=windamp*COS(val1*yr(i,j))
207# ifdef TL_IOMS
208 tl_sustr(i,j)=windamp*COS(val1*yr(i,j))
209# endif
210 END DO
211 END DO
212#elif defined FLT_TEST
213 DO j=JstrR,JendR
214 DO i=Istr,IendR
215 sustr(i,j)=1.0E-03_r8
216# ifdef TL_IOMS
217 tl_sustr(i,j)=1.0E-03_r8
218# endif
219 END DO
220 END DO
221#elif defined LAKE_SIGNELL
222 mxst=0.2500_r8 ! N/m2
223 ramp_u=15.0_r8 ! start ramp UP at RAMP_UP hours
224 ramp_time=10.0_r8 ! ramp from 0 to 1 over RAMP_TIME hours
225 ramp_d=50.0_r8 ! start ramp DOWN at RAMP_DOWN hours
226 DO j=JstrR,JendR
227 DO i=Istr,IendR
228 cff1=MIN((0.5_r8*(TANH((time(ng)/3600.0_r8-ramp_u)/ &
229 & (ramp_time/5.0_r8))+1.0_r8)), &
230 & (1.0_r8-(0.5_r8*(TANH((time(ng)/3600.0_r8-ramp_d)/ &
231 & (ramp_time/5.0_r8))+1.0_r8))))
232 sustr(i,j)=mxst/rho0*cff1
233# ifdef TL_IOMS
234 tl_sustr(i,j)=mxst/rho0*cff1
235# endif
236 END DO
237 END DO
238#elif defined LMD_TEST
239 IF (time(ng).le.57600.0_r8) THEN
240 windamp=-0.6_r8*SIN(pi*time(ng)/57600.0_r8)* &
241 & SIN(2.0_r8*pi/57600.0_r8)/rho0
242 ELSE
243 windamp=0.0_r8
244 END IF
245 DO j=JstrR,JendR
246 DO i=Istr,IendR
247 sustr(i,j)=windamp
248# ifdef TL_IOMS
249 tl_sustr(i,j)=windamp
250# endif
251 END DO
252 END DO
253#elif defined NJ_BIGHT
254!! windamp=0.086824313_r8
255!! winddir=0.5714286_r8
256!! if ((tdays(ng)-dstart).le.0.5_r8) then
257!! Ewind=windamp*winddir*SIN(pi*(tdays(ng)-dstart))/rho0
258!! Nwind=windamp*SIN(pi*(tdays(ng)-dstart))/rho0
259!! else
260!! Ewind=windamp*winddir/rho0
261!! Nwind=windamp/rho0
262!! endif
263 IF ((tdays(ng)-dstart).le.3.0_r8) THEN
264 winddir=60.0_r8
265 windamp=0.1_r8
266 ELSE IF (((tdays(ng)-dstart).gt.3.0_r8).and. &
267 & ((tdays(ng)-dstart).le.4.0_r8)) THEN
268 winddir= 60.0_r8*((tdays(ng)-dstart)-2.0_r8)- &
269 & 120.0_r8*((tdays(ng)-dstart)-2.0_r8)
270 windamp=0.0_r8
271 ELSE
272 winddir=-120.0_r8
273 windamp=0.0_r8
274 END IF
275 Ewind=windamp*COS(pi*winddir/180.0_r8)/rho0
276 Nwind=windamp*SIN(pi*winddir/180.0_r8)/rho0
277 DO j=JstrR,JendR
278 DO i=Istr,IendR
279 val1=0.5_r8*(angler(i-1,j)+angler(i,j))
280 sustr(i,j)=Ewind*COS(val1)+Nwind*SIN(val1)
281# ifdef TL_IOMS
282 tl_sustr(i,j)=Ewind*COS(val1)+Nwind*SIN(val1)
283# endif
284 END DO
285 END DO
286#elif defined SED_TOY
287 DO j=JstrR,JendR
288 DO i=Istr,IendR
289 cff=0.0001_r8
290 IF (time(ng).gt.3000.0_r8) THEN
291 cff=0.0_r8
292 END IF
293 sustr(i,j)=cff
294# ifdef TL_IOMS
295 tl_sustr(i,j)=cff
296# endif
297 END DO
298 END DO
299#elif defined SHOREFACE
300 DO j=JstrR,JendR
301 DO i=Istr,IendR
302 sustr(i,j)=0.0_r8
303# ifdef TL_IOMS
304 tl_sustr(i,j)=0.0_r8
305# endif
306 END DO
307 END DO
308#elif defined UPWELLING
309 IF ((tdays(ng)-dstart).le.2.0_r8) THEN
310 windamp=-0.1_r8*SIN(pi*(tdays(ng)-dstart)/4.0_r8)/rho0
311 ELSE
312 windamp=-0.1_r8/rho0
313 END IF
314 DO j=JstrR,JendR
315 DO i=Istr,IendR
316 sustr(i,j)=windamp
317# ifdef TL_IOMS
318 tl_sustr(i,j)=windamp
319# endif
320 END DO
321 END DO
322#elif defined WINDBASIN
323 IF ((tdays(ng)-dstart).le.2.0_r8) THEN
324 windamp=-0.1_r8*SIN(pi*(tdays(ng)-dstart)/4.0_r8)/rho0
325 ELSE
326 windamp=-0.1_r8/rho0
327 END IF
328 DO j=JstrR,JendR
329 DO i=Istr,IendR
330 sustr(i,j)=windamp
331# ifdef TL_IOMS
332 tl_sustr(i,j)=windamp
333# endif
334 END DO
335 END DO
336#else
337 DO j=JstrR,JendR
338 DO i=Istr,IendR
339 sustr(i,j)=0.0_r8
340# ifdef TL_IOMS
341 tl_sustr(i,j)=0.0_r8
342# endif
343 END DO
344 END DO
345#endif
346!
347!-----------------------------------------------------------------------
348! Set kinematic surface momentum flux (wind stress) component in the
349! ETA-direction (m2/s2) at horizontal V-points.
350!-----------------------------------------------------------------------
351!
352#if defined BL_TEST
353 DO j=JstrR,JendR
354 DO i=IstrR,IendR
355 svstr(i,j)=Nwind
356# ifdef TL_IOMS
357 tl_svstr(i,j)=Nwind
358# endif
359 END DO
360 END DO
361#elif defined LMD_TEST
362 IF (time(ng).le.57600.0_r8) THEN
363 windamp=-0.6_r8*SIN(pi*time(ng)/57600.0_r8)* &
364 & COS(2.0_r8*pi/57600.0_r8)/rho0
365 ELSE
366 windamp=0.0_r8
367 END IF
368 DO j=Jstr,JendR
369 DO i=IstrR,IendR
370 svstr(i,j)=windamp
371# ifdef TL_IOMS
372 tl_svstr(i,j)=windamp
373# endif
374 END DO
375 END DO
376#elif defined NJ_BIGHT
377 DO j=Jstr,JendR
378 DO i=IstrR,IendR
379 val1=0.5_r8*(angler(i,j)+angler(i,j-1))
380 svstr(i,j)=-Ewind*SIN(val1)+Nwind*COS(val1)
381# ifdef TL_IOMS
382 tl_svstr(i,j)=-Ewind*SIN(val1)+Nwind*COS(val1)
383# endif
384 END DO
385 END DO
386#elif defined SED_TOY
387 DO j=Jstr,JendR
388 DO i=IstrR,IendR
389 svstr(i,j)=0.0_r8
390# ifdef TL_IOMS
391 tl_svstr(i,j)=0.0_r8
392# endif
393 END DO
394 END DO
395#elif defined SHOREFACE
396 DO j=Jstr,JendR
397 DO i=IstrR,IendR
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#else
405 DO j=Jstr,JendR
406 DO i=IstrR,IendR
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#endif
414#if defined EW_PERIODIC || defined NS_PERIODIC
415 CALL exchange_u2d_tile (ng, Istr, Iend, Jstr, Jend, &
416 & LBi, UBi, LBj, UBj, &
417 & sustr)
418 CALL exchange_v2d_tile (ng, Istr, Iend, Jstr, Jend, &
419 & LBi, UBi, LBj, UBj, &
420 & svstr)
421# ifdef TL_IOMS
422 CALL exchange_u2d_tile (ng, Istr, Iend, Jstr, Jend, &
423 & LBi, UBi, LBj, UBj, &
424 & tl_sustr)
425 CALL exchange_v2d_tile (ng, Istr, Iend, Jstr, Jend, &
426 & LBi, UBi, LBj, UBj, &
427 & tl_svstr)
428# endif
429#endif
430#ifdef DISTRIBUTE
431 CALL mp_exchange2d (ng, model, 2, Istr, Iend, Jstr, Jend, &
432 & LBi, UBi, LBj, UBj, &
433 & NghostPoints, EWperiodic, NSperiodic, &
434 & sustr, svstr)
435# ifdef TL_IOMS
436 CALL mp_exchange2d (ng, model, 2, Istr, Iend, Jstr, Jend, &
437 & LBi, UBi, LBj, UBj, &
438 & NghostPoints, EWperiodic, NSperiodic, &
439 & tl_sustr, tl_svstr)
440# endif
441#endif
442 RETURN
443 END SUBROUTINE ana_smflux_tile