/[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 178 - (show annotations)
Fri Mar 11 18:47:26 2016 UTC (8 years, 2 months ago) by guez
File size: 7349 byte(s) 
Moved variables date0, deltat, datasz_max, ncvar_ids, point, buff_pos,
buffer, regular from module histcom_var to modules where they are
defined.

Removed procedure ioipslmpp, useless for a sequential program.

Added argument datasz_max to histwrite_real (to avoid circular
dependency with histwrite).

Removed useless variables and computations everywhere.

Changed real litteral constants from default kind to double precision
in lwb, lwu, lwvn, sw1s, swtt, swtt1, swu.

Removed unused arguments: paer of sw, sw1s, sw2s, swclr; pcldsw of
sw1s, sw2s; pdsig, prayl of swr; co2_ppm of clmain, clqh; tsol of
transp_lay; nsrf of screenp; kcrit and kknu of gwstress; pstd of
orosetup.

Added output of relative humidity.
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 jl, jk
50 REAL zcons1, ztmst, 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 ztmst = ptsphy
68 zpi = acos(-1.)
69
70 DO 1001 jl = 1, klon
71 zrho(jl,klev+1) = 0.0
72 pulow(jl) = 0.0
73 pvlow(jl) = 0.0
74 iknub(jl) = klev
75 iknul(jl) = klev
76 ll1(jl,klev+1) = .FALSE.
77 DO 1000 jk = 1, klev
78 pvom(jl,jk) = 0.0
79 pvol(jl,jk) = 0.0
80 pte(jl,jk) = 0.0
81 1000 CONTINUE
82 1001 CONTINUE
83
84
85 !* 2.1 DEFINE LOW LEVEL WIND, PROJECT WINDS IN PLANE OF
86 !* LOW LEVEL WIND, DETERMINE SECTOR IN WHICH TO TAKE
87 !* THE VARIANCE AND SET INDICATOR FOR CRITICAL LEVELS.
88
89
90
91 DO 2006 jk = klev, 1, -1
92 DO 2007 jl = 1, klon
93 IF (ktest(jl)==1) THEN
94 zhcrit(jl,jk) = amax1(ppic(jl)-pmea(jl),100.)
95 zhgeo = pgeom1(jl,jk)/rg
96 ll1(jl,jk) = (zhgeo>zhcrit(jl,jk))
97 IF (ll1(jl,jk) .NEQV. ll1(jl,jk+1)) THEN
98 iknub(jl) = jk
99 END IF
100 END IF
101 2007 CONTINUE
102 2006 CONTINUE
103
104 DO 2010 jl = 1, klon
105 IF (ktest(jl)==1) THEN
106 iknub(jl) = max(iknub(jl),klev/2)
107 iknul(jl) = max(iknul(jl),2*klev/3)
108 IF (iknub(jl)>nktopg) iknub(jl) = nktopg
109 IF (iknub(jl)==nktopg) iknul(jl) = klev
110 IF (iknub(jl)==iknul(jl)) iknub(jl) = iknul(jl) - 1
111 END IF
112 2010 CONTINUE
113
114 ! do 2011 jl=1,klon
115 ! IF(KTEST(JL).EQ.1) THEN
116 ! print *,' iknul= ',iknul(jl),' iknub=',iknub(jl)
117 ! ENDIF
118 !2011 continue
119
120 ! PRINT *,' DANS OROLIFT: 2010'
121
122 DO 223 jk = klev, 2, -1
123 DO 222 jl = 1, klon
124 zrho(jl,jk) = 2.*paphm1(jl,jk)*zcons1/(ptm1(jl,jk)+ptm1(jl,jk-1))
125 222 CONTINUE
126 223 CONTINUE
127 ! PRINT *,' DANS OROLIFT: 223'
128
129 !********************************************************************
130
131 !* DEFINE LOW LEVEL FLOW
132 ! -------------------
133 DO 2115 jk = klev, 1, -1
134 DO 2116 jl = 1, klon
135 IF (ktest(jl)==1) THEN
136 IF (jk>=iknub(jl) .AND. jk<=iknul(jl)) THEN
137 pulow(jl) = pulow(jl) + pum1(jl,jk)*(paphm1(jl,jk+1)-paphm1(jl, &
138 jk))
139 pvlow(jl) = pvlow(jl) + pvm1(jl,jk)*(paphm1(jl,jk+1)-paphm1(jl, &
140 jk))
141 zrho(jl,klev+1) = zrho(jl,klev+1) + zrho(jl,jk)*(paphm1(jl,jk+1) &
142 -paphm1(jl,jk))
143 END IF
144 END IF
145 2116 CONTINUE
146 2115 CONTINUE
147 DO 2110 jl = 1, klon
148 IF (ktest(jl)==1) THEN
149 pulow(jl) = pulow(jl)/(paphm1(jl,iknul(jl)+1)-paphm1(jl,iknub(jl)))
150 pvlow(jl) = pvlow(jl)/(paphm1(jl,iknul(jl)+1)-paphm1(jl,iknub(jl)))
151 zrho(jl,klev+1) = zrho(jl,klev+1)/(paphm1(jl,iknul(jl)+1)-paphm1(jl, &
152 iknub(jl)))
153 END IF
154 2110 CONTINUE
155
156 !* 3. COMPUTE MOUNTAIN LIFT
157
158 DO 301 jl = 1, klon
159 IF (ktest(jl)==1) THEN
160 ztau(jl,klev+1) = -gklift*zrho(jl,klev+1)*2.*romega*2*pvaror(jl)*sin &
161 (zpi/180.*plat(jl))*pvlow(jl)
162 ztav(jl,klev+1) = gklift*zrho(jl,klev+1)*2.*romega*2*pvaror(jl)* &
163 sin(zpi/180.*plat(jl))*pulow(jl)
164 ELSE
165 ztau(jl,klev+1) = 0.0
166 ztav(jl,klev+1) = 0.0
167 END IF
168 301 CONTINUE
169
170
171 !* 4. COMPUTE LIFT PROFILE
172 !* --------------------
173
174
175 DO jk = 1, klev
176 DO jl = 1, klon
177 IF (ktest(jl)==1) THEN
178 ztau(jl,jk) = ztau(jl,klev+1)*paphm1(jl,jk)/paphm1(jl,klev+1)
179 ztav(jl,jk) = ztav(jl,klev+1)*paphm1(jl,jk)/paphm1(jl,klev+1)
180 ELSE
181 ztau(jl,jk) = 0.0
182 ztav(jl,jk) = 0.0
183 END IF
184 end DO
185 end DO
186
187
188 !* 5. COMPUTE TENDENCIES.
189 !* -------------------
190 IF (lifthigh) THEN
191
192 ! PRINT *,' DANS OROLIFT: 500'
193
194 ! EXPLICIT SOLUTION AT ALL LEVELS
195
196 DO 524 jk = 1, klev
197 DO 523 jl = 1, klon
198 IF (ktest(jl)==1) THEN
199 zdelp = paphm1(jl,jk+1) - paphm1(jl,jk)
200 zdudt(jl) = -rg*(ztau(jl,jk+1)-ztau(jl,jk))/zdelp
201 zdvdt(jl) = -rg*(ztav(jl,jk+1)-ztav(jl,jk))/zdelp
202 END IF
203 523 CONTINUE
204 524 CONTINUE
205
206 ! PROJECT PERPENDICULARLY TO U NOT TO DESTROY ENERGY
207
208 DO jk = 1, klev
209 DO jl = 1, klon
210 IF (ktest(jl)==1) THEN
211
212 zslow = sqrt(pulow(jl)**2+pvlow(jl)**2)
213 zsqua = amax1(sqrt(pum1(jl,jk)**2+pvm1(jl,jk)**2),gvsec)
214 zscav = -zdudt(jl)*pvm1(jl,jk) + zdvdt(jl)*pum1(jl,jk)
215 IF (zsqua>gvsec) THEN
216 pvom(jl,jk) = -zscav*pvm1(jl,jk)/zsqua**2
217 pvol(jl,jk) = zscav*pum1(jl,jk)/zsqua**2
218 ELSE
219 pvom(jl,jk) = 0.0
220 pvol(jl,jk) = 0.0
221 END IF
222 zsqua = sqrt(pum1(jl,jk)**2+pum1(jl,jk)**2)
223 IF (zsqua<zslow) THEN
224 pvom(jl,jk) = zsqua/zslow*pvom(jl,jk)
225 pvol(jl,jk) = zsqua/zslow*pvol(jl,jk)
226 END IF
227
228 END IF
229 end DO
230 end DO
231
232 ! 6. LOW LEVEL LIFT, SEMI IMPLICIT:
233 ! ----------------------------------
234
235 ELSE
236
237 DO 601 jl = 1, klon
238 IF (ktest(jl)==1) THEN
239 DO jk = klev, iknub(jl), -1
240 zbet = gklift*2.*romega*sin(zpi/180.*plat(jl))*ztmst* &
241 (pgeom1(jl,iknub(jl)-1)-pgeom1(jl,jk))/ &
242 (pgeom1(jl,iknub(jl)-1)-pgeom1(jl,klev))
243 zdudt(jl) = -pum1(jl,jk)/ztmst/(1+zbet**2)
244 zdvdt(jl) = -pvm1(jl,jk)/ztmst/(1+zbet**2)
245 pvom(jl,jk) = zbet**2*zdudt(jl) - zbet*zdvdt(jl)
246 pvol(jl,jk) = zbet*zdudt(jl) + zbet**2*zdvdt(jl)
247 END DO
248 END IF
249 601 CONTINUE
250
251 END IF
252
253 RETURN
254 END
  ViewVC Help
Powered by ViewVC 1.1.21