source: LMDZ4/branches/LMDZ4V5.0-dev/libf/phylmd/cv3p1_closure.F @ 3603

Last change on this file since 3603 was 1322, checked in by Laurent Fairhead, 15 years ago

Improvements concerning wake parametrisation (from JYG, NR, IT, with more to come).
Alp_offset is read in form physiq.def file


Améliorations à la paramétrisation des poches froides (de JYG, NR, IT, d'autres
sont à venir)
Alp_offset est rajouté à la liste des paramètres lus dans physiq.def

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 17.9 KB
RevLine 
[1299]1!
2! $Id: cv3p1_closure.F 1322 2010-03-12 10:54:11Z musat $
3!
[879]4      SUBROUTINE cv3p1_closure(nloc,ncum,nd,icb,inb
5     :                      ,pbase,plcl,p,ph,tv,tvp,buoy
6     :                      ,Supmax,ok_inhib,Ale,Alp
7     o                      ,sig,w0,ptop2,cape,cin,m,iflag,coef
8     :                      ,Plim1,Plim2,asupmax,supmax0
9     :                      ,asupmaxmin,cbmf)
10
11*
12***************************************************************
13*                                                             *
14* CV3P1_CLOSURE                                               *
15*                  Ale & Alp Closure of Convect3              *
16*                                                             *
17* written by   :   Kerry Emanuel                              *
18* vectorization:   S. Bony                                    *
19* modified by :    Jean-Yves Grandpeix, 18/06/2003, 19.32.10  *
20*                  Julie Frohwirth,     14/10/2005  17.44.22  *
21***************************************************************
22*
23      implicit none
24
25#include "cvthermo.h"
26#include "cv3param.h"
27#include "YOMCST2.h"
28#include "YOMCST.h"
29#include "conema3.h"
[973]30#include "iniprint.h"
[879]31
32c input:
33      integer ncum, nd, nloc
34      integer icb(nloc), inb(nloc)
35      real pbase(nloc),plcl(nloc)
36      real p(nloc,nd), ph(nloc,nd+1)
37      real tv(nloc,nd),tvp(nloc,nd), buoy(nloc,nd)
38      real Supmax(nloc,nd)
39      logical ok_inhib ! enable convection inhibition by dryness
40      real Ale(nloc),Alp(nloc)
41
42c input/output:
43      real sig(nloc,nd), w0(nloc,nd), ptop2(nloc)
44
45c output:
46      real cape(nloc),cin(nloc)
47      real m(nloc,nd)
48      real Plim1(nloc),Plim2(nloc)
49      real asupmax(nloc,nd),supmax0(nloc)
50      real asupmaxmin(nloc)
51      integer iflag(nloc)
52c
53c local variables:
[973]54      integer il, i, j, k, icbmax, i0(nloc)
[879]55      real deltap, fac, w, amu
56      real rhodp
57      real Pbmxup
58      real dtmin(nloc,nd), sigold(nloc,nd)
59      real coefmix(nloc,nd)
60      real pzero(nloc),ptop2old(nloc)
61      real cina(nloc),cinb(nloc)
62      integer ibeg(nloc)
63      integer nsupmax(nloc)
64      real supcrit,temp(nloc,nd)
[1322]65      real P1(nloc),Pmin(nloc),plfc(nloc)
[879]66      real asupmax0(nloc)
67      logical ok(nloc)
68      real siglim(nloc,nd),wlim(nloc,nd),mlim(nloc,nd)
69      real wb2(nloc)
70      real cbmflim(nloc),cbmf1(nloc),cbmfmax(nloc),cbmf(nloc)
71      real cbmflast(nloc)
72      real coef(nloc)
73      real xp(nloc),xq(nloc),xr(nloc),discr(nloc),b3(nloc),b4(nloc)
74      real theta(nloc),bb(nloc)
75      real term1,term2,term3
76      real alp2(nloc) ! Alp with offset
77      real wb,sigmax
78      data wb /2./, sigmax /0.1/
[1299]79
80      CHARACTER (LEN=20) :: modname='cv3p1_closure'
81      CHARACTER (LEN=80) :: abort_message
[879]82c
83c      print *,' -> cv3p1_closure, Ale ',ale(1)
84c
85
86c -------------------------------------------------------
87c -- Initialization
88c -------------------------------------------------------
89
90c
91c
92      do il = 1,ncum
93       alp2(il) = max(alp(il),1.e-5)
[973]94cIM
95       alp2(il) = max(alp(il),1.e-12)
[879]96      enddo
97c
98      PBMXUP=50.    ! PBMXUP+PBCRIT = cloud depth above which mixed updraughts
99c                     exist (if any)
100
[973]101       if(prt_level.GE.20)
102     . print*,'cv3p1_param nloc ncum nd icb inb nl',nloc,ncum,nd,
103     . icb(nloc),inb(nloc),nl
[879]104      do k=1,nl
105       do il=1,ncum
106        m(il,k)=0.0
107       enddo
108      enddo
109
110c -------------------------------------------------------
111c -- Reset sig(i) and w0(i) for i>inb and i<icb
112c -------------------------------------------------------
113
114c update sig and w0 above LNB:
115
116      do 100 k=1,nl-1
117       do 110 il=1,ncum
118        if ((inb(il).lt.(nl-1)).and.(k.ge.(inb(il)+1)))then
119         sig(il,k)=beta*sig(il,k)
120     :            +2.*alpha*buoy(il,inb(il))*ABS(buoy(il,inb(il)))
121         sig(il,k)=AMAX1(sig(il,k),0.0)
122         w0(il,k)=beta*w0(il,k)
123        endif
124 110   continue
125 100  continue
126
[973]127c      if(prt.level.GE.20) print*,'cv3p1_param apres 100'
[879]128c compute icbmax:
129
130      icbmax=2
131      do 200 il=1,ncum
132        icbmax=MAX(icbmax,icb(il))
133 200  continue
[973]134!     if(prt.level.GE.20) print*,'cv3p1_param apres 200'
[879]135
136c update sig and w0 below cloud base:
137
138      do 300 k=1,icbmax
139       do 310 il=1,ncum
140        if (k.le.icb(il))then
141         sig(il,k)=beta*sig(il,k)-2.*alpha*buoy(il,icb(il))
142     $                                    *buoy(il,icb(il))
143         sig(il,k)=amax1(sig(il,k),0.0)
144         w0(il,k)=beta*w0(il,k)
145        endif
146310    continue
147300    continue
[973]148       if(prt_level.GE.20) print*,'cv3p1_param apres 300'
[879]149c -------------------------------------------------------------
150c -- Reset fractional areas of updrafts and w0 at initial time
151c -- and after 10 time steps of no convection
152c -------------------------------------------------------------
153
154      do 400 k=1,nl-1
155       do 410 il=1,ncum
156        if (sig(il,nd).lt.1.5.or.sig(il,nd).gt.12.0)then
157         sig(il,k)=0.0
158         w0(il,k)=0.0
159        endif
160 410   continue
161 400  continue
[973]162      if(prt_level.GE.20) print*,'cv3p1_param apres 400'
[879]163c
164c -------------------------------------------------------------
165Cjyg1
166C --  Calculate adiabatic ascent top pressure (ptop)
167c -------------------------------------------------------------
168C
169c
170cc 1. Start at first level where precipitations form
171      do il = 1,ncum
172        Pzero(il) = Plcl(il)-PBcrit
173      enddo
174c
175cc 2. Add offset
176      do il = 1,ncum
177        Pzero(il) = Pzero(il)-PBmxup
178      enddo
179      do il=1,ncum
180         ptop2old(il)=ptop2(il)
181      enddo
182c
183      do il = 1,ncum
184cCR:c est quoi ce 300??
185        P1(il) = Pzero(il)-300.
186      enddo
187
188c    compute asupmax=abs(supmax) up to lnm+1
189
190      DO il=1,ncum
191        ok(il)=.true.
192        nsupmax(il)=inb(il)
193      ENDDO
194
195      DO i = 1,nl
196        DO il = 1,ncum
197        IF (i .GT. icb(il) .AND. i .LE. inb(il)) THEN
198        IF (P(il,i) .LE. Pzero(il) .and.
199     $       supmax(il,i) .lt. 0 .and. ok(il)) THEN
200           nsupmax(il)=i
201           ok(il)=.false.
202        ENDIF    ! end IF (P(i) ...
203        ENDIF    ! end IF (icb+1 le i le inb)
204        ENDDO
205      ENDDO
206
[973]207      if(prt_level.GE.20) print*,'cv3p1_param apres 2.'
[879]208      DO i = 1,nl
209        DO il = 1,ncum
210        asupmax(il,i)=abs(supmax(il,i))
211        ENDDO
212      ENDDO
213
214c
215        DO il = 1,ncum
216        asupmaxmin(il)=10.
217        Pmin(il)=100.
[973]218!IM ??
219        asupmax0(il)=0.
[879]220        ENDDO
221
222cc 3.  Compute in which level is Pzero
223
[973]224cIM bug      i0 = 18
225       DO il = 1,ncum
226        i0(il) = nl
227       ENDDO
[879]228
229       DO i = 1,nl
230        DO il = 1,ncum
231         IF (i .GT. icb(il) .AND. i .LE. inb(il)) THEN
232           IF (P(il,i) .LE. Pzero(il) .AND. P(il,i) .GE. P1(il)) THEN
233            IF (Pzero(il) .GT. P(il,i) .AND.
234     $           Pzero(il) .LT. P(il,i-1)) THEN
[973]235             i0(il) = i
[879]236            ENDIF
237           ENDIF
238          ENDIF
239        ENDDO
240       ENDDO
[973]241       if(prt_level.GE.20) print*,'cv3p1_param apres 3.'
242
[879]243cc 4.  Compute asupmax at Pzero
244
245       DO i = 1,nl
246        DO il = 1,ncum
247         IF (i .GT. icb(il) .AND. i .LE. inb(il)) THEN
248           IF (P(il,i) .LE. Pzero(il) .AND. P(il,i) .GE. P1(il)) THEN
[973]249             asupmax0(il) =
250     $             ((Pzero(il)-P(il,i0(il)-1))*asupmax(il,i0(il))
251     $             -(Pzero(il)-P(il,i0(il)))*asupmax(il,i0(il)-1))
252     $             /(P(il,i0(il))-P(il,i0(il)-1))
[879]253           ENDIF
254         ENDIF
255        ENDDO
256       ENDDO
257
258
259      DO i = 1,nl
260        DO il = 1,ncum
261         IF (P(il,i) .EQ. Pzero(il)) THEN
262           asupmax(i,il) = asupmax0(il)
263         ENDIF
264        ENDDO
265      ENDDO
[973]266      if(prt_level.GE.20) print*,'cv3p1_param apres 4.'
[879]267
268cc 5. Compute asupmaxmin, minimum of asupmax
269
270      DO i = 1,nl
271        DO il = 1,ncum
272        IF (i .GT. icb(il) .AND. i .LE. inb(il)) THEN
273        IF (P(il,i) .LE. Pzero(il) .AND. P(il,i) .GE. P1(il)) THEN
274          IF (asupmax(il,i) .LT. asupmaxmin(il)) THEN
275            asupmaxmin(il)=asupmax(il,i)
276            Pmin(il)=P(il,i)
277          ENDIF
278        ENDIF
279        ENDIF
280        ENDDO
281      ENDDO
282
283      DO il = 1,ncum
[973]284!IM
285        if(prt_level.GE.20) THEN
286         print*,'cv3p1_closure il asupmax0 asupmaxmin',il,asupmax0(il),
287     $ asupmaxmin(il) ,Pzero(il),Pmin(il)
288        endif
[879]289          IF (asupmax0(il) .LT. asupmaxmin(il)) THEN
290             asupmaxmin(il) = asupmax0(il)
291             Pmin(il) = Pzero(il)
292          ENDIF
293      ENDDO
[973]294      if(prt_level.GE.20) print*,'cv3p1_param apres 5.'
[879]295
296c
297c   Compute Supmax at Pzero
298c
299      DO i = 1,nl
300        DO il = 1,ncum
301        IF (i .GT. icb(il) .AND. i .LE. inb(il)) THEN
302        IF (P(il,i) .LE. Pzero(il)) THEN
303         Supmax0(il) = ((P(il,i  )-Pzero(il))*aSupmax(il,i-1)
304     $             -(P(il,i-1)-Pzero(il))*aSupmax(il,i  ))
305     $             /(P(il,i)-P(il,i-1))
306         GO TO 425
307        ENDIF    ! end IF (P(i) ...
308        ENDIF    ! end IF (icb+1 le i le inb)
309        ENDDO
310      ENDDO
311
312425   continue
[973]313      if(prt_level.GE.20) print*,'cv3p1_param apres 425.'
[879]314
315cc 6. Calculate ptop2
316c
317      DO il = 1,ncum
318        IF (asupmaxmin(il) .LT. Supcrit1) THEN
319          Ptop2(il) = Pmin(il)
320        ENDIF
321
322        IF (asupmaxmin(il) .GT. Supcrit1
323     $ .AND. asupmaxmin(il) .LT. Supcrit2) THEN
324          Ptop2(il) = Ptop2old(il)
325        ENDIF
326
327        IF (asupmaxmin(il) .GT. Supcrit2) THEN
328            Ptop2(il) =  Ph(il,inb(il))
329        ENDIF
330      ENDDO
331c
[973]332      if(prt_level.GE.20) print*,'cv3p1_param apres 6.'
[879]333
334cc 7. Compute multiplying factor for adiabatic updraught mass flux
335c
336c
337      IF (ok_inhib) THEN
338c
339      DO i = 1,nl
340        DO il = 1,ncum
341         IF (i .le. nl) THEN
342         coefmix(il,i) = (min(ptop2(il),ph(il,i))-ph(il,i))
343     $                  /(ph(il,i+1)-ph(il,i))
344         coefmix(il,i) = min(coefmix(il,i),1.)
345         ENDIF
346        ENDDO
347      ENDDO
348c
349c
350      ELSE   ! when inhibition is not taken into account, coefmix=1
351c
352
353c
354      DO i = 1,nl
355        DO il = 1,ncum
356         IF (i .le. nl) THEN
357         coefmix(il,i) = 1.
358         ENDIF
359        ENDDO
360      ENDDO
361c
362      ENDIF  ! ok_inhib
[973]363      if(prt_level.GE.20) print*,'cv3p1_param apres 7.'
[879]364c -------------------------------------------------------------------
365c -------------------------------------------------------------------
366c
367
368Cjyg2
369C
370c==========================================================================
371C
372c
373c -------------------------------------------------------------
374c -- Calculate convective inhibition (CIN)
375c -------------------------------------------------------------
376
377c      do i=1,nloc
378c      print*,'avant cine p',pbase(i),plcl(i)
379c      enddo
[973]380c     do j=1,nd
381c     do i=1,nloc
[879]382c      print*,'avant cine t',tv(i),tvp(i)
[973]383c     enddo
384c     enddo
[879]385      CALL cv3_cine (nloc,ncum,nd,icb,inb
386     :                      ,pbase,plcl,p,ph,tv,tvp
[1322]387     :                      ,cina,cinb,plfc)
[879]388c
389      DO il = 1,ncum
390        cin(il) = cina(il)+cinb(il)
391      ENDDO
[973]392      if(prt_level.GE.20) print*,'cv3p1_param apres cv3_cine'
[879]393c -------------------------------------------------------------
394c --Update buoyancies to account for Ale
395c -------------------------------------------------------------
396c
397      CALL cv3_buoy (nloc,ncum,nd,icb,inb
398     :                      ,pbase,plcl,p,ph,Ale,Cin
399     :                      ,tv,tvp
400     :                      ,buoy )
[973]401      if(prt_level.GE.20) print*,'cv3p1_param apres cv3_buoy'
[879]402
403c -------------------------------------------------------------
404c -- Calculate convective available potential energy (cape),
405c -- vertical velocity (w), fractional area covered by
406c -- undilute updraft (sig), and updraft mass flux (m)
407c -------------------------------------------------------------
408
409      do 500 il=1,ncum
410       cape(il)=0.0
411 500  continue
412
413c compute dtmin (minimum buoyancy between ICB and given level k):
414
415      do k=1,nl
416       do il=1,ncum
417         dtmin(il,k)=100.0
418       enddo
419      enddo
420
421      do 550 k=1,nl
422       do 560 j=minorig,nl
423        do 570 il=1,ncum
424          if ( (k.ge.(icb(il)+1)).and.(k.le.inb(il)).and.
425     :         (j.ge.icb(il)).and.(j.le.(k-1)) )then
426           dtmin(il,k)=AMIN1(dtmin(il,k),buoy(il,j))
427          endif
428 570     continue
429 560   continue
430 550  continue
431
432c the interval on which cape is computed starts at pbase :
433
434      do 600 k=1,nl
435       do 610 il=1,ncum
436
437        if ((k.ge.(icb(il)+1)).and.(k.le.inb(il))) then
438
439         deltap = MIN(pbase(il),ph(il,k-1))-MIN(pbase(il),ph(il,k))
440         cape(il)=cape(il)+rrd*buoy(il,k-1)*deltap/p(il,k-1)
441         cape(il)=AMAX1(0.0,cape(il))
442         sigold(il,k)=sig(il,k)
443
444
445cjyg       Coefficient coefmix limits convection to levels where a sufficient
446c          fraction of mixed draughts are ascending.
447         siglim(il,k)=coefmix(il,k)*alpha1*dtmin(il,k)*ABS(dtmin(il,k))
448         siglim(il,k)=amax1(siglim(il,k),0.0)
449         siglim(il,k)=amin1(siglim(il,k),0.01)
450cc         fac=AMIN1(((dtcrit-dtmin(il,k))/dtcrit),1.0)
451         fac = 1.
452         wlim(il,k)=fac*SQRT(cape(il))
453         amu=siglim(il,k)*wlim(il,k)
454         rhodp = 0.007*p(il,k)*(ph(il,k)-ph(il,k+1))/tv(il,k)
455         mlim(il,k)=amu*rhodp
456c         print*, 'siglim ', k,siglim(1,k)
457        endif
458
459 610   continue
460 600  continue
[973]461      if(prt_level.GE.20) print*,'cv3p1_param apres 600'
[879]462
463      do 700 il=1,ncum
[973]464!IM beg
465        if(prt_level.GE.20) THEN
466         print*,'cv3p1_closure il icb mlim ph ph+1 ph+2',il,
467     $icb(il),mlim(il,icb(il)+1),ph(il,icb(il)),
468     $ph(il,icb(il)+1),ph(il,icb(il)+2)
469        endif
470
471        if (icb(il)+1.le.inb(il)) then
472!IM end
[879]473       mlim(il,icb(il))=0.5*mlim(il,icb(il)+1)
474     :             *(ph(il,icb(il))-ph(il,icb(il)+1))
475     :             /(ph(il,icb(il)+1)-ph(il,icb(il)+2))
[973]476!IM beg
477        endif !(icb(il.le.inb(il))) then
478!IM end
[879]479 700  continue
[973]480      if(prt_level.GE.20) print*,'cv3p1_param apres 700'
[879]481
482cjyg1
483c------------------------------------------------------------------------
484cc     Correct mass fluxes so that power used to overcome CIN does not
485cc     exceed Power Available for Lifting (PAL).
486c------------------------------------------------------------------------
487c
488      do il = 1,ncum
489       cbmflim(il) = 0.
490       cbmf(il) = 0.
491      enddo
492c
493cc 1. Compute cloud base mass flux of elementary system (Cbmf0=Cbmflim)
494c
495      do k= 1,nl
496       do il = 1,ncum
[1322]497!old       IF (k .ge. icb(il) .and. k .le. inb(il)) THEN
498!IM        IF (k .ge. icb(il)+1 .and. k .le. inb(il)) THEN
499       IF (k .ge. icb(il) .and. k .le. inb(il)         !cor jyg
500     $     .and. icb(il)+1 .le. inb(il)) THEN          !cor jyg
[879]501         cbmflim(il) = cbmflim(il)+MLIM(il,k)
502        ENDIF
503       enddo
504      enddo
[973]505      if(prt_level.GE.20) print*,'cv3p1_param apres cbmflim'
506
[879]507cc 1.5 Compute cloud base mass flux given by Alp closure (Cbmf1), maximum
508cc     allowed mass flux (Cbmfmax) and final target mass flux (Cbmf)
509cc     Cbmf is set to zero if Cbmflim (the mass flux of elementary cloud) is
510c--    exceedingly small.
511c
512      DO il = 1,ncum
513        wb2(il) = sqrt(2.*max(Ale(il)+cin(il),0.))
514      ENDDO
515c
516      DO il = 1,ncum
517       cbmf1(il) = alp2(il)/(0.5*wb*wb-Cin(il))
[973]518       if(cbmf1(il).EQ.0.AND.alp2(il).NE.0.) THEN
[1299]519        write(lunout,*)
520     &  'cv3p1_closure cbmf1=0 and alp NE 0 il alp2 alp cin ',il,
[973]521     . alp2(il),alp(il),cin(il)
[1299]522        abort_message = ''
523        CALL abort_gcm (modname,abort_message,1)
[973]524       endif
[879]525       cbmfmax(il) = sigmax*wb2(il)*100.*p(il,icb(il))
526     :              /(rrd*tv(il,icb(il)))
527      ENDDO
528c
529      DO il = 1,ncum
530       IF (cbmflim(il) .gt. 1.e-6) THEN
531cATTENTION TEST CR
532c         if (cbmfmax(il).lt.1.e-12) then
533        cbmf(il) = min(cbmf1(il),cbmfmax(il))
534c         else
535c         cbmf(il) = cbmf1(il)
536c         endif
537c        print*,'cbmf',cbmf1(il),cbmfmax(il)
538       ENDIF
539      ENDDO
[973]540      if(prt_level.GE.20) print*,'cv3p1_param apres cbmflim_testCR'
[879]541c
542cc 2. Compute coefficient and apply correction
543c
544      do il = 1,ncum
545       coef(il) = (cbmf(il)+1.e-10)/(cbmflim(il)+1.e-10)
546      enddo
[973]547      if(prt_level.GE.20) print*,'cv3p1_param apres coef_plantePLUS'
[879]548c
549      DO k = 1,nl
550        do il = 1,ncum
551         IF ( k .ge. icb(il)+1 .AND. k .le. inb(il)) THEN
[1322]552         amu=beta*sig(il,k)*w0(il,k)+
553     :   (1.-beta)*coef(il)*siglim(il,k)*wlim(il,k)
554         w0(il,k) = wlim(il,k)
555         w0(il,k) =max(w0(il,k),1.e-10)
556         sig(il,k)=amu/w0(il,k)
557         sig(il,k)=min(sig(il,k),1.)
[879]558cc         amu = 0.5*(SIG(il,k)+sigold(il,k))*W0(il,k)
559         M(il,k)=AMU*0.007*P(il,k)*(PH(il,k)-PH(il,k+1))/TV(il,k)
560         ENDIF
561        enddo
562      ENDDO
563cjyg2
564      DO il = 1,ncum
565       w0(il,icb(il))=0.5*w0(il,icb(il)+1)
566       m(il,icb(il))=0.5*m(il,icb(il)+1)
567     $       *(ph(il,icb(il))-ph(il,icb(il)+1))
568     $       /(ph(il,icb(il)+1)-ph(il,icb(il)+2))
569       sig(il,icb(il))=sig(il,icb(il)+1)
570       sig(il,icb(il)-1)=sig(il,icb(il))
571      ENDDO
[973]572      if(prt_level.GE.20) print*,'cv3p1_param apres w0_sig_M'
[879]573c
574cc 3. Compute final cloud base mass flux and set iflag to 3 if
575cc    cloud base mass flux is exceedingly small and is decreasing (i.e. if
576cc    the final mass flux (cbmflast) is greater than the target mass flux
577cc    (cbmf)).
578c
579      do il = 1,ncum
580       cbmflast(il) = 0.
581      enddo
582c
583      do k= 1,nl
584       do il = 1,ncum
585        IF (k .ge. icb(il) .and. k .le. inb(il)) THEN
[973]586 !IMpropo??      IF ((k.ge.(icb(il)+1)).and.(k.le.inb(il))) THEN
[879]587         cbmflast(il) = cbmflast(il)+M(il,k)
588        ENDIF
589       enddo
590      enddo
591c
592      do il = 1,ncum
593       IF (cbmflast(il) .lt. 1.e-6 .and.
594     $     cbmflast(il) .ge. cbmf(il)) THEN
595         iflag(il) = 3
596       ENDIF
597      enddo
598c
599      do k= 1,nl
600       do il = 1,ncum
601        IF (iflag(il) .ge. 3) THEN
602         M(il,k) = 0.
603         sig(il,k) = 0.
604         w0(il,k) = 0.
605        ENDIF
606       enddo
607      enddo
[973]608      if(prt_level.GE.20) print*,'cv3p1_param apres iflag'
[879]609c
610cc 4. Introduce a correcting factor for coef, in order to obtain an effective
611cc    sigdz larger in the present case (using cv3p1_closure) than in the old
612cc    closure (using cv3_closure).
[973]613      if (1.eq.0) then
614       do il = 1,ncum
615cc      coef(il) = 2.*coef(il)
616        coef(il) = 5.*coef(il)
617       enddo
618c version CVS du ..2008
619      else
620       if (iflag_cvl_sigd.eq.0) then
[879]621ctest pour verifier qu on fait la meme chose qu avant: sid constant
622        coef(1:ncum)=1.
[973]623       else
[879]624        coef(1:ncum) = min(2.*coef(1:ncum),5.)
625        coef(1:ncum) = max(2.*coef(1:ncum),0.2)
[973]626       endif
[879]627      endif
628c
[973]629      if(prt_level.GE.20) print*,'cv3p1_param FIN'
[879]630       return
631       end
632
633
Note: See TracBrowser for help on using the repository browser.