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

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

Parent Directory Parent Directory | Revision Log Revision Log

Revision 178 - (hide 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 guez 23 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 guez 38 USE suphec_m
23 guez 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 guez 178 INTEGER jl, jk
50     REAL zcons1, ztmst, zpi, zhgeo
51 guez 23 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 guez 108 DO jk = 1, klev
176     DO jl = 1, klon
177 guez 23 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 guez 108 end DO
185     end DO
186 guez 23
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 guez 108 DO jk = 1, klev
209     DO jl = 1, klon
210 guez 23 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 guez 108 end DO
230     end DO
231 guez 23
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