source: LMDZ4/branches/LMDZ4-dev-20091210/libf/phylmd/thermcell_dry.F90 @ 5444

Last change on this file since 5444 was 938, checked in by lmdzadmin, 17 years ago

Enleve prints par defaut
IM
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 7.0 KB
Line 
1       SUBROUTINE thermcell_dry(ngrid,nlay,zlev,pphi,ztv,alim_star,  &
2     &                            lalim,lmin,zmax,wmax,lev_out)
3
4!--------------------------------------------------------------------------
5!thermcell_dry: calcul de zmax et wmax du thermique sec
6!--------------------------------------------------------------------------
7       IMPLICIT NONE
8#include "YOMCST.h"       
9#include "iniprint.h"
10       INTEGER l,ig
11
12       INTEGER ngrid,nlay
13       REAL zlev(ngrid,nlay+1)
14       REAL pphi(ngrid,nlay)
15       REAl ztv(ngrid,nlay)
16       REAL alim_star(ngrid,nlay)
17       INTEGER lalim(ngrid)
18      integer lev_out                           ! niveau pour les print
19
20       REAL zmax(ngrid)
21       REAL wmax(ngrid)
22
23!variables locales
24       REAL zw2(ngrid,nlay+1)
25       REAL f_star(ngrid,nlay+1)
26       REAL ztva(ngrid,nlay+1)
27       REAL wmaxa(ngrid)
28       REAL wa_moy(ngrid,nlay+1)
29       REAL linter(ngrid),zlevinter(ngrid)
30       INTEGER lmix(ngrid),lmax(ngrid),lmin(ngrid)
31
32!initialisations
33       do ig=1,ngrid
34          do l=1,nlay+1
35             zw2(ig,l)=0.
36             wa_moy(ig,l)=0.
37          enddo
38       enddo
39       do ig=1,ngrid
40          do l=1,nlay
41             ztva(ig,l)=ztv(ig,l)
42          enddo
43       enddo
44       do ig=1,ngrid
45          wmax(ig)=0.
46          wmaxa(ig)=0.
47       enddo
48!calcul de la vitesse a partir de la CAPE en melangeant thetav
49
50!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
51! A eliminer
52! Ce if complique etait fait pour reperer la premiere couche instable
53! Ici, c'est lmin.
54!
55!       do l=1,nlay-2
56!         do ig=1,ngrid
57!            if (ztv(ig,l).gt.ztv(ig,l+1)  &
58!     &         .and.alim_star(ig,l).gt.1.e-10  &
59!     &         .and.zw2(ig,l).lt.1e-10) then
60!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
61
62
63! Calcul des F^*, integrale verticale de E^*
64       f_star(:,1)=0.
65       do l=1,nlay
66          f_star(:,l+1)=f_star(:,l)+alim_star(:,l)
67       enddo
68
69! niveau (reel) auquel zw2 s'annule FH :n'etait pas initialise
70       linter(:)=0.
71
72! couche la plus haute concernee par le thermique.
73       lmax(:)=1
74
75! Le niveau linter est une variable continue qui se trouve dans la couche
76! lmax
77
78       do l=1,nlay-2
79         do ig=1,ngrid
80            if (l.eq.lmin(ig).and.lalim(ig).gt.1) then
81
82!------------------------------------------------------------------------
83!  Calcul de la vitesse en haut de la premiere couche instable.
84!  Premiere couche du panache thermique
85!------------------------------------------------------------------------
86               zw2(ig,l+1)=2.*RG*(ztv(ig,l)-ztv(ig,l+1))/ztv(ig,l+1)  &
87     &                     *(zlev(ig,l+1)-zlev(ig,l))  &
88     &                     *0.4*pphi(ig,l)/(pphi(ig,l+1)-pphi(ig,l))
89
90!------------------------------------------------------------------------
91! Tant que la vitesse en bas de la couche et la somme du flux de masse
92! et de l'entrainement (c'est a dire le flux de masse en haut) sont
93! positifs, on calcul
94! 1. le flux de masse en haut  f_star(ig,l+1)
95! 2. la temperature potentielle virtuelle dans la couche ztva(ig,l)
96! 3. la vitesse au carré en haut zw2(ig,l+1)
97!------------------------------------------------------------------------
98
99!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
100!  A eliminer : dans cette version, si zw2 est > 0 on a un therique.
101!  et donc, au dessus, f_star(ig,l+1) est forcement suffisamment
102!  grand puisque on n'a pas de detrainement.
103!  f_star est une fonction croissante.
104!  c'est donc vraiment sur zw2 uniquement qu'il faut faire le test.
105!           else if ((zw2(ig,l).ge.1e-10).and.  &
106!    &               (f_star(ig,l)+alim_star(ig,l).gt.1.e-10)) then
107!              f_star(ig,l+1)=f_star(ig,l)+alim_star(ig,l)
108!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
109
110            else if (zw2(ig,l).ge.1e-10) then
111
112               ztva(ig,l)=(f_star(ig,l)*ztva(ig,l-1)+alim_star(ig,l)  &
113     &                    *ztv(ig,l))/f_star(ig,l+1)
114               zw2(ig,l+1)=zw2(ig,l)*(f_star(ig,l)/f_star(ig,l+1))**2+  &
115     &                     2.*RG*(ztva(ig,l)-ztv(ig,l))/ztv(ig,l)  &
116     &                     *(zlev(ig,l+1)-zlev(ig,l))
117            endif
118! determination de zmax continu par interpolation lineaire
119!------------------------------------------------------------------------
120
121            if (zw2(ig,l+1)>0. .and. zw2(ig,l+1).lt.1.e-10) then
122!               stop'On tombe sur le cas particulier de thermcell_dry'
123!               print*,'On tombe sur le cas particulier de thermcell_dry'
124                zw2(ig,l+1)=0.
125                linter(ig)=l+1
126                lmax(ig)=l
127            endif
128
129            if (zw2(ig,l+1).lt.0.) then
130               linter(ig)=(l*(zw2(ig,l+1)-zw2(ig,l))  &
131     &           -zw2(ig,l))/(zw2(ig,l+1)-zw2(ig,l))
132               zw2(ig,l+1)=0.
133               lmax(ig)=l
134            endif
135
136               wa_moy(ig,l+1)=sqrt(zw2(ig,l+1))
137
138            if (wa_moy(ig,l+1).gt.wmaxa(ig)) then
139!   lmix est le niveau de la couche ou w (wa_moy) est maximum
140               lmix(ig)=l+1
141               wmaxa(ig)=wa_moy(ig,l+1)
142            endif
143         enddo
144      enddo
145       if (prt_level.ge.1) print*,'fin calcul zw2'
146!
147!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
148! A eliminer :
149! Ce calcul de lmax est fait en meme temps que celui de linter, plus haut
150! Calcul de la couche correspondant a la hauteur du thermique
151!      do ig=1,ngrid
152!         lmax(ig)=lalim(ig)
153!      enddo
154!      do ig=1,ngrid
155!         do l=nlay,lalim(ig)+1,-1
156!            if (zw2(ig,l).le.1.e-10) then
157!               lmax(ig)=l-1
158!            endif
159!         enddo
160!      enddo
161!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
162
163!   
164! Determination de zw2 max
165      do ig=1,ngrid
166         wmax(ig)=0.
167      enddo
168
169      do l=1,nlay
170         do ig=1,ngrid
171            if (l.le.lmax(ig)) then
172                zw2(ig,l)=sqrt(zw2(ig,l))
173                wmax(ig)=max(wmax(ig),zw2(ig,l))
174            else
175                 zw2(ig,l)=0.
176            endif
177          enddo
178      enddo
179
180!   Longueur caracteristique correspondant a la hauteur des thermiques.
181      do  ig=1,ngrid
182         zmax(ig)=0.
183         zlevinter(ig)=zlev(ig,1)
184      enddo
185      do  ig=1,ngrid
186! calcul de zlevinter
187
188!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
189! FH A eliminer
190! Simplification
191!          zlevinter(ig)=(zlev(ig,lmax(ig)+1)-zlev(ig,lmax(ig)))*  &
192!     &    linter(ig)+zlev(ig,lmax(ig))-lmax(ig)*(zlev(ig,lmax(ig)+1)  &
193!     &    -zlev(ig,lmax(ig)))
194!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
195
196          zlevinter(ig)=zlev(ig,lmax(ig)) + &
197     &    (linter(ig)-lmax(ig))*(zlev(ig,lmax(ig)+1)-zlev(ig,lmax(ig)))
198           zmax(ig)=max(zmax(ig),zlevinter(ig)-zlev(ig,lmin(ig)))
199      enddo
200
201! Verification que lalim<=lmax
202      do ig=1,ngrid
203         if(lalim(ig)>lmax(ig)) then
204           if ( prt_level > 1 ) THEN
205            print*,'WARNING thermcell_dry ig=',ig,'  lalim=',lalim(ig),'  lmax(ig)=',lmax(ig)
206           endif
207           lmax(ig)=lalim(ig)
208         endif
209      enddo
210     
211      RETURN
212      END
Note: See TracBrowser for help on using the repository browser.