source: LMDZ4/branches/LMDZ4V5.0-dev/libf/phylmd/cv3_vertmix.F @ 1373

Last change on this file since 1373 was 879, checked in by Laurent Fairhead, 17 years ago

Suite de la bascule vers une physique avec thermiques, nouvelle convection, poche froide ...
LF

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 5.0 KB
Line 
1      SUBROUTINE cv3_vertmix(len,nd,iflag,plim1,plim2,p,ph,t,q,u,v
2     :                     ,w,wi,nk,tmix,thmix,qmix,qsmix
3     :                     ,umix,vmix,plcl)
4***************************************************************
5*                                                             *
6* CV3_VERTMIX   Brassage adiabatique d'une couche d'epaisseur *
7*               arbitraire.                                   *
8*                                                             *
9* written by   : Grandpeix Jean-Yves, 28/12/2001, 13.14.24    *
10* modified by :  Filiberti M-A 06/2005 vectorisation          *
11***************************************************************
12*
13       implicit none
14C==============================================================
15C
16C vertmix : determine theta et r du melange obtenu en brassant
17C adiabatiquement entre plim1 et plim2, avec une ponderation w.
18C
19C===============================================================
20
21#include "cvthermo.h"
22#include "YOETHF.h"
23#include "YOMCST.h"
24#include "FCTTRE.h"
25c input :
26      integer nd,len
27      integer nk(len),iflag(len)
28      real t(len,nd),q(len,nd),w(nd)
29      real u(len,nd),v(len,nd)
30      real p(len,nd),ph(len,nd+1)
31      real plim1(len),plim2(len)
32c output :
33      real tmix(len),thmix(len),qmix(len),wi(len,nd)
34      real umix(len),vmix(len)
35      real qsmix(len)
36      real plcl(len)
37c internal variables :
38      integer j1(len),j2(len),niflag7
39      real A,B
40      real ahm(len),dpw(len),coef(len)
41      real p1(len,nd),p2(len,nd)
42      real rdcp(len),a2(len),b2(len),pnk(len)
43      real rh(len),chi(len)
44      real cpn
45      real x,y,p0,p0m1,zdelta,zcor
46
47      integer i,j
48
49      do j = 1,nd
50        do i=1,len
51          if (plim1(i).le.ph(i,j)) j1(i) = j
52          if (plim2(i).ge.ph(i,j+1).and.plim2(i).lt.ph(i,j)) j2(i) = j
53        enddo
54      enddo
55c
56      do j=1,nd
57        do i = 1,len
58          wi(i,j) = 0.
59        enddo
60      enddo
61      do i = 1,len
62       ahm(i)=0.
63       qmix(i)=0.
64       umix(i)=0.
65       vmix(i)=0.
66       dpw(i) =0.
67       a2(i)=0.0
68       b2(i) = 0.
69       pnk(i) = p(i,nk(i))
70      enddo
71c
72      p0 = 1000.
73      p0m1 = 1./p0
74c
75      do i=1,len
76        coef(i) = 1./(plim1(i)-plim2(i))
77      end do
78c
79      do  j=1,nd
80        do i=1,len
81          if (j.ge.j1(i).and.j.le.j2(i)) then
82            p1(i,j) = min(ph(i,j),plim1(i))
83            p2(i,j) = max(ph(i,j+1),plim2(i))
84cCRtest:couplage thermiques: deja normalise
85c             wi(i,j) = w(j)
86c             print*,'wi',wi(i,j)
87            wi(i,j) = w(j)*(p1(i,j)-p2(i,j))*coef(i)
88            dpw(i) = dpw(i)+wi(i,j)
89          endif
90        end do
91      end do
92cCR:print
93c      do i=1,len
94c         print*,'plim',plim1(i),plim2(i)
95c      enddo
96      do  j=1,nd
97        do i=1,len
98          if (j.ge.j1(i).and.j.le.j2(i)) then
99            wi(i,j)=wi(i,j)/dpw(i)
100            ahm(i)=ahm(i)+(cpd*(1.-q(i,j))+q(i,j)*cpv)*t(i,j)*wi(i,j)
101            qmix(i)=qmix(i)+q(i,j)*wi(i,j)
102            umix(i)=umix(i)+u(i,j)*wi(i,j)
103            vmix(i)=vmix(i)+v(i,j)*wi(i,j)
104          endif
105        end do
106      end do
107c
108      do i=1,len
109         rdcp(i)=(rrd*(1.-qmix(i))+qmix(i)*rrv)/
110     :            (cpd*(1.-qmix(i))+qmix(i)*cpv)
111      end do
112c
113
114c
115      do 20 j=1,nd
116        do 18 i=1,len
117          if (j.ge.j1(i).and.j.le.j2(i)) then
118cc            x=(.5*(p1(i,j)+p2(i,j))*p0m1)**rdcp(i)
119            y=(.5*(p1(i,j)+p2(i,j))/pnk(i))**rdcp(i)
120cc            a2(i)=a2(i)+(cpd*(1.-qmix(i))+qmix(i)*cpv)*x*wi(i,j)
121            b2(i)=b2(i)+(cpd*(1.-qmix(i))+qmix(i)*cpv)*y*wi(i,j)
122          endif
123   18  continue
124   20  continue
125c
126       do i=1,len
127         tmix(i) = ahm(i)/b2(i)
128         thmix(i) =tmix(i)*(p0/pnk(i))**rdcp(i)
129c         print*,'thmix ahm',ahm(i),b2(i)
130c         print*,'thmix t',tmix(i),p0
131c         print*,'thmix p',pnk(i),rdcp(i)
132c         print*,'thmix',thmix(i)
133cc         thmix(i) = ahm(i)/a2(i)
134cc         tmix(i)= thmix(i)*(pnk(i)*p0m1)**rdcp(i)
135         zdelta=max(0.,sign(1.,rtt-tmix(i)))
136         qsmix(i)= r2es*FOEEW(tmix(i),zdelta)/(pnk(i)*100.)
137         qsmix(i)=min(0.5,qsmix(i))
138         zcor=1./(1.-retv*qsmix(i))
139         qsmix(i)=qsmix(i)*zcor
140       end do
141c
142!-------------------------------------------------------------------
143! --- Calculate lifted condensation level of air at parcel origin level
144! --- (Within 0.2% of formula of Bolton, MON. WEA. REV.,1980)
145!-------------------------------------------------------------------
146
147       A = 1669.0 ! convect3
148       B = 122.0  ! convect3
149
150
151       niflag7=0
152       do 260 i=1,len
153
154        if (iflag(i).ne.7) then ! modif sb Jun7th 2002
155c
156        rh(i)=qmix(i)/qsmix(i)
157        chi(i)=tmix(i)/(A-B*rh(i)-tmix(i)) ! convect3
158c   ATTENTION, la LIGNE DESSOUS A ETE RAJOUTEE ARBITRAIREMENT ET
159c  MASQUE UN PB POTENTIEL
160        chi(i)=max(chi(i),0.)
161        rh(i)=max(rh(i),0.)
162        plcl(i)=pnk(i)*(rh(i)**chi(i))
163        if(((plcl(i).lt.200.0).or.(plcl(i).ge.2000.0))
164     &   .and.(iflag(i).eq.0))iflag(i)=8
165
166        else
167
168          niflag7=niflag7+1
169          plcl(i)=plim2(i)
170c
171        endif ! iflag=7
172
173c      print*,'NIFLAG7  =',niflag7
174
175 260   continue
176
177      return
178      end
179
Note: See TracBrowser for help on using the repository browser.