source: LMDZ.3.3/branches/rel-LF/libf/dyn3d/redecoupe.F @ 115

Last change on this file since 115 was 2, checked in by lmdz, 25 years ago

Initial revision

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 18.9 KB
Line 
1      SUBROUTINE redecoupe(irec,massemn,pbarun,pbarvn,wn,tetan,phin,
2     s     nrec,avant,airefi,
3     s     zmfu, zmfd, zen_u, zde_u,zen_d, zde_d, coefkz,
4     s     yu1,yv1,ftsol,pctsrf,
5     s     frac_impa,frac_nucl,phisn)
6
7      IMPLICIT NONE
8
9#include "dimensions.h"
10#include "paramet.h"
11
12#include "comvert.h"
13#include "comconst.h"
14#include "comgeom2.h"
15
16#include "tracstoke.h"
17
18      integer irec,nrec,i,j
19
20      integer ig,l
21      integer imo,jmo,imn,jmn,ii,jj,ig
22      parameter (imn=iim,jmn=jjm,imo=imn/2,jmo=(jmn+1)/2)
23      integer ngrido,ngridn
24      parameter(ngrido=(jmo-1)*imo+2,ngridn=(jmn-1)*imn+2)
25
26
27      INTEGER nbsrf
28      PARAMETER (nbsrf=4) ! nombre de sous-fractions pour une maille
29
30      real zmfd(ngridn,llm),zde_d(ngridn,llm),zen_d(ngridn,llm)
31      real zmfu(ngridn,llm),zde_u(ngridn,llm),zen_u(ngridn,llm)
32      real coefkz(ngridn,llm)
33      real frac_impa(ngridn,llm),frac_nucl(ngridn,llm)
34      real yu1(ngridn), yv1(ngridn)
35      real ftsol(ngridn,nbsrf),pctsrf(ngridn,nbsrf)
36      integer imfu,imfd,ien_u,ide_u,
37     s      ien_d,ide_d,
38     s      icoefkz,izu1,izv1,
39     s      itsol,ipsf,
40     s      ilei, ilec
41      parameter(imfu=1,imfd=llm+1,ien_u=2*llm+1,ide_u=3*llm+1,
42     s      ien_d=4*llm+1,ide_d=5*llm+1,
43     s      icoefkz=6*llm+1,
44     s      ilei=7*llm+1,ilec=8*llm+1,
45     s      izu1=9*llm+1,izv1=9*llm+2,
46     s      itsol=9*llm+3,ipsf=9*llm+3+nbsrf)
47      logical avant
48
49
50      real massefi(ngridn,llm)
51
52      real massemn(imn+1,jmn+1,llm),tetan(imn+1,jmn+1,llm)
53      real pbarun(imn+1,jmn+1,llm),pbarvn(imn+1,jmn,llm)
54      real wn(imn+1,jmn+1,llm),phin(imn+1,jmn+1,llm)
55      real phisn(imn+1,jmn+1)
56
57      real massemo(imo+1,jmo+1,llm),tetao(imo+1,jmo+1,llm)
58      real pbaruo(imo+1,jmo+1,llm),pbarvo(imo+1,jmo,llm)
59      real wo(imo+1,jmo+1,llm),phio(imo+1,jmo+1,llm)
60      real phiso(imo+1,jmo+1)
61
62      real pbarvst(imo+1,jmo+1,llm)
63
64      real airefi(ngridn)
65
66      real xlecn(ngridn,9*llm+2+2*nbsrf),tmpn(imn+1,jmn+1)
67      real xleco(ngrido,9*llm+2+2*nbsrf),tmpo(imo+1,jmo+1)
68
69      real zcontrole(ngridn),zmass,tmpdyn(imn+1,jmn+1),zflux
70
71      real ziadvtrac,zrec,ziadvtrac2,zrec2
72      real zim,zjm,zlm,zklon,zklev
73
74      real zpi
75c  longitudes et latitudes lues
76      real rlonul(1:imo+1),rlatvl(1:jmo)
77      real rlonvl(1:imo+1),rlatul(1:jmo+1)
78c  longitudes et latitudes anciennes
79      real rlonuo(0:imo+1),rlatvo(0:jmo+1)
80c  longitudes et latitudes nouvelles
81      real rlonun(0:imn+1),rlatvn(0:jmn+1)
82      real aireo(imo+1,jmo+1)
83
84      integer ndecx(imo+1),ndecy(jmo+1)
85      real alphax(imn+1),alphay(jmn+1)
86      real alphaxo(imo+1)
87      real alpha(imn+1,jmn+1)
88
89      real aa,uu(0:imo+1),vv(imo+1,0:jmo+1)
90
91
92      integer iest(imo+1),iouest(imo+1)
93      integer jsud(jmo+1),jnord(jmo+1)
94
95      integer in,io,jn,jo,l,iin,jjn
96      integer i,j
97      real dlatm,dlatp,dlonm,dlonp
98
99
100      zpi=2.*asin(1.)
101
102
103c==================================================================
104c   Si le numero du record est 0 alors: INITIALISATION
105c==================================================================
106c
107      print*,'ENTREE DANS LECTFLUX'
108        print*,'IREC=',IREC
109      if(irec.eq.0) then
110
111        print*,'IREC==',0
112
113C test         call inigeom
114c==================================================================
115c   Definition des surdecoupages dans les deux directions
116c==================================================================
117
118      ndecx(1)=1
119      do io=2,imo
120         ndecx(io)=2
121      enddo
122      ndecx(imo+1)=1
123
124      ndecy(1)=1
125      do jo=2,jmo
126         ndecy(jo)=2
127      enddo
128      ndecy(jmo+1)=1
129
130      ii=0
131      do io=1,imo+1
132         ii=ii+ndecx(io)
133      enddo
134      if(ii.ne.iim) then
135         print*,'ii=',ii,'   iim=',iim
136         stop
137      endif
138
139      jj=0
140      do jo=1,jmo+1
141         jj=jj+ndecy(jo)
142      enddo
143      if(jj.ne.jjp1) then
144         print*,'jj=',jj,'   jjm=',jjm
145         stop
146      endif
147
148c==================================================================
149c   Calcul des jsud,... correspondant aux intersections des
150c   grilles.
151c==================================================================
152
153      iest(1)=0
154      do io=2,imo+1
155         iest(io)=iest(io-1)+ndecx(io-1)
156         iouest(io-1)=iest(io)
157      enddo
158      iouest(imo+1)=iest(imo+1)+ndecx(imo+1)
159
160      jnord(1)=0
161      do jo=2,jmo+1
162         jnord(jo)=jnord(jo-1)+ndecy(jo-1)
163         jsud(jo-1)=jnord(jo)
164      enddo
165      jsud(jmo+1)=jnord(jmo+1)+ndecy(jmo+1)
166
167c==================================================================
168c   ouverture des fichiers, lecture des entetes
169c==================================================================
170
171c   Fichier fluxmass
172#ifdef CRAY
173         CALL ASSIGN("assign -N ieee -F null f:fluxmass")
174#endif
175      open(47,file='fluxmass',form='unformatted',
176     s     access='direct'
177     s     ,recl=4*6*(imo+1)*(jmo+1)*llm)
178      read(47,rec=1) zrec,dtvr,ziadvtrac,zim,zjm,zlm
179     s   ,rlonul,rlonvl,rlatul,rlatvl,aireo
180     s    ,phiso
181
182      print*,'zrec,dtvr,ziadvtrac,zim,zjm,zlm'
183      print*,zrec,dtvr,ziadvtrac,zim,zjm,zlm
184      if((imo-nint(zim))*(jmo-nint(zjm)).ne.0) then
185        print*,'Modifier les dimensions dans redecoupe '
186        print*,'Mettre imo=',zim,'   jmo=',zjm
187        stop
188      endif
189
190
191c  Fichier physique
192c  Fichier lessivage (supprime les donnees utiles sont dans "physique")
193#ifdef CRAY
194         CALL ASSIGN("assign -N ieee -F null f:physique")
195#endif
196      open(49,file='physique',form='unformatted',
197     s     access='direct'
198     s     ,recl=4*ngrido*(9*llm+2+2*nbsrf))
199      read(49,rec=1) zrec2,ziadvtrac2,zklon,zklev
200      print*,'Entete du fichier physique'
201      print*,zrec2,ziadvtrac2,zklon,zklev
202
203      nrec=zrec
204      istdyn=ziadvtrac
205      istphy=ziadvtrac2
206
207c==================================================================
208c   Definition des anciennes latitudes et longitudes
209c   (qui pourraient etre relues plus tard)
210c==================================================================
211
212      rlonuo(0)=-zpi
213      do io=1,imo
214c        rlonuo(io)=2.*zpi/FLOAT(imo)*(io+0.5-0.5*FLOAT(imo)-1.)
215c        print*,'LON ',io,rlonuo(io),rlonul(io)
216         rlonuo(io)=rlonul(io)
217      enddo
218      rlonuo(imo+1)=zpi
219
220      rlatvo(0)=zpi/2.
221      do jo=1,jmo
222c        rlatvo(jo)=zpi/FLOAT(jmo)*(0.5*FLOAT(jmo)+1.-jo-0.5)
223c        print*,'LAT ',jo,rlatvo(jo),rlatvl(jo)
224         rlatvo(jo)=rlatvl(jo)
225      enddo
226      rlatvo(jmo+1)=-zpi/2.
227
228c     do jo=1,jmo+1
229c        do io=1,imo+1
230c           aireo(io,jo)=rad*rad
231c    s         *(rlonuo(io)-rlonuo(io-1))
232c    s         *(sin(rlatvo(jo-1))-sin(rlatvo(jo)))
233c           aireo(io,jo)=airel(io,jo)
234c        enddo
235c        aireo(1,jo)=aireo(1,jo)+aireo(imo+1,jo)
236c        aireo(imo+1,jo)=aireo(1,jo)
237c     enddo
238
239      do io=2,imo
240         alphaxo(io)=1.
241      enddo
242      alphaxo(1)=(rlonuo(1)-rlonuo(0))
243     s        /(rlonuo(1)-rlonuo(0)+rlonuo(imo+1)-rlonuo(imo))
244      alphaxo(imo+1)=1.-alphaxo(1)
245
246c==================================================================
247c    Definition des nouvelles latitudes et longitudes
248c==================================================================
249
250      rlonun(0)=-zpi
251      do io=1,imo+1
252         do iin=1,iouest(io)-iest(io)
253            in=iin+iest(io)
254            rlonun(in)=
255     s      rlonuo(io-1)+iin*(rlonuo(io)-rlonuo(io-1))
256     s      /ndecx(io)
257            alphax(in)=alphaxo(io)/ndecx(io)
258            print787,io,rlonuo(io-1)*180./zpi,in
259     s      ,iest(io),iouest(io),rlonun(in)*180./zpi,alphax(in)
260         enddo
261      enddo
262
263      rlatvn(0)=0.5*zpi
264      do jo=1,jmo+1
265         print*,'jo=',jo
266         do jjn=1,jsud(jo)-jnord(jo)
267            jn=jnord(jo)+jjn
268            rlatvn(jn)=rlatvo(jo-1)+jjn*(rlatvo(jo)-rlatvo(jo-1))
269     s      /ndecy(jo)
270            alphay(jn)=(sin(rlatvn(jn-1))-sin(rlatvn(jn)))
271     s                /(sin(rlatvo(jo-1))-sin(rlatvo(jo)))
272            print787,jo,rlatvo(jo-1)*180./zpi,jn
273     s      ,jnord(jo),jsud(jo),rlatvn(jn)*180./zpi,alphay(jn)
274         enddo
275      enddo
276
277787   format(i5,f10.2,3(i5),2(f10.2))
278      do in=1,imn
279         rlonu(in)=rlonun(in)
280         rlonv(in)=0.5*(rlonun(in)+rlonun(in-1))
281      enddo
282      rlonv(imn+1)=rlonv(1)+2.*zpi
283      rlonu(imn+1)=rlonu(1)+2.*zpi
284
285      do jn=1,jmn
286         rlatv(jn)=rlatvn(jn)
287      enddo
288      do jn=1,jmn+1
289         rlatu(jn)=0.5*(rlatvn(jn-1)+rlatvn(jn))
290      enddo
291
292      do jn=1,jmn+1
293         do in=1,imn
294            alpha(in,jn)=alphax(in)*alphay(jn)
295         enddo
296         alpha(imn+1,jn)=0.
297      enddo
298
299c     call dump2d(iip1,jjp1,alpha,'ALPHA   ')
300
301c      .  on a :  cu(i,j) = rad * COS(y) * dx/dX         .
302c      .          cv( j ) = rad          * dy/dY         .
303c   A 1 point scalaire P (i,j) de la grille, reguliere en (X,Y) , sont
304c   affectees 4 aires entourant P , calculees respectivement aux points
305c            ( i + 1/4, j - 1/4 )    :    aireij1 (i,j)
306c            ( i + 1/4, j + 1/4 )    :    aireij2 (i,j)
307c            ( i - 1/4, j + 1/4 )    :    aireij3 (i,j)
308c            ( i - 1/4, j - 1/4 )    :    aireij4 (i,j)
309c
310c                             . V
311c
312c                 aireij4 .        . aireij1
313c
314c                   U .       . P      . U
315c
316c                 aireij3 .        . aireij2
317c
318c                             . V
319
320
321      do j=1,jjp1
322         do i=1,iim
323            dlonp=rlonun(i)-rlonv(i)
324            dlonm=rlonv(i)-rlonun(i-1)
325            dlatp=sin(rlatvn(j-1))-sin(rlatu(j))
326            dlatm=sin(rlatu(j))-sin(rlatvn(j))
327            aireij1 ( i,j ) = rad*rad*dlatp*dlonp
328            aireij2 ( i,j ) = rad*rad*dlatm*dlonp
329            aireij3 ( i,j ) = rad*rad*dlatm*dlonm
330            aireij4 ( i,j ) = rad*rad*dlatp*dlonm
331      aire    ( i,j )  = aireij1(i,j) + aireij2(i,j) + aireij3(i,j) +
332     *                          aireij4(i,j)
333      alpha1  ( i,j )  = aireij1(i,j) / aire(i,j)
334      alpha2  ( i,j )  = aireij2(i,j) / aire(i,j)
335      alpha3  ( i,j )  = aireij3(i,j) / aire(i,j)
336      alpha4  ( i,j )  = aireij4(i,j) / aire(i,j)
337      alpha1p2( i,j )  = alpha1 (i,j) + alpha2 (i,j)
338      alpha1p4( i,j )  = alpha1 (i,j) + alpha4 (i,j)
339      alpha2p3( i,j )  = alpha2 (i,j) + alpha3 (i,j)
340      alpha3p4( i,j )  = alpha3 (i,j) + alpha4 (i,j)
341           enddo
342           aireij1(iip1,j)=aireij1(1,j)
343           aireij2(iip1,j)=aireij2(1,j)
344           aireij3(iip1,j)=aireij3(1,j)
345           aireij4(iip1,j)=aireij4(1,j)
346           aire(iip1,j)=aire(1,j)
347           alpha1(iip1,j)=alpha1(1,j)
348           alpha2(iip1,j)=alpha2(1,j)
349           alpha3(iip1,j)=alpha3(1,j)
350           alpha4(iip1,j)=alpha4(1,j)
351           alpha1p2(iip1,j)=alpha1p2(1,j)
352           alpha1p4(iip1,j)=alpha1p4(1,j)
353           alpha2p3(iip1,j)=alpha2p3(1,j)
354           alpha3p4(iip1,j)=alpha3p4(1,j)
355       enddo
356c     call dump2d(iip1,jjp1,aire,'AIRE   ')
357
358c     do jn=1,jjp1
359c        do in=1,iim
360c           aire(in,jn)=rad*rad*(sin(rlatvn(jn-1))-sin(rlatvn(jn)))
361c    s      *(rlonun(in)-rlonun(in-1))
362c           unsaire(in,jn)=1./aire(in,jn)
363c        enddo
364c        aire(iip1,jn)=aire(1,jn)
365c        unsaire(iip1,jn)=unsaire(1,jn)
366c     enddo
367c     call dump2d(iip1,jjp1,aire,'AIRE2   ')
368      DO 42 j = 1,jjp1
369      DO 41 i = 1,iim
370      unsaire(i,j) = 1./ aire(i,j)
371      aireu  (i,j) = aireij1(i,j) + aireij2(i,j) + aireij4(i+1,j) +
372     *                           aireij3(i+1,j)
373  41  CONTINUE
374      aireu  (iip1,j) = aireu  (1,j)
375      unsaire(iip1,j) = unsaire(1,j)
376  42  CONTINUE
377      DO 48 j = 1,jjm
378        DO i=1,iim
379         airev(i,j)     = aireij2(i,j)+ aireij3(i,j)+ aireij1(i,j+1) +
380     *                           aireij4(i,j+1)
381        ENDDO
382       airev   (iip1,j) = airev(1,j)
383  48  CONTINUE
384      apoln=0.
385      apols=0.
386      do i=1,iim
387         apoln=apoln+aire(i,1)
388         apols=apols+aire(i,jjp1)
389      enddo
390
391
392
393      do jn=1,jjp1
394         do in=1,iim
395            cu(in,jn)=rad*cos(rlatu(jn))*(rlonv(in+1)-rlonv(in))
396         enddo
397         cu(iip1,jn)=cu(1,jn)
398      enddo
399      do jn=1,jjm
400         do in=1,iim+1
401            cv(in,jn)=rad*(rlatu(jn)-rlatu(jn+1))
402         enddo
403      enddo
404
405
406c==================================================================
407c   Fin des initialisations
408      else ! irec=0
409c==================================================================
410
411
412c-----------------------------------------------------------------------
413c   Lecture des fichiers fluxmass et  physique:
414c   -----------------------------------------------------
415
416c==================================================================
417c  Variables dynamiques
418c==================================================================
419
420         read(47,rec=irec) massemo,pbaruo,pbarvst,wo,tetao,phio
421        do l=1,llm
422           do j=1,jmo
423              do i=1,imo+1
424                 pbarvo(i,j,l)=pbarvst(i,j,l)
425              enddo
426           enddo
427        enddo
428
429         do l=1,llm
430            do jo=1,jmo+1
431               do io=1,imo+1
432                  do jn=jnord(jo)+1,jsud(jo)
433                     do in=iest(io)+1,iouest(io)
434                        wn(in,jn,l)=alpha(in,jn)*wo(io,jo,l)
435                        massemn(in,jn,l)=alpha(in,jn)
436     s                     *massemo(io,jo,l)
437                        tetan(in,jn,l)=tetao(io,jo,l)
438                        phin(in,jn,l)=phio(io,jo,l)
439                     enddo
440                  enddo
441               enddo
442            enddo
443            do jn=1,jmn+1
444               wn(imn+1,jn,l)=wn(1,jn,l)
445               massemn(imn+1,jn,l)=massemn(1,jn,l)
446               tetan(imn+1,jn,l)=tetan(1,jn,l)
447               phin(imn+1,jn,l)=phin(1,jn,l)
448            enddo
449         enddo
450
451         do l=1,llm
452            do jo=1,jmo+1
453               uu(imo+1)=0.5*(pbaruo(imo,jo,l)+pbaruo(imo+1,jo,l))
454               uu(0)=uu(imo+1)
455               do io=1,imo
456                  uu(io)=pbaruo(io,jo,l)
457               enddo
458               do io=1,imo+1
459                  do jn=jnord(jo)+1,jsud(jo)
460                     aa=0.
461                     do in=iest(io)+1,iouest(io)
462                        aa=aa+alphax(in)
463                        pbarun(in,jn,l)=alphay(jn)*
464     s                    (uu(io-1)+aa*(uu(io)-uu(io-1)))
465                     enddo
466                  enddo
467               enddo
468            enddo
469            do jn=1,jmn+1
470               pbarun(imn+1,jn,l)=pbarun(1,jn,l)
471            enddo
472         enddo
473
474         do l=1,llm
475            do jo=1,jmo
476               do io=1,imo+1
477                  vv(io,jo)=pbarvo(io,jo,l)
478               enddo
479            enddo
480            do io=1,imo+1
481               vv(io,0)=0.
482               vv(io,jmo+1)=0.
483            enddo
484            do jo=1,jmo+1
485               do io=1,imo+1
486                  aa=0.
487c                 do jn=jnord(jo)+1,max(jsud(jo),jmo)
488                  do jn=jnord(jo)+1,min(jsud(jo),jmn)
489                     aa=aa+alphay(jn)
490                     do in=iest(io)+1,iouest(io)
491                        pbarvn(in,jn,l)=alphax(in)*
492     s                  (vv(io,jo-1)+aa*(vv(io,jo)-vv(io,jo-1)))
493                     enddo
494                  enddo
495               enddo
496            enddo
497            do jn=1,jmn
498               pbarvn(iip1,jn,l)=pbarvn(1,jn,l)
499            enddo
500         enddo
501
502c==================================================================
503c  Variables physiques
504c==================================================================
505         read(49,rec=irec) ((xleco(ig,l),ig=1,ngrido),
506     s                                    l=1,9*llm+2+2*nbsrf)
507
508c==================================================================
509c  Passage  a la nouvelle grille
510c==================================================================
511         do l=1,9*llm+2+2*nbsrf
512c   passage aa la grille dynamique ancienne
513            do io=1,imo+1
514               tmpo(io,1)=xleco(1,l)
515               tmpo(io,jmo+1)=xleco(ngrido,l)
516            enddo
517            do jo=2,jmo
518               do io=1,imo
519                  tmpo(io,jo)=xleco((jo-2)*imo+io+1,l)
520               enddo
521               tmpo(imo+1,jo)=tmpo(1,jo)
522            enddo
523c   passage a la grillle dynamique nouvelle
524            do jo=1,jmo+1
525               do io=1,imo+1
526                  do jn=jnord(jo)+1,jsud(jo)
527                     do in=iest(io)+1,iouest(io)
528                        tmpn(in,jn)=tmpo(io,jo)
529                     enddo
530                  enddo
531               enddo
532            enddo
533c   passage a la grille physique nouvelle
534            xlecn(1,l)=tmpn(1,1)
535            xlecn(ngridn,l)=tmpn(1,jmn+1)
536            do jn=2,jmn
537               do in=1,imn
538                  xlecn((jn-2)*imn+in+1,l)=tmpn(in,jn)
539               enddo
540            enddo
541         enddo
542
543c==================================================================
544
545       do l=1,llm
546          do ig=1,ngridn
547             coefkz(ig,l)=xlecn(ig,icoefkz+l-1)
548             frac_impa(ig,l)=xlecn(ig,ilei+l-1)
549             frac_nucl(ig,l)=xlecn(ig,ilec+l-1)
550          enddo
551       enddo
552       do l=1,nbsrf
553          do ig=1,ngridn
554             ftsol(ig,l)=xlecn(ig,itsol+l-1)
555             pctsrf(ig,l)=xlecn(ig,ipsf+l-1)
556          enddo
557       enddo
558       do ig=1,ngridn
559          yv1(ig)=xlecn(ig,izv1)
560          yu1(ig)=xlecn(ig,izu1)
561       enddo
562C
563      if(avant) then
564c   Simu directe
565       do l=1,llm
566          do ig=1,ngridn
567             zmfu(ig,l)=xlecn(ig,imfu+l-1)
568             zmfd(ig,l)=xlecn(ig,imfd+l-1)
569             zde_u(ig,l)=xlecn(ig,ide_u+l-1)
570             zen_u(ig,l)=xlecn(ig,ien_u+l-1)
571             zde_d(ig,l)=xlecn(ig,ide_d+l-1)
572             zen_d(ig,l)=xlecn(ig,ien_d+l-1)
573          enddo
574       enddo
575      else
576c   Simu retro
577       do l=1,llm
578          do ig=1,ngridn
579             zmfd(ig,l)=-xlecn(ig,imfu+l-1)
580             zmfu(ig,l)=-xlecn(ig,imfd+l-1)
581             zen_d(ig,l)=xlecn(ig,ide_u+l-1)
582             zde_d(ig,l)=xlecn(ig,ien_u+l-1)
583             zen_u(ig,l)=xlecn(ig,ide_d+l-1)
584             zde_u(ig,l)=xlecn(ig,ien_d+l-1)
585          enddo
586       enddo
587      endif
588
589c-----------------------------------------------------------------------
590c   PETIT CONTROLE SUR LES FLUX CONVECTIFS...
591c-----------------------------------------------------------------------
592
593         call gr_dyn_fi(llm,iip1,jjp1,ngridn,massemn,massefi)
594
595      print*,'Ap redec irec'
596         do ig=1,ngridn
597            zcontrole(ig)=1.
598         enddo
599c   zmass=(max(massemn(ig,l),massemn(ig,l-1))/airefi(ig)
600         do l=2,llm
601            do ig=1,ngridn
602               zmass=max(massefi(ig,l),massefi(ig,l-1))/airefi(ig)
603               zflux=max(abs(zmfu(ig,l)),abs(zmfd(ig,l)))*dtphys
604               if(zflux.gt.0.9*zmass) then
605                 zcontrole(ig)=min(zcontrole(ig),0.9*zmass/zflux)
606               endif
607            enddo
608         enddo
609
610         do ig=1,ngridn
611            if(zcontrole(ig).lt.0.99999) then
612               print*,'ATTENTION !!! on reduit les flux de masse '
613               print*,'convectifs au point ig=',ig
614            endif
615         enddo
616
617         call gr_fi_dyn(1,ngridn,iip1,jjp1,zcontrole,tmpdyn)
618
619         do l=1,llm
620            do ig=1,ngridn
621               zmfu(ig,l)=zmfu(ig,l)*zcontrole(ig)
622               zmfd(ig,l)=zmfd(ig,l)*zcontrole(ig)
623               zen_u(ig,l)=zen_u(ig,l)*zcontrole(ig)
624               zde_u(ig,l)=zde_u(ig,l)*zcontrole(ig)
625               zen_d(ig,l)=zen_d(ig,l)*zcontrole(ig)
626               zde_d(ig,l)=zde_d(ig,l)*zcontrole(ig)
627            enddo
628         enddo
629
630
631      endif ! irec=0
632
633
634      RETURN
635      END
636
637
Note: See TracBrowser for help on using the repository browser.