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

Annotation of /trunk/Sources/phylmd/CV30_routines/cv30_closure.f

Parent Directory Parent Directory | Revision Log Revision Log

Revision 201 - (hide annotations)
Mon Jun 6 17:42:15 2016 UTC (7 years, 11 months ago) by guez
File size: 3955 byte(s) 
Removed intermediary objects of cv_thermo_m, access suphec_m
directly. Procedure cv_thermo disappeared, all objects are named
constants.

In cv_driver and below, limited extents of arrays to what is needed.

lv, cpn and th in cv30_compress were set at level nl + 1 but lv1, cpn1
and th1 are not defined at this level. This did not lead to an error
because values at nl + 1 were not used.

Removed test on ok_sync in phystokenc because it is not read at run
time. Printing min and max of output NetCDF variables is heavy and
archaic.

Used histwrite_phy in phytrac.
1 guez 187 module cv30_closure_m
2 guez 47
3 guez 187 implicit none
4 guez 47
5 guez 187 contains
6 guez 47
7 guez 195 SUBROUTINE cv30_closure(icb, inb, pbase, p, ph, tv, buoy, sig, w0, cape, m)
8 guez 47
9 guez 195 ! CLOSURE
10 guez 187 ! Vectorization: S. Bony
11 guez 47
12 guez 187 use cv30_param_m, only: alpha, beta, dtcrit, minorig, nl
13 guez 195 USE dimphy, ONLY: klev, klon
14 guez 201 use suphec_m, only: rd
15 guez 47
16 guez 187 ! input:
17 guez 195 integer, intent(in):: icb(:), inb(:) ! (ncum)
18     real pbase(klon)
19 guez 196 real p(:, :) ! (klon, klev)
20     real, intent(in):: ph(:, :) ! (ncum, klev + 1)
21 guez 195 real tv(klon, klev), buoy(klon, klev)
22 guez 47
23 guez 187 ! input/output:
24 guez 195 real sig(klon, klev), w0(klon, klev)
25 guez 47
26 guez 187 ! output:
27 guez 195 real cape(klon)
28     real m(klon, klev)
29 guez 47
30 guez 187 ! Local:
31 guez 195 integer ncum
32 guez 187 integer i, j, k, icbmax
33     real deltap, fac, w, amu
34 guez 195 real dtmin(klon, klev), sigold(klon, klev)
35 guez 47
36 guez 187 !-------------------------------------------------------
37 guez 47
38 guez 195 ncum = size(icb)
39    
40 guez 187 ! Initialization
41 guez 47
42 guez 187 do k=1, nl
43     do i=1, ncum
44     m(i, k)=0.0
45     enddo
46     enddo
47 guez 47
48 guez 187 ! Reset sig(i) and w0(i) for i>inb and i<icb
49 guez 47
50 guez 187 ! update sig and w0 above LNB:
51 guez 47
52 guez 187 do k=1, nl-1
53     do i=1, ncum
54 guez 196 if ((inb(i) < (nl-1)).and.(k >= (inb(i) + 1)))then
55 guez 187 sig(i, k)=beta*sig(i, k) &
56 guez 196 + 2.*alpha*buoy(i, inb(i))*ABS(buoy(i, inb(i)))
57 guez 187 sig(i, k)=AMAX1(sig(i, k), 0.0)
58     w0(i, k)=beta*w0(i, k)
59     endif
60     end do
61     end do
62 guez 47
63 guez 187 ! compute icbmax:
64 guez 47
65 guez 187 icbmax=2
66     do i=1, ncum
67     icbmax=MAX(icbmax, icb(i))
68     end do
69 guez 47
70 guez 187 ! update sig and w0 below cloud base:
71 guez 47
72 guez 187 do k=1, icbmax
73     do i=1, ncum
74     if (k <= icb(i))then
75     sig(i, k)=beta*sig(i, k)-2.*alpha*buoy(i, icb(i))*buoy(i, icb(i))
76     sig(i, k)=amax1(sig(i, k), 0.0)
77     w0(i, k)=beta*w0(i, k)
78     endif
79     end do
80     end do
81 guez 47
82 guez 187 ! Reset fractional areas of updrafts and w0 at initial time
83     ! and after 10 time steps of no convection
84 guez 47
85 guez 187 do k=1, nl-1
86     do i=1, ncum
87 guez 195 if (sig(i, klev) < 1.5.or.sig(i, klev) > 12.0)then
88 guez 187 sig(i, k)=0.0
89     w0(i, k)=0.0
90     endif
91     end do
92     end do
93 guez 47
94 guez 187 ! Calculate convective available potential energy (cape),
95     ! vertical velocity (w), fractional area covered by
96     ! undilute updraft (sig), and updraft mass flux (m)
97 guez 47
98 guez 187 do i=1, ncum
99 guez 47 cape(i)=0.0
100 guez 187 end do
101 guez 47
102 guez 187 ! compute dtmin (minimum buoyancy between ICB and given level k):
103 guez 47
104 guez 187 do i=1, ncum
105     do k=1, nl
106     dtmin(i, k)=100.0
107 guez 47 enddo
108 guez 187 enddo
109 guez 47
110 guez 187 do i=1, ncum
111     do k=1, nl
112     do j=minorig, nl
113 guez 196 if ((k >= (icb(i) + 1)).and.(k <= inb(i)).and. &
114 guez 187 (j >= icb(i)).and.(j <= (k-1)))then
115     dtmin(i, k)=AMIN1(dtmin(i, k), buoy(i, j))
116     endif
117     end do
118     end do
119     end do
120 guez 47
121 guez 196 ! The interval on which cape is computed starts at pbase:
122 guez 47
123 guez 187 do k=1, nl
124     do i=1, ncum
125 guez 196 if ((k >= (icb(i) + 1)).and.(k <= inb(i))) then
126 guez 187 deltap = MIN(pbase(i), ph(i, k-1))-MIN(pbase(i), ph(i, k))
127 guez 201 cape(i)=cape(i) + rd*buoy(i, k-1)*deltap/p(i, k-1)
128 guez 187 cape(i)=AMAX1(0.0, cape(i))
129     sigold(i, k)=sig(i, k)
130 guez 47
131 guez 196 sig(i, k)=beta*sig(i, k) + alpha*dtmin(i, k)*ABS(dtmin(i, k))
132 guez 187 sig(i, k)=amax1(sig(i, k), 0.0)
133     sig(i, k)=amin1(sig(i, k), 0.01)
134     fac=AMIN1(((dtcrit-dtmin(i, k))/dtcrit), 1.0)
135 guez 196 w=(1.-beta)*fac*SQRT(cape(i)) + beta*w0(i, k)
136     amu=0.5*(sig(i, k) + sigold(i, k))*w
137     m(i, k)=amu*0.007*p(i, k)*(ph(i, k)-ph(i, k + 1))/tv(i, k)
138 guez 187 w0(i, k)=w
139     endif
140 guez 47
141 guez 187 end do
142     end do
143 guez 47
144 guez 187 do i=1, ncum
145 guez 196 w0(i, icb(i))=0.5*w0(i, icb(i) + 1)
146     m(i, icb(i))=0.5*m(i, icb(i) + 1) &
147     *(ph(i, icb(i))-ph(i, icb(i) + 1)) &
148     /(ph(i, icb(i) + 1)-ph(i, icb(i) + 2))
149     sig(i, icb(i))=sig(i, icb(i) + 1)
150 guez 187 sig(i, icb(i)-1)=sig(i, icb(i))
151     end do
152 guez 47
153 guez 187 end SUBROUTINE cv30_closure
154 guez 47
155 guez 187 end module cv30_closure_m
  ViewVC Help
Powered by ViewVC 1.1.21