phylmd/CV30_routines/cv30_mixing.f

module cv30_mixing_m

  implicit none

contains

  SUBROUTINE cv30_mixing(nloc, ncum, nd, na, icb, nk, inb, t, rr, rs, u, v, h, &
       lv, hp, ep, clw, m, sig, ment, qent, uent, vent, nent, sij, elij, &
       ments, qents)
    use cv30_param_m
    use cvthermo

    !
    ! a faire:
    ! - changer rr(il, 1) -> qnk(il)
    ! - vectorisation de la partie normalisation des flux (do 789)
    !

    ! inputs:
    integer, intent(in):: ncum, nd, na, nloc
    integer, intent(in):: icb(nloc), inb(nloc), nk(nloc)
    real sig(nloc, nd)
    real t(nloc, nd), rr(nloc, nd), rs(nloc, nd)
    real u(nloc, nd), v(nloc, nd)
    real lv(nloc, na), h(nloc, na), hp(nloc, na)
    real ep(nloc, na), clw(nloc, na)
    real m(nloc, na) ! input of convect3

    ! outputs:
    real ment(nloc, na, na), qent(nloc, na, na)
    real uent(nloc, na, na), vent(nloc, na, na)
    real sij(nloc, na, na), elij(nloc, na, na)
    real ments(nloc, nd, nd), qents(nloc, nd, nd)
    integer nent(nloc, nd)

    ! local variables:
    integer i, j, k, il, im, jm
    integer num1, num2
    real rti, bf2, anum, denom, dei, altem, cwat, stemp, qp
    real alt, smid, sjmin, sjmax, delp, delm
    real asij(nloc), smax(nloc), scrit(nloc)
    real asum(nloc, nd), bsum(nloc, nd), csum(nloc, nd)
    real wgh
    real zm(nloc, na)
    logical lwork(nloc)

    ! INITIALIZE VARIOUS ARRAYS USED IN THE COMPUTATIONS

    do j=1, nl
       do i=1, ncum
          nent(i, j)=0
          ! in convect3, m is computed in cv30_closure
          ! ori m(i, 1)=0.0
       end do
    end do

    do j=1, nl
       do k=1, nl
          do i=1, ncum
             qent(i, k, j)=rr(i, j)
             uent(i, k, j)=u(i, j)
             vent(i, k, j)=v(i, j)
             elij(i, k, j)=0.0
             !ym ment(i, k, j)=0.0
             !ym sij(i, k, j)=0.0
          end do
       end do
    end do

    !ym
    ment(1:ncum, 1:nd, 1:nd)=0.0
    sij(1:ncum, 1:nd, 1:nd)=0.0

    zm(:, :)=0.

    ! CALCULATE ENTRAINED AIR MASS FLUX (ment), TOTAL WATER MIXING
    ! RATIO (QENT), TOTAL CONDENSED WATER (elij), AND MIXING
    ! FRACTION (sij)

    do i=minorig+1, nl

       do j=minorig, nl
          do il=1, ncum
             if((i >= icb(il)).and.(i <= inb(il)).and. &
                  (j >= (icb(il)-1)).and.(j <= inb(il)))then

                rti=rr(il, 1)-ep(il, i)*clw(il, i)
                bf2=1.+lv(il, j)*lv(il, j)*rs(il, j)/(rrv*t(il, j)*t(il, j)*cpd)
                anum=h(il, j)-hp(il, i)+(cpv-cpd)*t(il, j)*(rti-rr(il, j))
                denom=h(il, i)-hp(il, i)+(cpd-cpv)*(rr(il, i)-rti)*t(il, j)
                dei=denom
                if(abs(dei) < 0.01)dei=0.01
                sij(il, i, j)=anum/dei
                sij(il, i, i)=1.0
                altem=sij(il, i, j)*rr(il, i)+(1.-sij(il, i, j))*rti-rs(il, j)
                altem=altem/bf2
                cwat=clw(il, j)*(1.-ep(il, j))
                stemp=sij(il, i, j)
                if((stemp < 0.0.or.stemp > 1.0.or.altem > cwat) &
                     .and.j > i)then
                   anum=anum-lv(il, j)*(rti-rs(il, j)-cwat*bf2)
                   denom=denom+lv(il, j)*(rr(il, i)-rti)
                   if(abs(denom) < 0.01)denom=0.01
                   sij(il, i, j)=anum/denom
                   altem=sij(il, i, j)*rr(il, i)+(1.-sij(il, i, j))*rti-rs(il, j)
                   altem=altem-(bf2-1.)*cwat
                end if
                if(sij(il, i, j) > 0.0.and.sij(il, i, j) < 0.95)then
                   qent(il, i, j)=sij(il, i, j)*rr(il, i)+(1.-sij(il, i, j))*rti
                   uent(il, i, j)=sij(il, i, j)*u(il, i)+(1.-sij(il, i, j))*u(il, nk(il))
                   vent(il, i, j)=sij(il, i, j)*v(il, i)+(1.-sij(il, i, j))*v(il, nk(il))
                   elij(il, i, j)=altem
                   elij(il, i, j)=amax1(0.0, elij(il, i, j))
                   ment(il, i, j)=m(il, i)/(1.-sij(il, i, j))
                   nent(il, i)=nent(il, i)+1
                end if
                sij(il, i, j)=amax1(0.0, sij(il, i, j))
                sij(il, i, j)=amin1(1.0, sij(il, i, j))
             endif ! new
          end do
       end do

       !  if no air can entrain at level i assume that updraft detrains 
       !  at that level and calculate detrained air flux and properties 

       do il=1, ncum
          if ((i >= icb(il)).and.(i <= inb(il)).and.(nent(il, i) == 0)) then
             !@ if(nent(il, i) == 0)then
             ment(il, i, i)=m(il, i)
             qent(il, i, i)=rr(il, nk(il))-ep(il, i)*clw(il, i)
             uent(il, i, i)=u(il, nk(il))
             vent(il, i, i)=v(il, nk(il))
             elij(il, i, i)=clw(il, i)
             !MAF sij(il, i, i)=1.0
             sij(il, i, i)=0.0
          end if
       end do
    end do

    ! NORMALIZE ENTRAINED AIR MASS FLUXES
    ! TO REPRESENT EQUAL PROBABILITIES OF MIXING

    asum = 0.
    csum = 0.

    do il=1, ncum
       lwork(il) = .FALSE.
    enddo

    DO i=minorig+1, nl

       num1=0
       do il=1, ncum
          if (i >= icb(il) .and. i <= inb(il)) num1=num1+1
       enddo
       if (num1 <= 0) cycle

       do il=1, ncum
          if (i >= icb(il) .and. i <= inb(il)) then
             lwork(il)=(nent(il, i) /= 0)
             qp=rr(il, 1)-ep(il, i)*clw(il, i)
             anum=h(il, i)-hp(il, i)-lv(il, i)*(qp-rs(il, i)) &
                  +(cpv-cpd)*t(il, i)*(qp-rr(il, i))
             denom=h(il, i)-hp(il, i)+lv(il, i)*(rr(il, i)-qp) &
                  +(cpd-cpv)*t(il, i)*(rr(il, i)-qp)
             if(abs(denom) < 0.01)denom=0.01
             scrit(il)=anum/denom
             alt=qp-rs(il, i)+scrit(il)*(rr(il, i)-qp)
             if(scrit(il) <= 0.0.or.alt <= 0.0)scrit(il)=1.0
             smax(il)=0.0
             asij(il)=0.0
          endif
       end do

       do j=nl, minorig, -1

          num2=0
          do il=1, ncum
             if (i >= icb(il) .and. i <= inb(il) .and. &
                  j >= (icb(il)-1) .and. j <= inb(il) &
                  .and. lwork(il)) num2=num2+1
          enddo
          if (num2 <= 0) cycle

          do il=1, ncum
             if (i >= icb(il) .and. i <= inb(il) .and. &
                  j >= (icb(il)-1) .and. j <= inb(il) &
                  .and. lwork(il)) then

                if(sij(il, i, j) > 1.0e-16.and.sij(il, i, j) < 0.95)then
                   wgh=1.0
                   if(j > i)then
                      sjmax=amax1(sij(il, i, j+1), smax(il))
                      sjmax=amin1(sjmax, scrit(il))
                      smax(il)=amax1(sij(il, i, j), smax(il))
                      sjmin=amax1(sij(il, i, j-1), smax(il))
                      sjmin=amin1(sjmin, scrit(il))
                      if(sij(il, i, j) < (smax(il)-1.0e-16))wgh=0.0
                      smid=amin1(sij(il, i, j), scrit(il))
                   else
                      sjmax=amax1(sij(il, i, j+1), scrit(il))
                      smid=amax1(sij(il, i, j), scrit(il))
                      sjmin=0.0
                      if(j > 1)sjmin=sij(il, i, j-1)
                      sjmin=amax1(sjmin, scrit(il))
                   endif
                   delp=abs(sjmax-smid)
                   delm=abs(sjmin-smid)
                   asij(il)=asij(il)+wgh*(delp+delm)
                   ment(il, i, j)=ment(il, i, j)*(delp+delm)*wgh
                endif
             endif
          end do

       end do

       do il=1, ncum
          if (i >= icb(il).and.i <= inb(il).and.lwork(il)) then
             asij(il)=amax1(1.0e-16, asij(il))
             asij(il)=1.0/asij(il)
             asum(il, i)=0.0
             bsum(il, i)=0.0
             csum(il, i)=0.0
          endif
       enddo

       do j=minorig, nl
          do il=1, ncum
             if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
                  .and. j >= (icb(il)-1) .and. j <= inb(il)) then
                ment(il, i, j)=ment(il, i, j)*asij(il)
             endif
          enddo
       end do

       do j=minorig, nl
          do il=1, ncum
             if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
                  .and. j >= (icb(il)-1) .and. j <= inb(il)) then
                asum(il, i)=asum(il, i)+ment(il, i, j)
                ment(il, i, j)=ment(il, i, j)*sig(il, j)
                bsum(il, i)=bsum(il, i)+ment(il, i, j)
             endif
          enddo
       end do

       do il=1, ncum
          if (i >= icb(il).and.i <= inb(il).and.lwork(il)) then
             bsum(il, i)=amax1(bsum(il, i), 1.0e-16)
             bsum(il, i)=1.0/bsum(il, i)
          endif
       enddo

       do j=minorig, nl
          do il=1, ncum
             if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
                  .and. j >= (icb(il)-1) .and. j <= inb(il)) then
                ment(il, i, j)=ment(il, i, j)*asum(il, i)*bsum(il, i)
             endif
          enddo
       end do

       do j=minorig, nl
          do il=1, ncum
             if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
                  .and. j >= (icb(il)-1) .and. j <= inb(il)) then
                csum(il, i)=csum(il, i)+ment(il, i, j)
             endif
          enddo
       end do

       do il=1, ncum
          if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
               .and. csum(il, i) < m(il, i)) then
             nent(il, i)=0
             ment(il, i, i)=m(il, i)
             qent(il, i, i)=rr(il, 1)-ep(il, i)*clw(il, i)
             uent(il, i, i)=u(il, nk(il))
             vent(il, i, i)=v(il, nk(il))
             elij(il, i, i)=clw(il, i)
             !MAF sij(il, i, i)=1.0
             sij(il, i, i)=0.0
          endif
       enddo ! il

    end DO

    ! MAF: renormalisation de MENT
    do jm=1, nd
       do im=1, nd
          do il=1, ncum
             zm(il, im)=zm(il, im)+(1.-sij(il, im, jm))*ment(il, im, jm)
          end do
       end do
    end do

    do jm=1, nd
       do im=1, nd
          do il=1, ncum
             if(zm(il, im) /= 0.) then
                ment(il, im, jm)=ment(il, im, jm)*m(il, im)/zm(il, im)
             endif
          end do
       end do
    end do

    do jm=1, nd
       do im=1, nd
          do il=1, ncum
             qents(il, im, jm)=qent(il, im, jm)
             ments(il, im, jm)=ment(il, im, jm)
          end do
       enddo
    enddo

  end SUBROUTINE cv30_mixing

end module cv30_mixing_m
1	guez	185	module cv30_mixing_m
2	guez	47
3	guez	97	implicit none
4	guez	47
5	guez	97	contains
6	guez	47
7	guez	185	SUBROUTINE cv30_mixing(nloc, ncum, nd, na, icb, nk, inb, t, rr, rs, u, v, h, &
8	guez	145	lv, hp, ep, clw, m, sig, ment, qent, uent, vent, nent, sij, elij, &
9			ments, qents)
10	guez	185	use cv30_param_m
11	guez	97	use cvthermo
12	guez	47
13	guez	189	!
14	guez	97	! a faire:
15	guez	145	! - changer rr(il, 1) -> qnk(il)
16			! - vectorisation de la partie normalisation des flux (do 789)
17	guez	189	!
18	guez	47
19	guez	97	! inputs:
20			integer, intent(in):: ncum, nd, na, nloc
21	guez	187	integer, intent(in):: icb(nloc), inb(nloc), nk(nloc)
22	guez	145	real sig(nloc, nd)
23			real t(nloc, nd), rr(nloc, nd), rs(nloc, nd)
24			real u(nloc, nd), v(nloc, nd)
25			real lv(nloc, na), h(nloc, na), hp(nloc, na)
26			real ep(nloc, na), clw(nloc, na)
27			real m(nloc, na) ! input of convect3
28	guez	47
29	guez	97	! outputs:
30	guez	145	real ment(nloc, na, na), qent(nloc, na, na)
31			real uent(nloc, na, na), vent(nloc, na, na)
32			real sij(nloc, na, na), elij(nloc, na, na)
33			real ments(nloc, nd, nd), qents(nloc, nd, nd)
34			integer nent(nloc, nd)
35	guez	47
36	guez	97	! local variables:
37			integer i, j, k, il, im, jm
38			integer num1, num2
39			real rti, bf2, anum, denom, dei, altem, cwat, stemp, qp
40			real alt, smid, sjmin, sjmax, delp, delm
41			real asij(nloc), smax(nloc), scrit(nloc)
42	guez	145	real asum(nloc, nd), bsum(nloc, nd), csum(nloc, nd)
43	guez	97	real wgh
44	guez	145	real zm(nloc, na)
45	guez	97	logical lwork(nloc)
46
47	guez	189	! INITIALIZE VARIOUS ARRAYS USED IN THE COMPUTATIONS
48	guez	97
49	guez	145	do j=1, nl
50			do i=1, ncum
51			nent(i, j)=0
52	guez	185	! in convect3, m is computed in cv30_closure
53	guez	145	! ori m(i, 1)=0.0
54	guez	97	end do
55			end do
56	guez	47
57	guez	145	do j=1, nl
58			do k=1, nl
59			do i=1, ncum
60			qent(i, k, j)=rr(i, j)
61			uent(i, k, j)=u(i, j)
62			vent(i, k, j)=v(i, j)
63			elij(i, k, j)=0.0
64			!ym ment(i, k, j)=0.0
65			!ym sij(i, k, j)=0.0
66	guez	97	end do
67			end do
68			end do
69	guez	47
70	guez	97	!ym
71	guez	145	ment(1:ncum, 1:nd, 1:nd)=0.0
72			sij(1:ncum, 1:nd, 1:nd)=0.0
73	guez	47
74	guez	145	zm(:, :)=0.
75	guez	47
76	guez	189	! CALCULATE ENTRAINED AIR MASS FLUX (ment), TOTAL WATER MIXING
77			! RATIO (QENT), TOTAL CONDENSED WATER (elij), AND MIXING
78			! FRACTION (sij)
79	guez	47
80	guez	145	do i=minorig+1, nl
81	guez	47
82	guez	145	do j=minorig, nl
83			do il=1, ncum
84			if((i >= icb(il)).and.(i <= inb(il)).and. &
85			(j >= (icb(il)-1)).and.(j <= inb(il)))then
86	guez	47
87	guez	145	rti=rr(il, 1)-ep(il, i)*clw(il, i)
88			bf2=1.+lv(il, j)lv(il, j)rs(il, j)/(rrvt(il, j)t(il, j)*cpd)
89			anum=h(il, j)-hp(il, i)+(cpv-cpd)t(il, j)(rti-rr(il, j))
90			denom=h(il, i)-hp(il, i)+(cpd-cpv)(rr(il, i)-rti)t(il, j)
91	guez	97	dei=denom
92	guez	145	if(abs(dei) < 0.01)dei=0.01
93			sij(il, i, j)=anum/dei
94			sij(il, i, i)=1.0
95			altem=sij(il, i, j)rr(il, i)+(1.-sij(il, i, j))rti-rs(il, j)
96	guez	97	altem=altem/bf2
97	guez	145	cwat=clw(il, j)*(1.-ep(il, j))
98			stemp=sij(il, i, j)
99			if((stemp < 0.0.or.stemp > 1.0.or.altem > cwat) &
100			.and.j > i)then
101			anum=anum-lv(il, j)(rti-rs(il, j)-cwatbf2)
102			denom=denom+lv(il, j)*(rr(il, i)-rti)
103			if(abs(denom) < 0.01)denom=0.01
104			sij(il, i, j)=anum/denom
105			altem=sij(il, i, j)rr(il, i)+(1.-sij(il, i, j))rti-rs(il, j)
106	guez	97	altem=altem-(bf2-1.)*cwat
107			end if
108	guez	145	if(sij(il, i, j) > 0.0.and.sij(il, i, j) < 0.95)then
109			qent(il, i, j)=sij(il, i, j)rr(il, i)+(1.-sij(il, i, j))rti
110			uent(il, i, j)=sij(il, i, j)u(il, i)+(1.-sij(il, i, j))u(il, nk(il))
111			vent(il, i, j)=sij(il, i, j)v(il, i)+(1.-sij(il, i, j))v(il, nk(il))
112			elij(il, i, j)=altem
113			elij(il, i, j)=amax1(0.0, elij(il, i, j))
114			ment(il, i, j)=m(il, i)/(1.-sij(il, i, j))
115			nent(il, i)=nent(il, i)+1
116	guez	97	end if
117	guez	145	sij(il, i, j)=amax1(0.0, sij(il, i, j))
118			sij(il, i, j)=amin1(1.0, sij(il, i, j))
119	guez	97	endif ! new
120			end do
121			end do
122	guez	47
123	guez	189	! if no air can entrain at level i assume that updraft detrains
124			! at that level and calculate detrained air flux and properties
125	guez	47
126	guez	145	do il=1, ncum
127			if ((i >= icb(il)).and.(i <= inb(il)).and.(nent(il, i) == 0)) then
128			!@ if(nent(il, i) == 0)then
129			ment(il, i, i)=m(il, i)
130			qent(il, i, i)=rr(il, nk(il))-ep(il, i)*clw(il, i)
131			uent(il, i, i)=u(il, nk(il))
132			vent(il, i, i)=v(il, nk(il))
133			elij(il, i, i)=clw(il, i)
134			!MAF sij(il, i, i)=1.0
135			sij(il, i, i)=0.0
136	guez	97	end if
137			end do
138			end do
139	guez	47
140	guez	189	! NORMALIZE ENTRAINED AIR MASS FLUXES
141			! TO REPRESENT EQUAL PROBABILITIES OF MIXING
142	guez	47
143	guez	186	asum = 0.
144			csum = 0.
145	guez	47
146	guez	145	do il=1, ncum
147	guez	47	lwork(il) = .FALSE.
148	guez	97	enddo
149	guez	47
150	guez	145	DO i=minorig+1, nl
151	guez	47
152	guez	97	num1=0
153	guez	145	do il=1, ncum
154			if (i >= icb(il) .and. i <= inb(il)) num1=num1+1
155	guez	97	enddo
156	guez	145	if (num1 <= 0) cycle
157	guez	47
158	guez	145	do il=1, ncum
159			if (i >= icb(il) .and. i <= inb(il)) then
160			lwork(il)=(nent(il, i) /= 0)
161			qp=rr(il, 1)-ep(il, i)*clw(il, i)
162			anum=h(il, i)-hp(il, i)-lv(il, i)*(qp-rs(il, i)) &
163			+(cpv-cpd)t(il, i)(qp-rr(il, i))
164			denom=h(il, i)-hp(il, i)+lv(il, i)*(rr(il, i)-qp) &
165			+(cpd-cpv)t(il, i)(rr(il, i)-qp)
166			if(abs(denom) < 0.01)denom=0.01
167	guez	97	scrit(il)=anum/denom
168	guez	145	alt=qp-rs(il, i)+scrit(il)*(rr(il, i)-qp)
169			if(scrit(il) <= 0.0.or.alt <= 0.0)scrit(il)=1.0
170	guez	97	smax(il)=0.0
171			asij(il)=0.0
172			endif
173			end do
174	guez	47
175	guez	145	do j=nl, minorig, -1
176	guez	47
177	guez	97	num2=0
178	guez	145	do il=1, ncum
179			if (i >= icb(il) .and. i <= inb(il) .and. &
180			j >= (icb(il)-1) .and. j <= inb(il) &
181			.and. lwork(il)) num2=num2+1
182	guez	97	enddo
183	guez	145	if (num2 <= 0) cycle
184	guez	47
185	guez	145	do il=1, ncum
186			if (i >= icb(il) .and. i <= inb(il) .and. &
187			j >= (icb(il)-1) .and. j <= inb(il) &
188			.and. lwork(il)) then
189	guez	47
190	guez	145	if(sij(il, i, j) > 1.0e-16.and.sij(il, i, j) < 0.95)then
191	guez	97	wgh=1.0
192	guez	145	if(j > i)then
193			sjmax=amax1(sij(il, i, j+1), smax(il))
194			sjmax=amin1(sjmax, scrit(il))
195			smax(il)=amax1(sij(il, i, j), smax(il))
196			sjmin=amax1(sij(il, i, j-1), smax(il))
197			sjmin=amin1(sjmin, scrit(il))
198			if(sij(il, i, j) < (smax(il)-1.0e-16))wgh=0.0
199			smid=amin1(sij(il, i, j), scrit(il))
200	guez	97	else
201	guez	145	sjmax=amax1(sij(il, i, j+1), scrit(il))
202			smid=amax1(sij(il, i, j), scrit(il))
203	guez	97	sjmin=0.0
204	guez	145	if(j > 1)sjmin=sij(il, i, j-1)
205			sjmin=amax1(sjmin, scrit(il))
206	guez	97	endif
207			delp=abs(sjmax-smid)
208			delm=abs(sjmin-smid)
209			asij(il)=asij(il)+wgh*(delp+delm)
210	guez	145	ment(il, i, j)=ment(il, i, j)(delp+delm)wgh
211	guez	97	endif
212			endif
213			end do
214	guez	47
215	guez	97	end do
216	guez	47
217	guez	145	do il=1, ncum
218			if (i >= icb(il).and.i <= inb(il).and.lwork(il)) then
219			asij(il)=amax1(1.0e-16, asij(il))
220	guez	97	asij(il)=1.0/asij(il)
221	guez	145	asum(il, i)=0.0
222			bsum(il, i)=0.0
223			csum(il, i)=0.0
224	guez	97	endif
225	guez	47	enddo
226
227	guez	145	do j=minorig, nl
228			do il=1, ncum
229			if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
230			.and. j >= (icb(il)-1) .and. j <= inb(il)) then
231			ment(il, i, j)=ment(il, i, j)*asij(il)
232	guez	97	endif
233			enddo
234			end do
235	guez	47
236	guez	145	do j=minorig, nl
237			do il=1, ncum
238			if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
239			.and. j >= (icb(il)-1) .and. j <= inb(il)) then
240			asum(il, i)=asum(il, i)+ment(il, i, j)
241			ment(il, i, j)=ment(il, i, j)*sig(il, j)
242			bsum(il, i)=bsum(il, i)+ment(il, i, j)
243	guez	97	endif
244			enddo
245			end do
246	guez	47
247	guez	145	do il=1, ncum
248			if (i >= icb(il).and.i <= inb(il).and.lwork(il)) then
249			bsum(il, i)=amax1(bsum(il, i), 1.0e-16)
250			bsum(il, i)=1.0/bsum(il, i)
251	guez	97	endif
252	guez	47	enddo
253
254	guez	145	do j=minorig, nl
255			do il=1, ncum
256			if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
257			.and. j >= (icb(il)-1) .and. j <= inb(il)) then
258			ment(il, i, j)=ment(il, i, j)asum(il, i)bsum(il, i)
259	guez	97	endif
260			enddo
261			end do
262
263	guez	145	do j=minorig, nl
264			do il=1, ncum
265			if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
266			.and. j >= (icb(il)-1) .and. j <= inb(il)) then
267			csum(il, i)=csum(il, i)+ment(il, i, j)
268	guez	97	endif
269			enddo
270			end do
271
272	guez	145	do il=1, ncum
273			if (i >= icb(il) .and. i <= inb(il) .and. lwork(il) &
274			.and. csum(il, i) < m(il, i)) then
275			nent(il, i)=0
276			ment(il, i, i)=m(il, i)
277			qent(il, i, i)=rr(il, 1)-ep(il, i)*clw(il, i)
278			uent(il, i, i)=u(il, nk(il))
279			vent(il, i, i)=v(il, nk(il))
280			elij(il, i, i)=clw(il, i)
281			!MAF sij(il, i, i)=1.0
282			sij(il, i, i)=0.0
283	guez	97	endif
284			enddo ! il
285	guez	47
286	guez	97	end DO
287	guez	145
288	guez	97	! MAF: renormalisation de MENT
289	guez	145	do jm=1, nd
290			do im=1, nd
291			do il=1, ncum
292			zm(il, im)=zm(il, im)+(1.-sij(il, im, jm))*ment(il, im, jm)
293	guez	97	end do
294			end do
295			end do
296	guez	145
297			do jm=1, nd
298			do im=1, nd
299			do il=1, ncum
300			if(zm(il, im) /= 0.) then
301			ment(il, im, jm)=ment(il, im, jm)*m(il, im)/zm(il, im)
302	guez	97	endif
303			end do
304	guez	47	end do
305	guez	97	end do
306	guez	145
307			do jm=1, nd
308			do im=1, nd
309			do il=1, ncum
310			qents(il, im, jm)=qent(il, im, jm)
311			ments(il, im, jm)=ment(il, im, jm)
312	guez	97	end do
313	guez	47	enddo
314	guez	97	enddo
315	guez	47
316	guez	185	end SUBROUTINE cv30_mixing
317	guez	97
318	guez	185	end module cv30_mixing_m