/[lmdze]/trunk/Sources/phylmd/CV30_routines/cv30_uncompress.f
ViewVC logotype

Contents of /trunk/Sources/phylmd/CV30_routines/cv30_uncompress.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 196 - (show annotations)
Mon May 23 13:50:39 2016 UTC (7 years, 11 months ago) by guez
File size: 2778 byte(s)
Removed argument icbmax of cv30_feed, not used in cv_driver (not used
in LMDZ either).

Clearer to have iflag1 = 0 in cv30_feed than in cv_driver. Clearer to
have iflag1 = 42 in cv30_uncompress than in cv_driver.

Removed argument iflag of cv30_compress. Since there is iflag1 = 42 in
cv30_uncompress, iflag1 that comes out of cv_driver is disconnected
from iflag1 computed by cv30_feed and cv30_trigger.

Removed some references to convect3 and convect4 in comments. This
program derives from the convect3 version, we do not need to know
about other versions.

Bug fix in cv30_undilute1: icbs1 was not made >= 2.

1 module cv30_uncompress_m
2
3 implicit none
4
5 contains
6
7 SUBROUTINE cv30_uncompress(idcum, iflag, precip, VPrecip, sig, w0, ft, fq, &
8 fu, fv, inb, Ma, upwd, dnwd, dnwd0, qcondc, cape, da, phi, mp, iflag1, &
9 precip1, VPrecip1, sig1, w01, ft1, fq1, fu1, fv1, inb1, Ma1, upwd1, &
10 dnwd1, dnwd01, qcondc1, cape1, da1, phi1, mp1)
11
12 ! UNCOMPRESS THE FIELDS
13
14 USE cv30_param_m, ONLY: nl
15 use dimphy, only: klon, klev
16
17 integer, intent(in):: idcum(:) ! (ncum)
18 integer, intent(in):: iflag(:) ! (ncum)
19 real, intent(in):: precip(klon)
20 real, intent(in):: VPrecip(klon, klev+1)
21 real, intent(in):: sig(klon, klev), w0(klon, klev)
22 real, intent(in), dimension(klon, klev):: ft, fq, fu, fv
23 integer, intent(in):: inb(klon)
24 real, intent(in):: Ma(klon, klev)
25 real, intent(in):: upwd(klon, klev), dnwd(klon, klev), dnwd0(klon, klev)
26 real, intent(in):: qcondc(klon, klev)
27 real, intent(in):: cape(klon)
28 real, intent(in):: da(klon, klev), phi(klon, klev, klev), mp(klon, klev)
29
30 ! outputs:
31 integer, intent(out):: iflag1(:) ! (klon)
32 real precip1(klon)
33 real VPrecip1(klon, klev+1)
34 real sig1(klon, klev), w01(klon, klev)
35 real ft1(klon, klev), fq1(klon, klev), fu1(klon, klev), fv1(klon, klev)
36 integer, intent(inout):: inb1(klon)
37 real Ma1(klon, klev)
38 real upwd1(klon, klev), dnwd1(klon, klev), dnwd01(klon, klev)
39 real qcondc1(klon, klev)
40 real cape1(klon)
41 real, intent(inout):: da1(klon, klev), phi1(klon, klev, klev)
42 real, intent(inout):: mp1(klon, klev)
43
44 ! Local:
45 integer ncum, i, k, j
46
47 !-------------------------------------------------------------------
48
49 ncum = size(idcum)
50 iflag1 = 42 ! for non convective points
51
52 do i=1, ncum
53 precip1(idcum(i))=precip(i)
54 iflag1(idcum(i))=iflag(i)
55 inb1(idcum(i))=inb(i)
56 cape1(idcum(i))=cape(i)
57 end do
58
59 do k=1, nl
60 do i=1, ncum
61 VPrecip1(idcum(i), k)=VPrecip(i, k)
62 sig1(idcum(i), k)=sig(i, k)
63 w01(idcum(i), k)=w0(i, k)
64 ft1(idcum(i), k)=ft(i, k)
65 fq1(idcum(i), k)=fq(i, k)
66 fu1(idcum(i), k)=fu(i, k)
67 fv1(idcum(i), k)=fv(i, k)
68 Ma1(idcum(i), k)=Ma(i, k)
69 upwd1(idcum(i), k)=upwd(i, k)
70 dnwd1(idcum(i), k)=dnwd(i, k)
71 dnwd01(idcum(i), k)=dnwd0(i, k)
72 qcondc1(idcum(i), k)=qcondc(i, k)
73 da1(idcum(i), k)=da(i, k)
74 mp1(idcum(i), k)=mp(i, k)
75 end do
76 end do
77
78 do i=1, ncum
79 sig1(idcum(i), klev)=sig(i, klev)
80 end do
81
82 do j=1, klev
83 do k=1, klev
84 do i=1, ncum
85 phi1(idcum(i), k, j)=phi(i, k, j)
86 end do
87 end do
88 end do
89
90 end SUBROUTINE cv30_uncompress
91
92 end module cv30_uncompress_m

  ViewVC Help
Powered by ViewVC 1.1.21