1 |
! $Header: $ |
2 |
! $Name: $ |
3 |
|
4 |
module lscale_cond_mod |
5 |
|
6 |
!----------------------------------------------------------------------- |
7 |
!use fms_mod, only: file_exist, error_mesg, open_file, & |
8 |
! check_nml_error, mpp_pe, FATAL, & |
9 |
! close_file |
10 |
use simple_sat_vapor_pres_mod, only: escomp, descomp |
11 |
|
12 |
use constants_mod, only: HLv,HLs,Cp_air,Grav,rdgas,rvgas |
13 |
use gcm_params_mod, only: gcm_LEN_MBUF, gcm_SQZ_R, gcm_stdMsgUnit |
14 |
|
15 |
implicit none |
16 |
private |
17 |
!----------------------------------------------------------------------- |
18 |
! ---- public interfaces ---- |
19 |
|
20 |
public lscale_cond, lscale_cond_init |
21 |
|
22 |
!----------------------------------------------------------------------- |
23 |
! ---- version number ---- |
24 |
|
25 |
character(len=128) :: version = '$Id: lscale_cond_mod.F90,v 1.2 2013/02/06 04:07:52 jmc Exp $' |
26 |
character(len=128) :: tag = '$Name: $' |
27 |
|
28 |
!----------------------------------------------------------------------- |
29 |
! ---- local/private data ---- |
30 |
|
31 |
real, parameter :: d622 = rdgas/rvgas |
32 |
real, parameter :: d378 = 1.-d622 |
33 |
|
34 |
logical :: do_init=.true. |
35 |
|
36 |
!----------------------------------------------------------------------- |
37 |
! --- namelist ---- |
38 |
|
39 |
real :: hc=1.00 |
40 |
logical :: do_evap=.false. |
41 |
|
42 |
namelist /lscale_cond_nml/ hc, do_evap |
43 |
|
44 |
!----------------------------------------------------------------------- |
45 |
! description of namelist variables |
46 |
! |
47 |
! hc = relative humidity at which large scale condensation |
48 |
! occurs, where 0 <= hc <= 1 (default: hc=1.) |
49 |
! |
50 |
! do_evap = flag for the re-evaporation of moisture in |
51 |
! sub-saturated layers below, if do_evap=.true. then |
52 |
! re-evaporation is performed (default: do_evap=.false.) |
53 |
! |
54 |
!----------------------------------------------------------------------- |
55 |
|
56 |
contains |
57 |
|
58 |
!####################################################################### |
59 |
|
60 |
! subroutine lscale_cond (tin, qin, pfull, phalf, coldT, & |
61 |
! rain, snow, tdel, qdel, mask, conv) |
62 |
subroutine lscale_cond (tin, qin, pfull, phalf, coldT, & |
63 |
rain, snow, tdel, qdel, qsat, & |
64 |
myThid, mask, conv ) |
65 |
|
66 |
|
67 |
!----------------------------------------------------------------------- |
68 |
! |
69 |
! large scale condensation |
70 |
! |
71 |
!----------------------------------------------------------------------- |
72 |
! |
73 |
! input: tin temperature at full model levels |
74 |
! qin specific humidity of water vapor at full |
75 |
! model levels |
76 |
! pfull pressure at full model levels |
77 |
! phalf pressure at half (interface) model levels |
78 |
! coldT should precipitation be snow at this point? |
79 |
! optional: |
80 |
! mask optional mask (0 or 1.) |
81 |
! conv logical flag; if true then no large-scale |
82 |
! adjustment is performed at that grid-point or |
83 |
! model level |
84 |
! |
85 |
! output: rain liquid precipitation (kg/m2) |
86 |
! snow frozen precipitation (kg/m2) |
87 |
! tdel temperature tendency at full model levels |
88 |
! qdel specific humidity tendency (of water vapor) at |
89 |
! full model levels |
90 |
! qsat saturated specific humidity |
91 |
! |
92 |
!----------------------------------------------------------------------- |
93 |
!--------------------- interface arguments ----------------------------- |
94 |
|
95 |
real , intent(in) , dimension(:,:,:) :: tin, qin, pfull, phalf |
96 |
logical , intent(in) , dimension(:,:):: coldT |
97 |
real , intent(out), dimension(:,:) :: rain,snow |
98 |
real , intent(out), dimension(:,:,:) :: tdel, qdel, qsat |
99 |
integer, intent(in) :: myThid |
100 |
real , intent(in) , dimension(:,:,:), optional :: mask |
101 |
logical, intent(in) , dimension(:,:,:), optional :: conv |
102 |
!----------------------------------------------------------------------- |
103 |
!---------------------- local data ------------------------------------- |
104 |
|
105 |
logical,dimension(size(tin,1),size(tin,2),size(tin,3)) :: do_adjust |
106 |
real,dimension(size(tin,1),size(tin,2),size(tin,3)) :: & |
107 |
esat, desat, dqsat, pmes, pmass |
108 |
real,dimension(size(tin,1),size(tin,2)) :: hlcp, precip |
109 |
integer k, kx |
110 |
!----------------------------------------------------------------------- |
111 |
! computation of precipitation by condensation processes |
112 |
!----------------------------------------------------------------------- |
113 |
|
114 |
! if (do_init) call error_mesg ('lscale_cond', & |
115 |
! 'lscale_cond_init has not been called.', FATAL) |
116 |
|
117 |
kx=size(tin,3) |
118 |
|
119 |
!----- compute proper latent heat -------------------------------------- |
120 |
WHERE (coldT) |
121 |
hlcp = HLs/Cp_air |
122 |
ELSEWHERE |
123 |
hlcp = HLv/Cp_air |
124 |
END WHERE |
125 |
|
126 |
!----- saturation vapor pressure (esat) & specific humidity (qsat) ----- |
127 |
|
128 |
call escomp (tin,esat) |
129 |
call descomp (tin,desat) |
130 |
|
131 |
|
132 |
esat(:,:,:)=esat(:,:,:)*hc |
133 |
|
134 |
where (pfull(:,:,:) > d378*esat(:,:,:)) |
135 |
pmes(:,:,:)=1.0/(pfull(:,:,:)-d378*esat(:,:,:)) |
136 |
qsat(:,:,:)=d622*esat(:,:,:)*pmes(:,:,:) |
137 |
qsat(:,:,:)=max(0.0,qsat(:,:,:)) |
138 |
dqsat(:,:,:)=d622*pfull(:,:,:)*desat(:,:,:)*pmes(:,:,:)*pmes(:,:,:) |
139 |
elsewhere |
140 |
pmes(:,:,:)=0.0 |
141 |
qsat(:,:,:)=0.0 |
142 |
dqsat(:,:,:)=0.0 |
143 |
endwhere |
144 |
|
145 |
!--------- do adjustment where greater than saturated value ------------ |
146 |
|
147 |
if (present(conv)) then |
148 |
!!!! do_adjust(:,:,:)=(.not.conv(:,:,:) .and. qin(:,:,:) > qsat(:,:,:)) |
149 |
do_adjust(:,:,:)=(.not.conv(:,:,:) .and. & |
150 |
(qin(:,:,:) - qsat(:,:,:))*qsat(:,:,:) > 0.0) |
151 |
else |
152 |
!!!! do_adjust(:,:,:)=(qin(:,:,:) > qsat(:,:,:)) |
153 |
do_adjust(:,:,:)=( (qin(:,:,:) - qsat(:,:,:))*qsat(:,:,:) > 0.0) |
154 |
endif |
155 |
|
156 |
if (present(mask)) then |
157 |
do_adjust(:,:,:)=do_adjust(:,:,:) .and. (mask(:,:,:) > 0.5) |
158 |
end if |
159 |
|
160 |
!----------- compute adjustments to temp and spec humidity ------------- |
161 |
do k = 1,kx |
162 |
where (do_adjust(:,:,k)) |
163 |
qdel(:,:,k)=(qsat(:,:,k)-qin(:,:,k))/(1.0+hlcp(:,:)*dqsat(:,:,k)) |
164 |
tdel(:,:,k)=-hlcp(:,:)*qdel(:,:,k) |
165 |
elsewhere |
166 |
qdel(:,:,k)=0.0 |
167 |
tdel(:,:,k)=0.0 |
168 |
endwhere |
169 |
end do |
170 |
|
171 |
|
172 |
!------------ pressure mass of each layer ------------------------------ |
173 |
|
174 |
do k=1,kx |
175 |
pmass(:,:,k)=(phalf(:,:,k+1)-phalf(:,:,k))/Grav |
176 |
enddo |
177 |
|
178 |
!------------ re-evaporation of precipitation in dry layer below ------- |
179 |
|
180 |
if (do_evap) then |
181 |
if (present(mask)) then |
182 |
! call precip_evap (pmass,tin,qin,qsat,dqsat,hlcp,tdel,qdel,mask) |
183 |
call precip_evap (pmass,tin,qin,qsat,dqsat,hlcp,tdel,qdel, & |
184 |
myThid, mask ) |
185 |
else |
186 |
call precip_evap (pmass,tin,qin,qsat,dqsat,hlcp,tdel,qdel, & |
187 |
myThid ) |
188 |
endif |
189 |
endif |
190 |
|
191 |
!------------ integrate precip ----------------------------------------- |
192 |
|
193 |
precip(:,:)=0.0 |
194 |
do k=1,kx |
195 |
precip(:,:)=precip(:,:)-pmass(:,:,k)*qdel(:,:,k) |
196 |
enddo |
197 |
precip(:,:)=max(precip(:,:),0.0) |
198 |
|
199 |
!assign precip to snow or rain |
200 |
WHERE (coldT) |
201 |
snow = precip |
202 |
rain = 0. |
203 |
ELSEWHERE |
204 |
rain = precip |
205 |
snow = 0. |
206 |
END WHERE |
207 |
|
208 |
!----------------------------------------------------------------------- |
209 |
|
210 |
end subroutine lscale_cond |
211 |
|
212 |
!####################################################################### |
213 |
|
214 |
! subroutine precip_evap (pmass, tin, qin, qsat, dqsat, hlcp, & |
215 |
! tdel, qdel, mask) |
216 |
subroutine precip_evap (pmass, tin, qin, qsat, dqsat, hlcp, & |
217 |
tdel, qdel, myThid, mask ) |
218 |
|
219 |
!----------------------------------------------------------------------- |
220 |
! performs re-evaporation of falling precipitation |
221 |
!----------------------------------------------------------------------- |
222 |
real, intent(in), dimension(:,:,:) :: pmass, tin, qin, qsat, dqsat |
223 |
real, intent(in), dimension(:,:) :: hlcp |
224 |
real, intent(inout), dimension(:,:,:) :: tdel, qdel |
225 |
integer, intent(in) :: myThid |
226 |
real, intent(in), dimension(:,:,:), optional :: mask |
227 |
!----------------------------------------------------------------------- |
228 |
real, dimension(size(tin,1),size(tin,2)) :: exq, def |
229 |
|
230 |
integer k |
231 |
!----------------------------------------------------------------------- |
232 |
exq(:,:)=0.0 |
233 |
|
234 |
do k=1,size(tin,3) |
235 |
|
236 |
where (qdel(:,:,k) < 0.0) exq(:,:) = exq(:,:) - & |
237 |
qdel(:,:,k)*pmass(:,:,k) |
238 |
|
239 |
if (present(mask)) exq(:,:) = exq(:,:)*mask(:,:,k) |
240 |
|
241 |
! ---- evaporate precip where needed ------ |
242 |
|
243 |
where ( (qdel(:,:,k) >= 0.0) .and. (exq(:,:) > 0.0) ) |
244 |
exq(:,:) = exq(:,:) / pmass(:,:,k) |
245 |
def(:,:) = (qsat(:,:,k)-qin(:,:,k))/(1.+hlcp(:,:)*dqsat(:,:,k)) |
246 |
def(:,:) = min(max(def(:,:),0.0),exq(:,:)) |
247 |
qdel(:,:,k) = qdel(:,:,k) + def(:,:) |
248 |
tdel(:,:,k) = tdel(:,:,k) - def(:,:)*hlcp(:,:) |
249 |
exq(:,:) = (exq(:,:)-def(:,:))*pmass(:,:,k) |
250 |
endwhere |
251 |
|
252 |
enddo |
253 |
|
254 |
!----------------------------------------------------------------------- |
255 |
|
256 |
end subroutine precip_evap |
257 |
|
258 |
!####################################################################### |
259 |
|
260 |
subroutine lscale_cond_init (myThid) |
261 |
|
262 |
!----------------------------------------------------------------------- |
263 |
! |
264 |
! initialization for large scale condensation |
265 |
! |
266 |
!----------------------------------------------------------------------- |
267 |
|
268 |
integer, intent(in) :: myThid |
269 |
|
270 |
!----------------------------------------------------------------------- |
271 |
|
272 |
! integer unit,io,ierr |
273 |
integer :: iUnit |
274 |
CHARACTER*(gcm_LEN_MBUF) :: msgBuf |
275 |
|
276 |
|
277 |
!----------- read namelist --------------------------------------------- |
278 |
|
279 |
! _BARRIER |
280 |
! _BEGIN_MASTER(myThid) |
281 |
CALL BARRIER(myThid) |
282 |
IF ( myThid.EQ.1 ) THEN |
283 |
|
284 |
WRITE(msgBuf,'(A)') 'LSCALE_COND_INIT: opening data.atm_gray' |
285 |
CALL PRINT_MESSAGE( msgBuf, gcm_stdMsgUnit, gcm_SQZ_R, myThid ) |
286 |
CALL OPEN_COPY_DATA_FILE( & |
287 |
'data.atm_gray', 'LSCALE_COND_INIT', & |
288 |
iUnit, & |
289 |
myThid ) |
290 |
! Read parameters from open data file |
291 |
READ(UNIT=iUnit,NML=lscale_cond_nml) |
292 |
WRITE(msgBuf,'(A)') & |
293 |
'LSCALE_COND_INIT: finished reading data.atm_gray' |
294 |
CALL PRINT_MESSAGE( msgBuf, gcm_stdMsgUnit, gcm_SQZ_R, myThid ) |
295 |
! Close the open data file |
296 |
CLOSE(iUnit) |
297 |
|
298 |
! if (file_exist('input.nml')) then |
299 |
! unit = open_file (file='input.nml', action='read') |
300 |
! ierr=1; do while (ierr /= 0) |
301 |
! read (unit, nml=lscale_cond_nml, iostat=io, end=10) |
302 |
! ierr = check_nml_error (io,'lscale_cond_nml') |
303 |
! enddo |
304 |
! 10 call close_file (unit) |
305 |
! endif |
306 |
|
307 |
!---------- output namelist -------------------------------------------- |
308 |
|
309 |
! unit = open_file (file='logfile.out', action='append') |
310 |
! if ( mpp_pe() == 0 ) then |
311 |
! write (unit,'(/,80("="),/(a))') trim(version), trim(tag) |
312 |
! write (unit,nml=lscale_cond_nml) |
313 |
! endif |
314 |
! call close_file (unit) |
315 |
|
316 |
do_init=.false. |
317 |
|
318 |
ENDIF |
319 |
CALL BARRIER(myThid) |
320 |
|
321 |
end subroutine lscale_cond_init |
322 |
|
323 |
!####################################################################### |
324 |
|
325 |
end module lscale_cond_mod |
326 |
|