phylmd/CV30_routines/cv30_mixing.f90

module cv30_mixing_m

  implicit none

contains

  SUBROUTINE cv30_mixing(icb, inb, t, rr, rs, u, v, h, lv, hp, ep, clw, m, &
       sig, ment, qent, uent, vent, nent, sij, elij, ments, qents)

    ! MIXING

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

    use cv30_param_m, only: minorig, nl
    USE dimphy, ONLY: klev, klon
    use suphec_m, only: rcpd, rcpv, rv

    ! inputs:

    integer, intent(in):: icb(:) ! (ncum) {2 <= icb <= nl - 3}

    integer, intent(in):: inb(:) ! (ncum)
    ! first model level above the level of neutral buoyancy of the
    ! parcel (1 <= inb <= nl - 1)

    real, intent(in):: t(klon, klev), rr(klon, klev), rs(klon, klev)
    real u(klon, klev), v(klon, klev)
    real, intent(in):: h(klon, klev)
    real, intent(in):: lv(:, :) ! (klon, klev)
    real, intent(in):: hp(klon, klev)
    real ep(klon, klev), clw(klon, klev)
    real m(klon, klev) ! input of convect3
    real sig(klon, klev)

    ! outputs:
    real ment(klon, klev, klev), qent(klon, klev, klev)
    real uent(klon, klev, klev), vent(klon, klev, klev)
    integer, intent(out):: nent(:, 2:) ! (ncum, 2:nl - 1)
    real sij(klon, klev, klev), elij(klon, klev, klev)
    real ments(klon, klev, klev), qents(klon, klev, klev)

    ! Local:
    integer ncum, 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(klon), smax(klon), scrit(klon)
    real asum(klon, klev), bsum(klon, klev), csum(klon, klev)
    real wgh
    real zm(klon, klev)
    logical lwork(klon)

    !-------------------------------------------------------------------------

    ncum = size(icb)

    ! INITIALIZE VARIOUS ARRAYS USED IN THE COMPUTATIONS

    nent = 0

    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
          end do
       end do
    end do

    ment(1:ncum, 1:klev, 1:klev) = 0.0
    sij(1:ncum, 1:klev, 1:klev) = 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) / (rv &
                     * t(il, j) * t(il, j) * rcpd)
                anum = h(il, j) - hp(il, i) + (rcpv - rcpd) * t(il, j) * (rti &
                     - rr(il, j))
                denom = h(il, i) - hp(il, i) + (rcpd - rcpv) * (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, minorig)
                   vent(il, i, j) = sij(il, i, j) * v(il, i) + (1. &
                        - sij(il, i, j)) * v(il, minorig)
                   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
          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)) then
             if (nent(il, i) == 0) then
                ment(il, i, i) = m(il, i)
                qent(il, i, i) = rr(il, minorig) - ep(il, i) * clw(il, i)
                uent(il, i, i) = u(il, minorig)
                vent(il, i, i) = v(il, minorig)
                elij(il, i, i) = clw(il, i)
                sij(il, i, i) = 0.0
             end if
          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)) &
                  + (rcpv - rcpd) * t(il, i) * (qp - rr(il, i))
             denom = h(il, i) - hp(il, i) + lv(il, i) * (rr(il, i) - qp) &
                  + (rcpd - rcpv) * 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, minorig)
             vent(il, i, i) = v(il, minorig)
             elij(il, i, i) = clw(il, i)
             sij(il, i, i) = 0.0
          endif
       enddo ! il

    end DO

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