| 1 |
|
| 2 |
SUBROUTINE cv3_closure(nloc,ncum,nd,icb,inb & |
| 3 |
,pbase,p,ph,tv,buoy & |
| 4 |
,sig,w0,cape,m) |
| 5 |
use cvparam3 |
| 6 |
use cvthermo |
| 7 |
implicit none |
| 8 |
|
| 9 |
!=================================================================== |
| 10 |
! --- CLOSURE OF CONVECT3 |
| 11 |
! |
| 12 |
! vectorization: S. Bony |
| 13 |
!=================================================================== |
| 14 |
|
| 15 |
|
| 16 |
! input: |
| 17 |
integer ncum, nd, nloc |
| 18 |
integer icb(nloc), inb(nloc) |
| 19 |
real pbase(nloc) |
| 20 |
real p(nloc,nd), ph(nloc,nd+1) |
| 21 |
real tv(nloc,nd), buoy(nloc,nd) |
| 22 |
|
| 23 |
! input/output: |
| 24 |
real sig(nloc,nd), w0(nloc,nd) |
| 25 |
|
| 26 |
! output: |
| 27 |
real cape(nloc) |
| 28 |
real m(nloc,nd) |
| 29 |
|
| 30 |
! local variables: |
| 31 |
integer i, j, k, icbmax |
| 32 |
real deltap, fac, w, amu |
| 33 |
real dtmin(nloc,nd), sigold(nloc,nd) |
| 34 |
|
| 35 |
|
| 36 |
! ------------------------------------------------------- |
| 37 |
! -- Initialization |
| 38 |
! ------------------------------------------------------- |
| 39 |
|
| 40 |
do k=1,nl |
| 41 |
do i=1,ncum |
| 42 |
m(i,k)=0.0 |
| 43 |
enddo |
| 44 |
enddo |
| 45 |
|
| 46 |
! ------------------------------------------------------- |
| 47 |
! -- Reset sig(i) and w0(i) for i>inb and i<icb |
| 48 |
! ------------------------------------------------------- |
| 49 |
|
| 50 |
! update sig and w0 above LNB: |
| 51 |
|
| 52 |
do 100 k=1,nl-1 |
| 53 |
do 110 i=1,ncum |
| 54 |
if ((inb(i).lt.(nl-1)).and.(k.ge.(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 |
110 continue |
| 61 |
100 continue |
| 62 |
|
| 63 |
! compute icbmax: |
| 64 |
|
| 65 |
icbmax=2 |
| 66 |
do 200 i=1,ncum |
| 67 |
icbmax=MAX(icbmax,icb(i)) |
| 68 |
200 continue |
| 69 |
|
| 70 |
! update sig and w0 below cloud base: |
| 71 |
|
| 72 |
do 300 k=1,icbmax |
| 73 |
do 310 i=1,ncum |
| 74 |
if (k.le.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 |
310 continue |
| 80 |
300 continue |
| 81 |
|
| 82 |
!! if(inb.lt.(nl-1))then |
| 83 |
!! do 85 i=inb+1,nl-1 |
| 84 |
!! sig(i)=beta*sig(i)+2.*alpha*buoy(inb)* |
| 85 |
!! 1 abs(buoy(inb)) |
| 86 |
!! sig(i)=amax1(sig(i),0.0) |
| 87 |
!! w0(i)=beta*w0(i) |
| 88 |
!! 85 continue |
| 89 |
!! end if |
| 90 |
|
| 91 |
!! do 87 i=1,icb |
| 92 |
!! sig(i)=beta*sig(i)-2.*alpha*buoy(icb)*buoy(icb) |
| 93 |
!! sig(i)=amax1(sig(i),0.0) |
| 94 |
!! w0(i)=beta*w0(i) |
| 95 |
!! 87 continue |
| 96 |
|
| 97 |
! ------------------------------------------------------------- |
| 98 |
! -- Reset fractional areas of updrafts and w0 at initial time |
| 99 |
! -- and after 10 time steps of no convection |
| 100 |
! ------------------------------------------------------------- |
| 101 |
|
| 102 |
do 400 k=1,nl-1 |
| 103 |
do 410 i=1,ncum |
| 104 |
if (sig(i,nd).lt.1.5.or.sig(i,nd).gt.12.0)then |
| 105 |
sig(i,k)=0.0 |
| 106 |
w0(i,k)=0.0 |
| 107 |
endif |
| 108 |
410 continue |
| 109 |
400 continue |
| 110 |
|
| 111 |
! ------------------------------------------------------------- |
| 112 |
! -- Calculate convective available potential energy (cape), |
| 113 |
! -- vertical velocity (w), fractional area covered by |
| 114 |
! -- undilute updraft (sig), and updraft mass flux (m) |
| 115 |
! ------------------------------------------------------------- |
| 116 |
|
| 117 |
do 500 i=1,ncum |
| 118 |
cape(i)=0.0 |
| 119 |
500 continue |
| 120 |
|
| 121 |
! compute dtmin (minimum buoyancy between ICB and given level k): |
| 122 |
|
| 123 |
do i=1,ncum |
| 124 |
do k=1,nl |
| 125 |
dtmin(i,k)=100.0 |
| 126 |
enddo |
| 127 |
enddo |
| 128 |
|
| 129 |
do 550 i=1,ncum |
| 130 |
do 560 k=1,nl |
| 131 |
do 570 j=minorig,nl |
| 132 |
if ( (k.ge.(icb(i)+1)).and.(k.le.inb(i)).and. & |
| 133 |
(j.ge.icb(i)).and.(j.le.(k-1)) )then |
| 134 |
dtmin(i,k)=AMIN1(dtmin(i,k),buoy(i,j)) |
| 135 |
endif |
| 136 |
570 continue |
| 137 |
560 continue |
| 138 |
550 continue |
| 139 |
|
| 140 |
! the interval on which cape is computed starts at pbase : |
| 141 |
|
| 142 |
do 600 k=1,nl |
| 143 |
do 610 i=1,ncum |
| 144 |
|
| 145 |
if ((k.ge.(icb(i)+1)).and.(k.le.inb(i))) then |
| 146 |
|
| 147 |
deltap = MIN(pbase(i),ph(i,k-1))-MIN(pbase(i),ph(i,k)) |
| 148 |
cape(i)=cape(i)+rrd*buoy(i,k-1)*deltap/p(i,k-1) |
| 149 |
cape(i)=AMAX1(0.0,cape(i)) |
| 150 |
sigold(i,k)=sig(i,k) |
| 151 |
|
| 152 |
! dtmin(i,k)=100.0 |
| 153 |
! do 97 j=icb(i),k-1 ! mauvaise vectorisation |
| 154 |
! dtmin(i,k)=AMIN1(dtmin(i,k),buoy(i,j)) |
| 155 |
! 97 continue |
| 156 |
|
| 157 |
sig(i,k)=beta*sig(i,k)+alpha*dtmin(i,k)*ABS(dtmin(i,k)) |
| 158 |
sig(i,k)=amax1(sig(i,k),0.0) |
| 159 |
sig(i,k)=amin1(sig(i,k),0.01) |
| 160 |
fac=AMIN1(((dtcrit-dtmin(i,k))/dtcrit),1.0) |
| 161 |
w=(1.-beta)*fac*SQRT(cape(i))+beta*w0(i,k) |
| 162 |
amu=0.5*(sig(i,k)+sigold(i,k))*w |
| 163 |
m(i,k)=amu*0.007*p(i,k)*(ph(i,k)-ph(i,k+1))/tv(i,k) |
| 164 |
w0(i,k)=w |
| 165 |
endif |
| 166 |
|
| 167 |
610 continue |
| 168 |
600 continue |
| 169 |
|
| 170 |
do 700 i=1,ncum |
| 171 |
w0(i,icb(i))=0.5*w0(i,icb(i)+1) |
| 172 |
m(i,icb(i))=0.5*m(i,icb(i)+1) & |
| 173 |
*(ph(i,icb(i))-ph(i,icb(i)+1)) & |
| 174 |
/(ph(i,icb(i)+1)-ph(i,icb(i)+2)) |
| 175 |
sig(i,icb(i))=sig(i,icb(i)+1) |
| 176 |
sig(i,icb(i)-1)=sig(i,icb(i)) |
| 177 |
700 continue |
| 178 |
|
| 179 |
|
| 180 |
!! cape=0.0 |
| 181 |
!! do 98 i=icb+1,inb |
| 182 |
!! deltap = min(pbase,ph(i-1))-min(pbase,ph(i)) |
| 183 |
!! cape=cape+rrd*buoy(i-1)*deltap/p(i-1) |
| 184 |
!! dcape=rrd*buoy(i-1)*deltap/p(i-1) |
| 185 |
!! dlnp=deltap/p(i-1) |
| 186 |
!! cape=amax1(0.0,cape) |
| 187 |
!! sigold=sig(i) |
| 188 |
|
| 189 |
!! dtmin=100.0 |
| 190 |
!! do 97 j=icb,i-1 |
| 191 |
!! dtmin=amin1(dtmin,buoy(j)) |
| 192 |
!! 97 continue |
| 193 |
|
| 194 |
!! sig(i)=beta*sig(i)+alpha*dtmin*abs(dtmin) |
| 195 |
!! sig(i)=amax1(sig(i),0.0) |
| 196 |
!! sig(i)=amin1(sig(i),0.01) |
| 197 |
!! fac=amin1(((dtcrit-dtmin)/dtcrit),1.0) |
| 198 |
!! w=(1.-beta)*fac*sqrt(cape)+beta*w0(i) |
| 199 |
!! amu=0.5*(sig(i)+sigold)*w |
| 200 |
!! m(i)=amu*0.007*p(i)*(ph(i)-ph(i+1))/tv(i) |
| 201 |
!! w0(i)=w |
| 202 |
!! 98 continue |
| 203 |
!! w0(icb)=0.5*w0(icb+1) |
| 204 |
!! m(icb)=0.5*m(icb+1)*(ph(icb)-ph(icb+1))/(ph(icb+1)-ph(icb+2)) |
| 205 |
!! sig(icb)=sig(icb+1) |
| 206 |
!! sig(icb-1)=sig(icb) |
| 207 |
|
| 208 |
return |
| 209 |
end |
Parent Directory
| 
Revision Log