/[MITgcm]/MITgcm/pkg/exf/exf_wind.F
ViewVC logotype

Contents of /MITgcm/pkg/exf/exf_wind.F

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


Revision 1.4 - (show annotations) (download)
Wed Jan 10 21:44:38 2007 UTC (17 years, 5 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint58u_post, checkpoint58v_post
Changes since 1.3: +37 -43 lines
when wind-stress is given as input (#undef ALLOW_ATM_WIND),
 by-pass turbulent momentum calculation.

1 C $Header: /u/gcmpack/MITgcm/pkg/exf/exf_wind.F,v 1.3 2006/06/03 21:36:15 mlosch Exp $
2 C $Name: $
3
4 #include "EXF_OPTIONS.h"
5
6 subroutine exf_wind(mytime, myiter, mythid)
7
8 c ==================================================================
9 c SUBROUTINE exf_wind
10 c ==================================================================
11 c
12 c o Prepare wind speed and stress fields
13 c
14 c ==================================================================
15 c SUBROUTINE exf_wind
16 c ==================================================================
17
18 implicit none
19
20 c == global variables ==
21
22 #include "SIZE.h"
23 #include "EEPARAMS.h"
24 #include "PARAMS.h"
25 c #include "DYNVARS.h"
26 c #include "GRID.h"
27
28 #include "exf_param.h"
29 #include "exf_fields.h"
30 #include "exf_constants.h"
31
32 #ifdef ALLOW_AUTODIFF_TAMC
33 #include "tamc.h"
34 #endif
35
36 c == routine arguments ==
37
38 integer mythid
39 integer myiter
40 _RL mytime
41
42 c == local variables ==
43
44 integer bi,bj
45 integer i,j
46
47 _RL ustmp
48 c _RL ustar
49 c _RL tmp1, tmp2, tmp3, tmp4, tmp5
50
51 c == external functions ==
52
53 integer ilnblnk
54 external ilnblnk
55
56 c == end of interface ==
57
58 c-- Use atmospheric state to compute surface fluxes.
59
60 c Loop over tiles.
61 do bj = mybylo(mythid),mybyhi(mythid)
62 do bi = mybxlo(mythid),mybxhi(mythid)
63 do j = 1,sny
64 do i = 1,snx
65 c
66 c-- Initialise
67 us(i,j,bi,bj) = 0. _d 0
68 cw(i,j,bi,bj) = 0. _d 0
69 sw(i,j,bi,bj) = 0. _d 0
70 sh(i,j,bi,bj) = 0. _d 0
71 c
72 #ifdef ALLOW_ATM_WIND
73 c Wind speed and direction.
74 ustmp = uwind(i,j,bi,bj)*uwind(i,j,bi,bj) +
75 & vwind(i,j,bi,bj)*vwind(i,j,bi,bj)
76 if ( ustmp .ne. 0. _d 0 ) then
77 us(i,j,bi,bj) = sqrt(ustmp)
78 cw(i,j,bi,bj) = uwind(i,j,bi,bj)/us(i,j,bi,bj)
79 sw(i,j,bi,bj) = vwind(i,j,bi,bj)/us(i,j,bi,bj)
80 else
81 us(i,j,bi,bj) = 0. _d 0
82 cw(i,j,bi,bj) = 0. _d 0
83 sw(i,j,bi,bj) = 0. _d 0
84 endif
85 #else /* ifndef ALLOW_ATM_WIND */
86 c
87 c The variables us, sh and rdn have to be computed from
88 c given wind stresses inverting relationship for neutral
89 c drag coeff. cdn.
90 c The inversion is based on linear and quadratic form of
91 c cdn(umps); ustar can be directly computed from stress;
92 C- no need for wind-stress inversion since everything
93 C (ustar, ... etc ...) is derived directly from wind-stress
94
95 c ustmp = ustress(i,j,bi,bj)*ustress(i,j,bi,bj) +
96 c & vstress(i,j,bi,bj)*vstress(i,j,bi,bj)
97 c if ( ustmp .ne. 0. _d 0 ) then
98 c ustar = sqrt(ustmp/atmrho)
99 c cw(i,j,bi,bj) = ustress(i,j,bi,bj)/sqrt(ustmp)
100 c sw(i,j,bi,bj) = vstress(i,j,bi,bj)/sqrt(ustmp)
101 c else
102 c ustar = 0. _d 0
103 c cw(i,j,bi,bj) = 0. _d 0
104 c sw(i,j,bi,bj) = 0. _d 0
105 c endif
106 c
107 c if ( ustar .eq. 0. _d 0 ) then
108 c us(i,j,bi,bj) = 0. _d 0
109 c else if ( ustar .lt. ustofu11 ) then
110 c tmp1 = -cquadrag_2/cquadrag_1/2
111 c tmp2 = sqrt(tmp1*tmp1 + ustar*ustar/cquadrag_1)
112 c us(i,j,bi,bj) = sqrt(tmp1 + tmp2)
113 c else
114 c tmp3 = clindrag_2/clindrag_1/3
115 c tmp4 = ustar*ustar/clindrag_1/2 - tmp3**3
116 c tmp5 = sqrt(ustar*ustar/clindrag_1*
117 c & (ustar*ustar/clindrag_1/4 - tmp3**3))
118 c us(i,j,bi,bj) = (tmp4 + tmp5)**(1/3) +
119 c & tmp3**2 * (tmp4 + tmp5)**(-1/3) - tmp3
120 c endif
121 #endif /* ifndef ALLOW_ATM_WIND */
122
123 c-- set lower limit
124 sh(i,j,bi,bj) = max(us(i,j,bi,bj),umin)
125
126 c-- if wspeed available, overwrite sh,
127 c-- otherwise fill wspeed array for later use
128 if ( wspeedfile .NE. ' ' ) then
129 us(i,j,bi,bj) = wspeed(i,j,bi,bj)
130 sh(i,j,bi,bj) = max(wspeed(i,j,bi,bj),umin)
131 else
132 wspeed(i,j,bi,bj) = sh(i,j,bi,bj)
133 endif
134
135 enddo
136 enddo
137 enddo
138 enddo
139
140 return
141 end

  ViewVC Help
Powered by ViewVC 1.1.22