/[lmdze]/trunk/phylmd/CV_routines/cv_yield.f
ViewVC logotype

Annotation of /trunk/phylmd/CV_routines/cv_yield.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 103 - (hide annotations)
Fri Aug 29 13:00:05 2014 UTC (9 years, 9 months ago) by guez
File size: 13321 byte(s)
Renamed module cvparam to cv_param. Deleted procedure
cv_param. Changed variables of module cv_param into parameters.

In procedures cv_driver, cv_uncompress and cv3_uncompress, removed
some arguments giving dimensions and used module variables klon and
klev instead.

In procedures gradiv2, laplacien_gam and laplacien, changed
declarations of local variables because klevel is not always klev.

Removed code for nudging surface pressure.

Removed arguments pim and pjm of tau2alpha. Added assignment of false
to variable first.

Replaced real argument del of procedures foeew and FOEDE by logical
argument.

1 guez 52 SUBROUTINE cv_yield(nloc,ncum,nd,nk,icb,inb,delt &
2     ,t,q,u,v,gz,p,ph,h,hp,lv,cpn &
3     ,ep,clw,frac,m,mp,qp,up,vp &
4     ,wt,water,evap &
5     ,ment,qent,uent,vent,nent,elij &
6     ,tv,tvp &
7     ,iflag,wd,qprime,tprime &
8     ,precip,cbmf,ft,fq,fu,fv,Ma,qcondc)
9     use cvthermo
10 guez 103 use cv_param
11 guez 52 implicit none
12    
13    
14     ! inputs
15 guez 97 integer, intent(in):: ncum, nd, nloc
16 guez 52 integer nk(nloc), icb(nloc), inb(nloc)
17     integer nent(nloc,nd)
18     real, intent(in):: delt
19     real t(nloc,nd), q(nloc,nd), u(nloc,nd), v(nloc,nd)
20     real gz(nloc,nd)
21     real p(nloc,nd), ph(nloc,nd+1), h(nloc,nd)
22     real hp(nloc,nd), lv(nloc,nd)
23     real cpn(nloc,nd), ep(nloc,nd), clw(nloc,nd), frac(nloc)
24     real m(nloc,nd), mp(nloc,nd), qp(nloc,nd)
25     real up(nloc,nd), vp(nloc,nd)
26     real wt(nloc,nd), water(nloc,nd), evap(nloc,nd)
27     real ment(nloc,nd,nd), qent(nloc,nd,nd), elij(nloc,nd,nd)
28     real uent(nloc,nd,nd), vent(nloc,nd,nd)
29     real tv(nloc,nd), tvp(nloc,nd)
30    
31     ! outputs
32     integer iflag(nloc) ! also an input
33     real cbmf(nloc) ! also an input
34     real wd(nloc), tprime(nloc), qprime(nloc)
35     real precip(nloc)
36     real ft(nloc,nd), fq(nloc,nd), fu(nloc,nd), fv(nloc,nd)
37     real Ma(nloc,nd)
38     real qcondc(nloc,nd)
39    
40     ! local variables
41     integer i,j,ij,k,num1
42     real dpinv,cpinv,awat,fqold,ftold,fuold,fvold,delti
43     real work(nloc), am(nloc),amp1(nloc),ad(nloc)
44     real ents(nloc), uav(nloc),vav(nloc),lvcp(nloc,nd)
45     real qcond(nloc,nd), nqcond(nloc,nd), wa(nloc,nd) ! cld
46     real siga(nloc,nd), ax(nloc,nd), mac(nloc,nd) ! cld
47    
48    
49     ! -- initializations:
50    
51     delti = 1.0/delt
52    
53     do 160 i=1,ncum
54     precip(i)=0.0
55     wd(i)=0.0
56     tprime(i)=0.0
57     qprime(i)=0.0
58     do 170 k=1,nl+1
59     ft(i,k)=0.0
60     fu(i,k)=0.0
61     fv(i,k)=0.0
62     fq(i,k)=0.0
63     lvcp(i,k)=lv(i,k)/cpn(i,k)
64     qcondc(i,k)=0.0 ! cld
65     qcond(i,k)=0.0 ! cld
66     nqcond(i,k)=0.0 ! cld
67     170 continue
68     160 continue
69    
70     !
71     ! *** Calculate surface precipitation in mm/day ***
72     !
73     do 1190 i=1,ncum
74     if(iflag(i).le.1)then
75     precip(i) = wt(i,1)*sigd*water(i,1)*86400/g
76     endif
77     1190 continue
78     !
79     !
80     ! *** Calculate downdraft velocity scale and surface temperature and ***
81     ! *** water vapor fluctuations ***
82     !
83     do i=1,ncum
84     wd(i)=betad*abs(mp(i,icb(i)))*0.01*rrd*t(i,icb(i)) &
85     /(sigd*p(i,icb(i)))
86     qprime(i)=0.5*(qp(i,1)-q(i,1))
87     tprime(i)=lv0*qprime(i)/cpd
88     enddo
89     !
90     ! *** Calculate tendencies of lowest level potential temperature ***
91     ! *** and mixing ratio ***
92     !
93     do 1200 i=1,ncum
94     work(i)=0.01/(ph(i,1)-ph(i,2))
95     am(i)=0.0
96     1200 continue
97     do 1220 k=2,nl
98     do 1210 i=1,ncum
99     if((nk(i).eq.1).and.(k.le.inb(i)).and.(nk(i).eq.1))then
100     am(i)=am(i)+m(i,k)
101     endif
102     1210 continue
103     1220 continue
104     do 1240 i=1,ncum
105     if((g*work(i)*am(i)).ge.delti)iflag(i)=1
106     ft(i,1)=ft(i,1)+g*work(i)*am(i)*(t(i,2)-t(i,1) &
107     +(gz(i,2)-gz(i,1))/cpn(i,1))
108     ft(i,1)=ft(i,1)-lvcp(i,1)*sigd*evap(i,1)
109     ft(i,1)=ft(i,1)+sigd*wt(i,2)*(cl-cpd)*water(i,2)*(t(i,2) &
110     -t(i,1))*work(i)/cpn(i,1)
111     fq(i,1)=fq(i,1)+g*mp(i,2)*(qp(i,2)-q(i,1))* &
112     work(i)+sigd*evap(i,1)
113     fq(i,1)=fq(i,1)+g*am(i)*(q(i,2)-q(i,1))*work(i)
114     fu(i,1)=fu(i,1)+g*work(i)*(mp(i,2)*(up(i,2)-u(i,1)) &
115     +am(i)*(u(i,2)-u(i,1)))
116     fv(i,1)=fv(i,1)+g*work(i)*(mp(i,2)*(vp(i,2)-v(i,1)) &
117     +am(i)*(v(i,2)-v(i,1)))
118     1240 continue
119     do 1260 j=2,nl
120     do 1250 i=1,ncum
121     if(j.le.inb(i))then
122     fq(i,1)=fq(i,1) &
123     +g*work(i)*ment(i,j,1)*(qent(i,j,1)-q(i,1))
124     fu(i,1)=fu(i,1) &
125     +g*work(i)*ment(i,j,1)*(uent(i,j,1)-u(i,1))
126     fv(i,1)=fv(i,1) &
127     +g*work(i)*ment(i,j,1)*(vent(i,j,1)-v(i,1))
128     endif
129     1250 continue
130     1260 continue
131     !
132     ! *** Calculate tendencies of potential temperature and mixing ratio ***
133     ! *** at levels above the lowest level ***
134     !
135     ! *** First find the net saturated updraft and downdraft mass fluxes ***
136     ! *** through each level ***
137     !
138     do 1500 i=2,nl+1
139     !
140     num1=0
141     do 1265 ij=1,ncum
142     if(i.le.inb(ij))num1=num1+1
143     1265 continue
144     if(num1.le.0)go to 1500
145     !
146     call zilch(amp1,ncum)
147     call zilch(ad,ncum)
148     !
149     do 1280 k=i+1,nl+1
150     do 1270 ij=1,ncum
151     if((i.ge.nk(ij)).and.(i.le.inb(ij)) &
152     .and.(k.le.(inb(ij)+1)))then
153     amp1(ij)=amp1(ij)+m(ij,k)
154     endif
155     1270 continue
156     1280 continue
157     !
158     do 1310 k=1,i
159     do 1300 j=i+1,nl+1
160     do 1290 ij=1,ncum
161     if((j.le.(inb(ij)+1)).and.(i.le.inb(ij)))then
162     amp1(ij)=amp1(ij)+ment(ij,k,j)
163     endif
164     1290 continue
165     1300 continue
166     1310 continue
167     do 1340 k=1,i-1
168     do 1330 j=i,nl+1
169     do 1320 ij=1,ncum
170     if((i.le.inb(ij)).and.(j.le.inb(ij)))then
171     ad(ij)=ad(ij)+ment(ij,j,k)
172     endif
173     1320 continue
174     1330 continue
175     1340 continue
176     !
177     do 1350 ij=1,ncum
178     if(i.le.inb(ij))then
179     dpinv=0.01/(ph(ij,i)-ph(ij,i+1))
180     cpinv=1.0/cpn(ij,i)
181     !
182     ft(ij,i)=ft(ij,i) &
183     +g*dpinv*(amp1(ij)*(t(ij,i+1)-t(ij,i) &
184     +(gz(ij,i+1)-gz(ij,i))*cpinv) &
185     -ad(ij)*(t(ij,i)-t(ij,i-1)+(gz(ij,i)-gz(ij,i-1))*cpinv)) &
186     -sigd*lvcp(ij,i)*evap(ij,i)
187     ft(ij,i)=ft(ij,i)+g*dpinv*ment(ij,i,i)*(hp(ij,i)-h(ij,i)+ &
188     t(ij,i)*(cpv-cpd)*(q(ij,i)-qent(ij,i,i)))*cpinv
189     ft(ij,i)=ft(ij,i)+sigd*wt(ij,i+1)*(cl-cpd)*water(ij,i+1)* &
190     (t(ij,i+1)-t(ij,i))*dpinv*cpinv
191     fq(ij,i)=fq(ij,i)+g*dpinv*(amp1(ij)*(q(ij,i+1)-q(ij,i))- &
192     ad(ij)*(q(ij,i)-q(ij,i-1)))
193     fu(ij,i)=fu(ij,i)+g*dpinv*(amp1(ij)*(u(ij,i+1)-u(ij,i))- &
194     ad(ij)*(u(ij,i)-u(ij,i-1)))
195     fv(ij,i)=fv(ij,i)+g*dpinv*(amp1(ij)*(v(ij,i+1)-v(ij,i))- &
196     ad(ij)*(v(ij,i)-v(ij,i-1)))
197     endif
198     1350 continue
199     do 1370 k=1,i-1
200     do 1360 ij=1,ncum
201     if(i.le.inb(ij))then
202     awat=elij(ij,k,i)-(1.-ep(ij,i))*clw(ij,i)
203     awat=max(awat,0.0)
204     fq(ij,i)=fq(ij,i) &
205     +g*dpinv*ment(ij,k,i)*(qent(ij,k,i)-awat-q(ij,i))
206     fu(ij,i)=fu(ij,i) &
207     +g*dpinv*ment(ij,k,i)*(uent(ij,k,i)-u(ij,i))
208     fv(ij,i)=fv(ij,i) &
209     +g*dpinv*ment(ij,k,i)*(vent(ij,k,i)-v(ij,i))
210     ! (saturated updrafts resulting from mixing) ! cld
211     qcond(ij,i)=qcond(ij,i)+(elij(ij,k,i)-awat) ! cld
212     nqcond(ij,i)=nqcond(ij,i)+1. ! cld
213     endif
214     1360 continue
215     1370 continue
216     do 1390 k=i,nl+1
217     do 1380 ij=1,ncum
218     if((i.le.inb(ij)).and.(k.le.inb(ij)))then
219     fq(ij,i)=fq(ij,i) &
220     +g*dpinv*ment(ij,k,i)*(qent(ij,k,i)-q(ij,i))
221     fu(ij,i)=fu(ij,i) &
222     +g*dpinv*ment(ij,k,i)*(uent(ij,k,i)-u(ij,i))
223     fv(ij,i)=fv(ij,i) &
224     +g*dpinv*ment(ij,k,i)*(vent(ij,k,i)-v(ij,i))
225     endif
226     1380 continue
227     1390 continue
228     do 1400 ij=1,ncum
229     if(i.le.inb(ij))then
230     fq(ij,i)=fq(ij,i) &
231     +sigd*evap(ij,i)+g*(mp(ij,i+1)* &
232     (qp(ij,i+1)-q(ij,i)) &
233     -mp(ij,i)*(qp(ij,i)-q(ij,i-1)))*dpinv
234     fu(ij,i)=fu(ij,i) &
235     +g*(mp(ij,i+1)*(up(ij,i+1)-u(ij,i))-mp(ij,i)* &
236     (up(ij,i)-u(ij,i-1)))*dpinv
237     fv(ij,i)=fv(ij,i) &
238     +g*(mp(ij,i+1)*(vp(ij,i+1)-v(ij,i))-mp(ij,i)* &
239     (vp(ij,i)-v(ij,i-1)))*dpinv
240     ! (saturated downdrafts resulting from mixing) ! cld
241     do k=i+1,inb(ij) ! cld
242     qcond(ij,i)=qcond(ij,i)+elij(ij,k,i) ! cld
243     nqcond(ij,i)=nqcond(ij,i)+1. ! cld
244     enddo ! cld
245     ! (particular case: no detraining level is found) ! cld
246     if (nent(ij,i).eq.0) then ! cld
247     qcond(ij,i)=qcond(ij,i)+(1.-ep(ij,i))*clw(ij,i) ! cld
248     nqcond(ij,i)=nqcond(ij,i)+1. ! cld
249     endif ! cld
250     if (nqcond(ij,i).ne.0.) then ! cld
251     qcond(ij,i)=qcond(ij,i)/nqcond(ij,i) ! cld
252     endif ! cld
253     endif
254     1400 continue
255     1500 continue
256     !
257     ! *** Adjust tendencies at top of convection layer to reflect ***
258     ! *** actual position of the level zero cape ***
259     !
260     do 503 ij=1,ncum
261     fqold=fq(ij,inb(ij))
262     fq(ij,inb(ij))=fq(ij,inb(ij))*(1.-frac(ij))
263     fq(ij,inb(ij)-1)=fq(ij,inb(ij)-1) &
264     +frac(ij)*fqold*((ph(ij,inb(ij))-ph(ij,inb(ij)+1))/ &
265     (ph(ij,inb(ij)-1)-ph(ij,inb(ij))))*lv(ij,inb(ij)) &
266     /lv(ij,inb(ij)-1)
267     ftold=ft(ij,inb(ij))
268     ft(ij,inb(ij))=ft(ij,inb(ij))*(1.-frac(ij))
269     ft(ij,inb(ij)-1)=ft(ij,inb(ij)-1) &
270     +frac(ij)*ftold*((ph(ij,inb(ij))-ph(ij,inb(ij)+1))/ &
271     (ph(ij,inb(ij)-1)-ph(ij,inb(ij))))*cpn(ij,inb(ij)) &
272     /cpn(ij,inb(ij)-1)
273     fuold=fu(ij,inb(ij))
274     fu(ij,inb(ij))=fu(ij,inb(ij))*(1.-frac(ij))
275     fu(ij,inb(ij)-1)=fu(ij,inb(ij)-1) &
276     +frac(ij)*fuold*((ph(ij,inb(ij))-ph(ij,inb(ij)+1))/ &
277     (ph(ij,inb(ij)-1)-ph(ij,inb(ij))))
278     fvold=fv(ij,inb(ij))
279     fv(ij,inb(ij))=fv(ij,inb(ij))*(1.-frac(ij))
280     fv(ij,inb(ij)-1)=fv(ij,inb(ij)-1) &
281     +frac(ij)*fvold*((ph(ij,inb(ij))-ph(ij,inb(ij)+1))/ &
282     (ph(ij,inb(ij)-1)-ph(ij,inb(ij))))
283     503 continue
284     !
285     ! *** Very slightly adjust tendencies to force exact ***
286     ! *** enthalpy, momentum and tracer conservation ***
287     !
288     do 682 ij=1,ncum
289     ents(ij)=0.0
290     uav(ij)=0.0
291     vav(ij)=0.0
292     do 681 i=1,inb(ij)
293     ents(ij)=ents(ij) &
294     +(cpn(ij,i)*ft(ij,i)+lv(ij,i)*fq(ij,i))*(ph(ij,i)-ph(ij,i+1))
295     uav(ij)=uav(ij)+fu(ij,i)*(ph(ij,i)-ph(ij,i+1))
296     vav(ij)=vav(ij)+fv(ij,i)*(ph(ij,i)-ph(ij,i+1))
297     681 continue
298     682 continue
299     do 683 ij=1,ncum
300     ents(ij)=ents(ij)/(ph(ij,1)-ph(ij,inb(ij)+1))
301     uav(ij)=uav(ij)/(ph(ij,1)-ph(ij,inb(ij)+1))
302     vav(ij)=vav(ij)/(ph(ij,1)-ph(ij,inb(ij)+1))
303     683 continue
304     do 642 ij=1,ncum
305     do 641 i=1,inb(ij)
306     ft(ij,i)=ft(ij,i)-ents(ij)/cpn(ij,i)
307     fu(ij,i)=(1.-cu)*(fu(ij,i)-uav(ij))
308     fv(ij,i)=(1.-cu)*(fv(ij,i)-vav(ij))
309     641 continue
310     642 continue
311     !
312     do 1810 k=1,nl+1
313     do 1800 i=1,ncum
314     if((q(i,k)+delt*fq(i,k)).lt.0.0)iflag(i)=10
315     1800 continue
316     1810 continue
317     !
318     !
319     do 1900 i=1,ncum
320     if(iflag(i).gt.2)then
321     precip(i)=0.0
322     cbmf(i)=0.0
323     endif
324     1900 continue
325     do 1920 k=1,nl
326     do 1910 i=1,ncum
327     if(iflag(i).gt.2)then
328     ft(i,k)=0.0
329     fq(i,k)=0.0
330     fu(i,k)=0.0
331     fv(i,k)=0.0
332     qcondc(i,k)=0.0 ! cld
333     endif
334     1910 continue
335     1920 continue
336    
337     do k=1,nl+1
338     do i=1,ncum
339     Ma(i,k) = 0.
340     enddo
341     enddo
342     do k=nl,1,-1
343     do i=1,ncum
344     Ma(i,k) = Ma(i,k+1)+m(i,k)
345     enddo
346     enddo
347    
348     !
349     ! *** diagnose the in-cloud mixing ratio ***
350     ! *** of condensed water ***
351     !
352     DO ij=1,ncum
353     do i=1,nd
354     mac(ij,i)=0.0
355     wa(ij,i)=0.0
356     siga(ij,i)=0.0
357     enddo
358     do i=nk(ij),inb(ij)
359     do k=i+1,inb(ij)+1
360     mac(ij,i)=mac(ij,i)+m(ij,k)
361     enddo
362     enddo
363     do i=icb(ij),inb(ij)-1
364     ax(ij,i)=0.
365     do j=icb(ij),i
366     ax(ij,i)=ax(ij,i)+rrd*(tvp(ij,j)-tv(ij,j)) &
367     *(ph(ij,j)-ph(ij,j+1))/p(ij,j)
368     enddo
369     if (ax(ij,i).gt.0.0) then
370     wa(ij,i)=sqrt(2.*ax(ij,i))
371     endif
372     enddo
373     do i=1,nl
374     if (wa(ij,i).gt.0.0) &
375     siga(ij,i)=mac(ij,i)/wa(ij,i) &
376     *rrd*tvp(ij,i)/p(ij,i)/100./delta
377     siga(ij,i) = min(siga(ij,i),1.0)
378     qcondc(ij,i)=siga(ij,i)*clw(ij,i)*(1.-ep(ij,i)) &
379     + (1.-siga(ij,i))*qcond(ij,i)
380     enddo
381     ENDDO
382    
383     return
384     end

  ViewVC Help
Powered by ViewVC 1.1.21