1 |
|
module cv30_closure_m |
2 |
|
|
3 |
SUBROUTINE cv3_closure(nloc,ncum,nd,icb,inb & |
implicit none |
4 |
,pbase,p,ph,tv,buoy & |
|
5 |
,sig,w0,cape,m) |
contains |
6 |
use cv3_param_m |
|
7 |
use cvthermo |
SUBROUTINE cv30_closure(icb, inb, pbase, p, ph, tv, buoy, sig, w0, cape, m) |
8 |
implicit none |
|
9 |
|
! CLOSURE |
10 |
!=================================================================== |
! Vectorization: S. Bony |
11 |
! --- CLOSURE OF CONVECT3 |
|
12 |
! |
use cv30_param_m, only: alpha, beta, dtcrit, minorig, nl |
13 |
! vectorization: S. Bony |
use cv_thermo_m, only: rrd |
14 |
!=================================================================== |
USE dimphy, ONLY: klev, klon |
15 |
|
|
16 |
|
! input: |
17 |
! input: |
integer, intent(in):: icb(:), inb(:) ! (ncum) |
18 |
integer, intent(in):: ncum, nd, nloc |
real pbase(klon) |
19 |
integer icb(nloc), inb(nloc) |
real p(:, :) ! (klon, klev) |
20 |
real pbase(nloc) |
real, intent(in):: ph(:, :) ! (ncum, klev + 1) |
21 |
real p(nloc,nd), ph(nloc,nd+1) |
real tv(klon, klev), buoy(klon, klev) |
22 |
real tv(nloc,nd), buoy(nloc,nd) |
|
23 |
|
! input/output: |
24 |
! input/output: |
real sig(klon, klev), w0(klon, klev) |
25 |
real sig(nloc,nd), w0(nloc,nd) |
|
26 |
|
! output: |
27 |
! output: |
real cape(klon) |
28 |
real cape(nloc) |
real m(klon, klev) |
29 |
real m(nloc,nd) |
|
30 |
|
! Local: |
31 |
! local variables: |
integer ncum |
32 |
integer i, j, k, icbmax |
integer i, j, k, icbmax |
33 |
real deltap, fac, w, amu |
real deltap, fac, w, amu |
34 |
real dtmin(nloc,nd), sigold(nloc,nd) |
real dtmin(klon, klev), sigold(klon, klev) |
35 |
|
|
36 |
|
!------------------------------------------------------- |
37 |
! ------------------------------------------------------- |
|
38 |
! -- Initialization |
ncum = size(icb) |
39 |
! ------------------------------------------------------- |
|
40 |
|
! Initialization |
41 |
do k=1,nl |
|
42 |
do i=1,ncum |
do k=1, nl |
43 |
m(i,k)=0.0 |
do i=1, ncum |
44 |
|
m(i, k)=0.0 |
45 |
enddo |
enddo |
46 |
enddo |
enddo |
47 |
|
|
48 |
|
! Reset sig(i) and w0(i) for i>inb and i<icb |
49 |
|
|
50 |
! ------------------------------------------------------- |
! update sig and w0 above LNB: |
|
! -- Reset sig(i) and w0(i) for i>inb and i<icb |
|
|
! ------------------------------------------------------- |
|
|
|
|
|
! update sig and w0 above LNB: |
|
|
|
|
|
do 100 k=1,nl-1 |
|
|
do 110 i=1,ncum |
|
|
if ((inb(i).lt.(nl-1)).and.(k.ge.(inb(i)+1)))then |
|
|
sig(i,k)=beta*sig(i,k) & |
|
|
+2.*alpha*buoy(i,inb(i))*ABS(buoy(i,inb(i))) |
|
|
sig(i,k)=AMAX1(sig(i,k),0.0) |
|
|
w0(i,k)=beta*w0(i,k) |
|
|
endif |
|
|
110 continue |
|
|
100 continue |
|
|
|
|
|
! compute icbmax: |
|
|
|
|
|
icbmax=2 |
|
|
do 200 i=1,ncum |
|
|
icbmax=MAX(icbmax,icb(i)) |
|
|
200 continue |
|
|
|
|
|
! update sig and w0 below cloud base: |
|
|
|
|
|
do 300 k=1,icbmax |
|
|
do 310 i=1,ncum |
|
|
if (k.le.icb(i))then |
|
|
sig(i,k)=beta*sig(i,k)-2.*alpha*buoy(i,icb(i))*buoy(i,icb(i)) |
|
|
sig(i,k)=amax1(sig(i,k),0.0) |
|
|
w0(i,k)=beta*w0(i,k) |
|
|
endif |
|
|
310 continue |
|
|
300 continue |
|
|
|
|
|
!! if(inb.lt.(nl-1))then |
|
|
!! do 85 i=inb+1,nl-1 |
|
|
!! sig(i)=beta*sig(i)+2.*alpha*buoy(inb)* |
|
|
!! 1 abs(buoy(inb)) |
|
|
!! sig(i)=amax1(sig(i),0.0) |
|
|
!! w0(i)=beta*w0(i) |
|
|
!! 85 continue |
|
|
!! end if |
|
|
|
|
|
!! do 87 i=1,icb |
|
|
!! sig(i)=beta*sig(i)-2.*alpha*buoy(icb)*buoy(icb) |
|
|
!! sig(i)=amax1(sig(i),0.0) |
|
|
!! w0(i)=beta*w0(i) |
|
|
!! 87 continue |
|
|
|
|
|
! ------------------------------------------------------------- |
|
|
! -- Reset fractional areas of updrafts and w0 at initial time |
|
|
! -- and after 10 time steps of no convection |
|
|
! ------------------------------------------------------------- |
|
|
|
|
|
do 400 k=1,nl-1 |
|
|
do 410 i=1,ncum |
|
|
if (sig(i,nd).lt.1.5.or.sig(i,nd).gt.12.0)then |
|
|
sig(i,k)=0.0 |
|
|
w0(i,k)=0.0 |
|
|
endif |
|
|
410 continue |
|
|
400 continue |
|
|
|
|
|
! ------------------------------------------------------------- |
|
|
! -- Calculate convective available potential energy (cape), |
|
|
! -- vertical velocity (w), fractional area covered by |
|
|
! -- undilute updraft (sig), and updraft mass flux (m) |
|
|
! ------------------------------------------------------------- |
|
51 |
|
|
52 |
do 500 i=1,ncum |
do k=1, nl-1 |
53 |
|
do i=1, ncum |
54 |
|
if ((inb(i) < (nl-1)).and.(k >= (inb(i) + 1)))then |
55 |
|
sig(i, k)=beta*sig(i, k) & |
56 |
|
+ 2.*alpha*buoy(i, inb(i))*ABS(buoy(i, inb(i))) |
57 |
|
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 |
|
|
63 |
|
! compute icbmax: |
64 |
|
|
65 |
|
icbmax=2 |
66 |
|
do i=1, ncum |
67 |
|
icbmax=MAX(icbmax, icb(i)) |
68 |
|
end do |
69 |
|
|
70 |
|
! update sig and w0 below cloud base: |
71 |
|
|
72 |
|
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 |
|
|
82 |
|
! Reset fractional areas of updrafts and w0 at initial time |
83 |
|
! and after 10 time steps of no convection |
84 |
|
|
85 |
|
do k=1, nl-1 |
86 |
|
do i=1, ncum |
87 |
|
if (sig(i, klev) < 1.5.or.sig(i, klev) > 12.0)then |
88 |
|
sig(i, k)=0.0 |
89 |
|
w0(i, k)=0.0 |
90 |
|
endif |
91 |
|
end do |
92 |
|
end do |
93 |
|
|
94 |
|
! Calculate convective available potential energy (cape), |
95 |
|
! vertical velocity (w), fractional area covered by |
96 |
|
! undilute updraft (sig), and updraft mass flux (m) |
97 |
|
|
98 |
|
do i=1, ncum |
99 |
cape(i)=0.0 |
cape(i)=0.0 |
100 |
500 continue |
end do |
101 |
|
|
102 |
! compute dtmin (minimum buoyancy between ICB and given level k): |
! compute dtmin (minimum buoyancy between ICB and given level k): |
103 |
|
|
104 |
do i=1,ncum |
do i=1, ncum |
105 |
do k=1,nl |
do k=1, nl |
106 |
dtmin(i,k)=100.0 |
dtmin(i, k)=100.0 |
107 |
enddo |
enddo |
108 |
enddo |
enddo |
109 |
|
|
110 |
do 550 i=1,ncum |
do i=1, ncum |
111 |
do 560 k=1,nl |
do k=1, nl |
112 |
do 570 j=minorig,nl |
do j=minorig, nl |
113 |
if ( (k.ge.(icb(i)+1)).and.(k.le.inb(i)).and. & |
if ((k >= (icb(i) + 1)).and.(k <= inb(i)).and. & |
114 |
(j.ge.icb(i)).and.(j.le.(k-1)) )then |
(j >= icb(i)).and.(j <= (k-1)))then |
115 |
dtmin(i,k)=AMIN1(dtmin(i,k),buoy(i,j)) |
dtmin(i, k)=AMIN1(dtmin(i, k), buoy(i, j)) |
116 |
|
endif |
117 |
|
end do |
118 |
|
end do |
119 |
|
end do |
120 |
|
|
121 |
|
! The interval on which cape is computed starts at pbase: |
122 |
|
|
123 |
|
do k=1, nl |
124 |
|
do i=1, ncum |
125 |
|
if ((k >= (icb(i) + 1)).and.(k <= inb(i))) then |
126 |
|
deltap = MIN(pbase(i), ph(i, k-1))-MIN(pbase(i), ph(i, k)) |
127 |
|
cape(i)=cape(i) + rrd*buoy(i, k-1)*deltap/p(i, k-1) |
128 |
|
cape(i)=AMAX1(0.0, cape(i)) |
129 |
|
sigold(i, k)=sig(i, k) |
130 |
|
|
131 |
|
sig(i, k)=beta*sig(i, k) + alpha*dtmin(i, k)*ABS(dtmin(i, k)) |
132 |
|
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 |
|
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 |
|
w0(i, k)=w |
139 |
endif |
endif |
|
570 continue |
|
|
560 continue |
|
|
550 continue |
|
|
|
|
|
! the interval on which cape is computed starts at pbase : |
|
|
|
|
|
do 600 k=1,nl |
|
|
do 610 i=1,ncum |
|
|
|
|
|
if ((k.ge.(icb(i)+1)).and.(k.le.inb(i))) then |
|
|
|
|
|
deltap = MIN(pbase(i),ph(i,k-1))-MIN(pbase(i),ph(i,k)) |
|
|
cape(i)=cape(i)+rrd*buoy(i,k-1)*deltap/p(i,k-1) |
|
|
cape(i)=AMAX1(0.0,cape(i)) |
|
|
sigold(i,k)=sig(i,k) |
|
|
|
|
|
! dtmin(i,k)=100.0 |
|
|
! do 97 j=icb(i),k-1 ! mauvaise vectorisation |
|
|
! dtmin(i,k)=AMIN1(dtmin(i,k),buoy(i,j)) |
|
|
! 97 continue |
|
|
|
|
|
sig(i,k)=beta*sig(i,k)+alpha*dtmin(i,k)*ABS(dtmin(i,k)) |
|
|
sig(i,k)=amax1(sig(i,k),0.0) |
|
|
sig(i,k)=amin1(sig(i,k),0.01) |
|
|
fac=AMIN1(((dtcrit-dtmin(i,k))/dtcrit),1.0) |
|
|
w=(1.-beta)*fac*SQRT(cape(i))+beta*w0(i,k) |
|
|
amu=0.5*(sig(i,k)+sigold(i,k))*w |
|
|
m(i,k)=amu*0.007*p(i,k)*(ph(i,k)-ph(i,k+1))/tv(i,k) |
|
|
w0(i,k)=w |
|
|
endif |
|
|
|
|
|
610 continue |
|
|
600 continue |
|
|
|
|
|
do 700 i=1,ncum |
|
|
w0(i,icb(i))=0.5*w0(i,icb(i)+1) |
|
|
m(i,icb(i))=0.5*m(i,icb(i)+1) & |
|
|
*(ph(i,icb(i))-ph(i,icb(i)+1)) & |
|
|
/(ph(i,icb(i)+1)-ph(i,icb(i)+2)) |
|
|
sig(i,icb(i))=sig(i,icb(i)+1) |
|
|
sig(i,icb(i)-1)=sig(i,icb(i)) |
|
|
700 continue |
|
|
|
|
|
|
|
|
!! cape=0.0 |
|
|
!! do 98 i=icb+1,inb |
|
|
!! deltap = min(pbase,ph(i-1))-min(pbase,ph(i)) |
|
|
!! cape=cape+rrd*buoy(i-1)*deltap/p(i-1) |
|
|
!! dcape=rrd*buoy(i-1)*deltap/p(i-1) |
|
|
!! dlnp=deltap/p(i-1) |
|
|
!! cape=amax1(0.0,cape) |
|
|
!! sigold=sig(i) |
|
|
|
|
|
!! dtmin=100.0 |
|
|
!! do 97 j=icb,i-1 |
|
|
!! dtmin=amin1(dtmin,buoy(j)) |
|
|
!! 97 continue |
|
|
|
|
|
!! sig(i)=beta*sig(i)+alpha*dtmin*abs(dtmin) |
|
|
!! sig(i)=amax1(sig(i),0.0) |
|
|
!! sig(i)=amin1(sig(i),0.01) |
|
|
!! fac=amin1(((dtcrit-dtmin)/dtcrit),1.0) |
|
|
!! w=(1.-beta)*fac*sqrt(cape)+beta*w0(i) |
|
|
!! amu=0.5*(sig(i)+sigold)*w |
|
|
!! m(i)=amu*0.007*p(i)*(ph(i)-ph(i+1))/tv(i) |
|
|
!! w0(i)=w |
|
|
!! 98 continue |
|
|
!! w0(icb)=0.5*w0(icb+1) |
|
|
!! m(icb)=0.5*m(icb+1)*(ph(icb)-ph(icb+1))/(ph(icb+1)-ph(icb+2)) |
|
|
!! sig(icb)=sig(icb+1) |
|
|
!! sig(icb-1)=sig(icb) |
|
140 |
|
|
141 |
return |
end do |
142 |
end |
end do |
143 |
|
|
144 |
|
do i=1, ncum |
145 |
|
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 |
|
sig(i, icb(i)-1)=sig(i, icb(i)) |
151 |
|
end do |
152 |
|
|
153 |
|
end SUBROUTINE cv30_closure |
154 |
|
|
155 |
|
end module cv30_closure_m |