/[lmdze]/trunk/Sources/phylmd/Orography/orolift.f
ViewVC logotype

Contents of /trunk/Sources/phylmd/Orography/orolift.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 108 - (show annotations)
Tue Sep 16 14:00:41 2014 UTC (9 years, 7 months ago) by guez
Original Path: trunk/phylmd/Orography/orolift.f
File size: 7571 byte(s)
Imported writefield from LMDZ. Close at the end of gcm the files which
were created by writefiled (not done in LMDZ).

Removed procedures for the output of Grads files. Removed calls to
dump2d. In guide, replaced calls to wrgrads by calls to writefield.

In vlspltqs, removed redundant programming of saturation
pressure. Call foeew from module FCTTRE instead.

Bug fix in interpre: size of w exceeding size of correponding actual
argument wg in advtrac.

In leapfrog, call guide until the end of the run, instead of six hours
before the end.

Bug fix in readsulfate_preind: type of arguments.

1 SUBROUTINE orolift(nlon,nlev,ktest,ptsphy,paphm1,pgeom1,ptm1,pum1,pvm1, &
2 plat,pmea,pvaror,ppic &
3 ,pulow,pvlow,pvom,pvol,pte)
4
5
6 !**** *OROLIFT: SIMULATE THE GEOSTROPHIC LIFT.
7
8 ! PURPOSE.
9 ! --------
10
11 !** INTERFACE.
12 ! ----------
13 ! CALLED FROM *lift_noro
14 ! ----------
15
16 ! AUTHOR.
17 ! -------
18 ! F.LOTT LMD 22/11/95
19
20 USE dimens_m
21 USE dimphy
22 USE suphec_m
23 USE yoegwd
24 IMPLICIT NONE
25
26
27 !-----------------------------------------------------------------------
28
29 !* 0.1 ARGUMENTS
30 ! ---------
31
32
33 INTEGER nlon, nlev
34 REAL pte(nlon,nlev), pvol(nlon,nlev), pvom(nlon,nlev), pulow(nlon), &
35 pvlow(nlon)
36 REAL pum1(nlon,nlev), pvm1(nlon,nlev), ptm1(nlon,nlev)
37 REAL, INTENT (IN) :: plat(nlon)
38 REAL pmea(nlon)
39 REAL, INTENT (IN) :: pvaror(nlon)
40 REAL ppic(nlon), pgeom1(nlon,nlev), paphm1(nlon,nlev+1)
41
42 INTEGER ktest(nlon)
43 REAL, INTENT (IN) :: ptsphy
44 !-----------------------------------------------------------------------
45
46 !* 0.2 LOCAL ARRAYS
47 ! ------------
48 LOGICAL lifthigh
49 INTEGER klevm1, jl, ilevh, jk
50 REAL zcons1, ztmst, zrtmst, zpi, zhgeo
51 REAL zdelp, zslow, zsqua, zscav, zbet
52 INTEGER iknub(klon), iknul(klon)
53 LOGICAL ll1(klon,klev+1)
54
55 REAL ztau(klon,klev+1), ztav(klon,klev+1), zrho(klon,klev+1)
56 REAL zdudt(klon), zdvdt(klon)
57 REAL zhcrit(klon,klev)
58 !-----------------------------------------------------------------------
59
60 !* 1.1 INITIALIZATIONS
61 ! ---------------
62
63 lifthigh = .FALSE.
64
65 IF (nlon/=klon .OR. nlev/=klev) STOP
66 zcons1 = 1./rd
67 klevm1 = klev - 1
68 ztmst = ptsphy
69 zrtmst = 1./ztmst
70 zpi = acos(-1.)
71
72 DO 1001 jl = 1, klon
73 zrho(jl,klev+1) = 0.0
74 pulow(jl) = 0.0
75 pvlow(jl) = 0.0
76 iknub(jl) = klev
77 iknul(jl) = klev
78 ilevh = klev/3
79 ll1(jl,klev+1) = .FALSE.
80 DO 1000 jk = 1, klev
81 pvom(jl,jk) = 0.0
82 pvol(jl,jk) = 0.0
83 pte(jl,jk) = 0.0
84 1000 CONTINUE
85 1001 CONTINUE
86
87
88 !* 2.1 DEFINE LOW LEVEL WIND, PROJECT WINDS IN PLANE OF
89 !* LOW LEVEL WIND, DETERMINE SECTOR IN WHICH TO TAKE
90 !* THE VARIANCE AND SET INDICATOR FOR CRITICAL LEVELS.
91
92
93
94 DO 2006 jk = klev, 1, -1
95 DO 2007 jl = 1, klon
96 IF (ktest(jl)==1) THEN
97 zhcrit(jl,jk) = amax1(ppic(jl)-pmea(jl),100.)
98 zhgeo = pgeom1(jl,jk)/rg
99 ll1(jl,jk) = (zhgeo>zhcrit(jl,jk))
100 IF (ll1(jl,jk) .NEQV. ll1(jl,jk+1)) THEN
101 iknub(jl) = jk
102 END IF
103 END IF
104 2007 CONTINUE
105 2006 CONTINUE
106
107 DO 2010 jl = 1, klon
108 IF (ktest(jl)==1) THEN
109 iknub(jl) = max(iknub(jl),klev/2)
110 iknul(jl) = max(iknul(jl),2*klev/3)
111 IF (iknub(jl)>nktopg) iknub(jl) = nktopg
112 IF (iknub(jl)==nktopg) iknul(jl) = klev
113 IF (iknub(jl)==iknul(jl)) iknub(jl) = iknul(jl) - 1
114 END IF
115 2010 CONTINUE
116
117 ! do 2011 jl=1,klon
118 ! IF(KTEST(JL).EQ.1) THEN
119 ! print *,' iknul= ',iknul(jl),' iknub=',iknub(jl)
120 ! ENDIF
121 !2011 continue
122
123 ! PRINT *,' DANS OROLIFT: 2010'
124
125 DO 223 jk = klev, 2, -1
126 DO 222 jl = 1, klon
127 zrho(jl,jk) = 2.*paphm1(jl,jk)*zcons1/(ptm1(jl,jk)+ptm1(jl,jk-1))
128 222 CONTINUE
129 223 CONTINUE
130 ! PRINT *,' DANS OROLIFT: 223'
131
132 !********************************************************************
133
134 !* DEFINE LOW LEVEL FLOW
135 ! -------------------
136 DO 2115 jk = klev, 1, -1
137 DO 2116 jl = 1, klon
138 IF (ktest(jl)==1) THEN
139 IF (jk>=iknub(jl) .AND. jk<=iknul(jl)) THEN
140 pulow(jl) = pulow(jl) + pum1(jl,jk)*(paphm1(jl,jk+1)-paphm1(jl, &
141 jk))
142 pvlow(jl) = pvlow(jl) + pvm1(jl,jk)*(paphm1(jl,jk+1)-paphm1(jl, &
143 jk))
144 zrho(jl,klev+1) = zrho(jl,klev+1) + zrho(jl,jk)*(paphm1(jl,jk+1) &
145 -paphm1(jl,jk))
146 END IF
147 END IF
148 2116 CONTINUE
149 2115 CONTINUE
150 DO 2110 jl = 1, klon
151 IF (ktest(jl)==1) THEN
152 pulow(jl) = pulow(jl)/(paphm1(jl,iknul(jl)+1)-paphm1(jl,iknub(jl)))
153 pvlow(jl) = pvlow(jl)/(paphm1(jl,iknul(jl)+1)-paphm1(jl,iknub(jl)))
154 zrho(jl,klev+1) = zrho(jl,klev+1)/(paphm1(jl,iknul(jl)+1)-paphm1(jl, &
155 iknub(jl)))
156 END IF
157 2110 CONTINUE
158
159
160 200 CONTINUE
161
162 !***********************************************************
163
164 !* 3. COMPUTE MOUNTAIN LIFT
165
166 300 CONTINUE
167
168 DO 301 jl = 1, klon
169 IF (ktest(jl)==1) THEN
170 ztau(jl,klev+1) = -gklift*zrho(jl,klev+1)*2.*romega*2*pvaror(jl)*sin &
171 (zpi/180.*plat(jl))*pvlow(jl)
172 ztav(jl,klev+1) = gklift*zrho(jl,klev+1)*2.*romega*2*pvaror(jl)* &
173 sin(zpi/180.*plat(jl))*pulow(jl)
174 ELSE
175 ztau(jl,klev+1) = 0.0
176 ztav(jl,klev+1) = 0.0
177 END IF
178 301 CONTINUE
179
180
181 !* 4. COMPUTE LIFT PROFILE
182 !* --------------------
183
184
185 400 CONTINUE
186
187 DO jk = 1, klev
188 DO jl = 1, klon
189 IF (ktest(jl)==1) THEN
190 ztau(jl,jk) = ztau(jl,klev+1)*paphm1(jl,jk)/paphm1(jl,klev+1)
191 ztav(jl,jk) = ztav(jl,klev+1)*paphm1(jl,jk)/paphm1(jl,klev+1)
192 ELSE
193 ztau(jl,jk) = 0.0
194 ztav(jl,jk) = 0.0
195 END IF
196 end DO
197 end DO
198
199
200 !* 5. COMPUTE TENDENCIES.
201 !* -------------------
202 IF (lifthigh) THEN
203
204 500 CONTINUE
205 ! PRINT *,' DANS OROLIFT: 500'
206
207 ! EXPLICIT SOLUTION AT ALL LEVELS
208
209 DO 524 jk = 1, klev
210 DO 523 jl = 1, klon
211 IF (ktest(jl)==1) THEN
212 zdelp = paphm1(jl,jk+1) - paphm1(jl,jk)
213 zdudt(jl) = -rg*(ztau(jl,jk+1)-ztau(jl,jk))/zdelp
214 zdvdt(jl) = -rg*(ztav(jl,jk+1)-ztav(jl,jk))/zdelp
215 END IF
216 523 CONTINUE
217 524 CONTINUE
218
219 ! PROJECT PERPENDICULARLY TO U NOT TO DESTROY ENERGY
220
221 DO jk = 1, klev
222 DO jl = 1, klon
223 IF (ktest(jl)==1) THEN
224
225 zslow = sqrt(pulow(jl)**2+pvlow(jl)**2)
226 zsqua = amax1(sqrt(pum1(jl,jk)**2+pvm1(jl,jk)**2),gvsec)
227 zscav = -zdudt(jl)*pvm1(jl,jk) + zdvdt(jl)*pum1(jl,jk)
228 IF (zsqua>gvsec) THEN
229 pvom(jl,jk) = -zscav*pvm1(jl,jk)/zsqua**2
230 pvol(jl,jk) = zscav*pum1(jl,jk)/zsqua**2
231 ELSE
232 pvom(jl,jk) = 0.0
233 pvol(jl,jk) = 0.0
234 END IF
235 zsqua = sqrt(pum1(jl,jk)**2+pum1(jl,jk)**2)
236 IF (zsqua<zslow) THEN
237 pvom(jl,jk) = zsqua/zslow*pvom(jl,jk)
238 pvol(jl,jk) = zsqua/zslow*pvol(jl,jk)
239 END IF
240
241 END IF
242 end DO
243 end DO
244
245 ! 6. LOW LEVEL LIFT, SEMI IMPLICIT:
246 ! ----------------------------------
247
248 ELSE
249
250 DO 601 jl = 1, klon
251 IF (ktest(jl)==1) THEN
252 DO jk = klev, iknub(jl), -1
253 zbet = gklift*2.*romega*sin(zpi/180.*plat(jl))*ztmst* &
254 (pgeom1(jl,iknub(jl)-1)-pgeom1(jl,jk))/ &
255 (pgeom1(jl,iknub(jl)-1)-pgeom1(jl,klev))
256 zdudt(jl) = -pum1(jl,jk)/ztmst/(1+zbet**2)
257 zdvdt(jl) = -pvm1(jl,jk)/ztmst/(1+zbet**2)
258 pvom(jl,jk) = zbet**2*zdudt(jl) - zbet*zdvdt(jl)
259 pvol(jl,jk) = zbet*zdudt(jl) + zbet**2*zdvdt(jl)
260 END DO
261 END IF
262 601 CONTINUE
263
264 END IF
265
266 RETURN
267 END

  ViewVC Help
Powered by ViewVC 1.1.21