SURFEX  V8_0
Surfex V8_0 release
 All Classes Files Functions Variables
ch_dep_isba.F90
Go to the documentation of this file.
1 !SFX_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
2 !SFX_LIC This is part of the SURFEX software governed by the CeCILL-C licence
3 !SFX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
4 !SFX_LIC for details. version 1.
5 ! #########
6  SUBROUTINE ch_dep_isba (PUSTAR, PHU, PPSN, &
7  pveg, plai, psand, pclay, presa, &
8  prs, pz0, pta, ppa, ptrad, pno, prock, &
9  hsv, psoilrc_so2, psoilrc_o3, pdep )
10 !###########################################################
11 !!
12 !! PURPOSE
13 !! -------
14 !!
15 !! Compute dry deposition velocity for chemical species on nature area
16 !!
17 !! AUTHOR
18 !! ------
19 !! P.Tulet * Laboratoire d'Aerologie*
20 !!
21 !! MODIFICATIONS
22 !! -------------
23 !! Original 20/02/97
24 !! Modification 21/07/00 (Guenais/Tulet) add deposition on town and
25 !! vegetation class
26 !! Modification 18/01/01 (Solmon/Tulet) patch dry deposition
27 !! Modification 18/07/03 (Tulet) surface externalization
28 !! Modification 01/2004 (Tulet Masson) removes patch calculation
29 !! Modification 03/2006 (Le Moigne) pb in where test with some
30 !! compilation options
31 !!
32 !-------------------------------------------------------------------------------
33 !
34 !* 0. DECLARATIONS
35 ! ------------
36 !
37 USE modd_isba_par
38 USE modd_data_cover_par
39 USE modd_csts
40 USE modd_ch_isba, ONLY : xrcclayso2, xrcclayo3, xrcsandso2, xrcsando3, &
41  xrcsnowso2, xrcsnowo3, xlandrext
42 !
43 USE modd_ch_surf
44 USE modd_surf_par, ONLY: xundef
45 !
46 USE yomhook ,ONLY : lhook, dr_hook
47 USE parkind1 ,ONLY : jprb
48 !
49 IMPLICIT NONE
50 !
51 !* 0.1 Declarations of dummy arguments :
52 !
53 REAL, DIMENSION(:), INTENT(IN) :: pustar ! friction velocity
54 REAL, DIMENSION(:), INTENT(IN) :: phu ! soil humidity
55 REAL, DIMENSION(:), INTENT(IN) :: ppsn ! fraction of the grid covered
56  ! by snow
57 REAL, DIMENSION(:), INTENT(IN) :: prs ! stomatal resistance
58 REAL, DIMENSION(:), INTENT(IN) :: pz0 ! vegetation roughness length
59 REAL, DIMENSION(:), INTENT(IN) :: pveg ! vegetation fraction
60 REAL, DIMENSION(:), INTENT(IN) :: plai ! Leaf area index
61 REAL, DIMENSION(:,:), INTENT(IN) :: psand ! Sand fraction
62 REAL, DIMENSION(:,:), INTENT(IN) :: pclay ! Clay fraction
63 REAL, DIMENSION(:), INTENT(IN) :: presa ! aerodynamical resistance
64 REAL, DIMENSION(:), INTENT(IN) :: pta ! air temperature forcing (K)
65 REAL, DIMENSION(:), INTENT(IN) :: ppa ! surface atmospheric pressure
66 REAL, DIMENSION(:), INTENT(IN) :: ptrad ! radiative temperature (K)
67 REAL, DIMENSION(:), INTENT(IN) :: psoilrc_so2 ! bare soil resistance for SO2
68 REAL, DIMENSION(:), INTENT(IN) :: psoilrc_o3 ! bare soil resistance for O3
69 REAL, DIMENSION(:,:), INTENT(OUT) :: pdep ! deposition dry velocity (m/s)
70 REAL, DIMENSION(:), INTENT(IN) :: pno, prock ! fractions of bare soil, rock
71  CHARACTER(LEN=6), DIMENSION(:), INTENT(IN) :: hsv ! name of chemical
72  ! species
73 !
74 !* 0.2 Declarations of local variables :
75 !
76 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zdiffmolval
77 ! Molecular diffusivity
78 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zscmdt
79 ! Sc(:)hmidt number
80 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: znatrb
81 ! nature quasi-laminar resistances
82 !
83 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zhenryvalcor
84 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zstomrc
85 ! stomatal surface resistance
86 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zmesorc
87 ! mesophyl resistance
88 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zextrc
89 ! leaf uptake external surface resistance
90 !
91 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zsoilrc
92 ! bare soil surface resistance
93 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: znatrc
94 ! nature surface resistances where vegetation is
95 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zsnowrc
96 ! snow surface resistance
97 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zclayrc
98 ! clay surface resistance
99 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zsandrc
100 ! sand surface resistance
101 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zbarerc
102 ! nature surface resistances for bare soils
103 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zrockrc
104 ! nature surface resistances for rocks
105 !
106 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zres_vegtype
107 REAL , DIMENSION(SIZE(PTRAD,1),size(HSV,1)) :: zres_snowtype
108 !
109 ! final nature resistance by vegtype
110 REAL, DIMENSION(SIZE(PTRAD,1)) :: ztype1_sand, ztype1_clay, ztype1_snow ! Type soil 1
111 REAL, DIMENSION(SIZE(PTRAD,1)) :: zustar
112 REAL, DIMENSION(SIZE(PTRAD,1)) :: zdiffmolh2o
113 ! final nature resistance
114 REAL, DIMENSION(SIZE(PTRAD,1)) :: zlandext
115 ! computed Rext from Wesely tabulations (89)
116 REAL, DIMENSION(SIZE(PTRAD,1)) :: zincrc
117 ! in-canopy transport resistance
118 REAL, DIMENSION(SIZE(PTRAD,1)) :: zcoef1, zcoef2, zcoef3, zcoef4, zcoef5, zinv1, zinv2
119 REAL, DIMENSION(SIZE(PTRAD,1)) :: ztcor
120 !
121 REAL, DIMENSION(size(HSV,1)) :: zvar1, zvar2, zfact1
122 !
123 REAL :: ztype2_sand, ztype2_clay, ztype2_snow ! Type soil 2
124 !
125 INTEGER :: jsv, ji
126 !
127 REAL(KIND=JPRB) :: zhook_handle
128 !
129 !============================================================================
130 !
131 ! Primilary
132 ! ---------
133 !Default values
134 !--------------
135 !
136 IF (lhook) CALL dr_hook('CH_DEP_ISBA',0,zhook_handle)
137 !
138 ! Default type soil
139 ! TYPE1 = RCCLAY(SO2) = 1000
140 ! TYPE2 = RCCLAY(O3) = 100
141 IF (xrcclayo3.NE.xundef) THEN
142  ztype2_clay = xrcclayo3
143 ELSE
144  ztype2_clay = 100.
145 ENDIF
146 !
147 ! TYPE1 = RCSAND(SO2) = 1000
148 ! TYPE2 = RCSAND(O3) = 200
149 IF (xrcsando3.NE.xundef) THEN
150  ztype2_sand = xrcsando3
151 ELSE
152  ztype2_sand = 200.
153 ENDIF
154 !
155 ! TYPE1 = RCSNOW(SO2)
156 ! TYPE2 = RCSNOW(O3)
157 IF (xrcsnowo3 /=xundef) THEN
158  ztype2_snow = xrcsnowo3
159 ELSE
160  ztype2_snow = 2000.
161 ENDIF
162 !
163 DO ji = 1, SIZE(pveg)
164  !
165  IF (xrcclayso2.NE.xundef) THEN
166  ztype1_clay(ji) = xrcclayso2
167  ELSE
168  ztype1_clay(ji) = 1000.
169  ENDIF
170  !
171  IF (xrcsandso2.NE.xundef) THEN
172  ztype1_sand(ji) = xrcsandso2
173  ELSE
174  ztype1_sand(ji) = 1000.
175  ENDIF
176  !
177  IF (xrcsnowso2/=xundef) THEN
178  ztype1_snow(ji) = xrcsnowso2
179  ELSEIF (ptrad(ji) > 275.) THEN
180  ztype1_snow(ji) = 540.
181  ELSE
182  ztype1_snow(ji) = 70. * (275. - ptrad(ji))
183  ENDIF
184  !
185  !
186  zustar(ji) = max(pustar(ji), 1e-9)
187  !
188  zcoef5(ji) = 1./(xkarman*zustar(ji))
189  !
190  ! 3.2.5 In-canopy transport resistance
191  ! ------------------------------
192  !
193  IF (pveg(ji) > 0.) THEN
194  zincrc(ji) = 14. * plai(ji) * 4. * pz0(ji) / zustar(ji)
195  ELSE
196  zincrc(ji) = 1e-4
197  ENDIF
198  !
199  !
200  ! computed Rext from Wesely tabulations (89)
201  !
202  IF ( xlandrext.NE.xundef ) THEN
203  ! user value
204  zlandext(ji) = xlandrext
205  ELSEIF (plai(ji) /= xundef) THEN
206  ! computed value
207  zlandext(ji) = 6000. - 4000. * tanh( 1.6 * (plai(ji) - 1.6) )
208  ELSE
209  zlandext(ji) = 9999.
210  END IF
211  !
212  !
213  zcoef1(ji) = 1./298. - 1./pta(ji)
214  !
215  zdiffmolh2o(ji) = 2.22e-05 + 1.46e-07 * (pta(ji) * (ppa(ji)/xp00)**(xrd/xcpd) - 273.)
216  zcoef2(ji) = prs(ji) * zdiffmolh2o(ji)
217  !
218  zcoef3(ji) = 1./zlandext(ji)
219  !
220  IF ( ptrad(ji) < 271.) THEN
221  zcoef4(ji) = 1000. * exp(-ptrad(ji) + 269.)
222  ELSE
223  zcoef4(ji) = 0.
224  ENDIF
225  !
226  ztcor(ji) = min(2.5e3, zcoef4(ji))
227  !
228  !
229  zinv1(ji) = 1.e-5/psoilrc_so2(ji)
230  !
231  zinv2(ji) = 1./psoilrc_o3(ji)
232  !
233 ENDDO
234 !
235 !
236 DO jsv = 1, SIZE(hsv,1)
237  !
238  zvar1(jsv) = xsrealreactval(jsv) / 3000.
239  zvar2(jsv) = xsrealreactval(jsv) * 100.
240  !
241  zfact1(jsv) = 1.46e-07 * sqrt(18. / xsrealmassmolval(jsv))
242  !
243 ENDDO
244 !
245 !============================================================================
246 !
247 DO jsv = 1, SIZE(hsv,1)
248  !
249  DO ji = 1, SIZE(pta)
250  !
251  ! 2.0 Quasi-laminar resistance
252  ! ------------------------
253  !
254  ! compute molecular diffusivity for each species (Langevin, 1905)
255  ! ----------------------------------------------
256  !
257  zdiffmolval(ji,jsv) = 2.22e-05 + (pta(ji) - 273.0) * zfact1(jsv)
258  !
259  ! computation of Rb for each cover type
260  ! -------------------------------------
261  !
262  zscmdt(ji,jsv) = 0.15e-4 / zdiffmolval(ji,jsv)
263  znatrb(ji,jsv) = ((zscmdt(ji,jsv)/0.72)**(2./3.)) * zcoef5(ji)
264  !
265  IF (plai(ji)/=xundef) znatrb(ji,jsv) = 2. * znatrb(ji,jsv)
266  !
267  ENDDO
268  !
269 ENDDO
270 !
271 DO jsv = 1, SIZE(hsv,1)
272  !
273  DO ji = 1, SIZE(pta)
274  !
275  !============================================================================
276  !
277  ! 3.0 Surface resistance on NATURE
278  ! --------------------------------
279  !
280  ! 3.0.1 Stomatal resistance
281  ! -------------------
282  !
283  !ZEXTRC_O3(:) = 1./(1./(3.*ZLANDEXT(:) + 1./3000.))
284  !
285  zhenryvalcor(ji,jsv) = xsrealhenryval(jsv,1) * exp(xsrealhenryval(jsv,2) * zcoef1(ji))
286  !
287  IF (prs(ji)>0.) THEN
288  !
289  zstomrc(ji,jsv) = zcoef2(ji) / zdiffmolval(ji,jsv)
290  !
291  ! 3.2.2 Mesophyl resistance
292  ! -------------------
293  !
294  zmesorc(ji,jsv) = 1. / ( zhenryvalcor(ji,jsv)/3000. + zvar2(jsv) )
295  !
296  ELSE
297  !
298  zstomrc(ji,jsv) = 9999.
299  !
300  zmesorc(ji,jsv) = 9999.
301  !
302  ENDIF
303  !
304  ! 3.2.4 External leaf uptake resistance (Wesely, 1989)
305  ! -------------------------------
306  !
307  IF (phu(ji) >= 1.) THEN ! for dew-wetted surface
308  !
309  ! compute Rext for any species exept O3
310  ! taking acount of (Walmsley, Wesely, 95, technical note, Atm Env vol 30)
311  zextrc(ji,jsv) = 1./( zcoef3(ji) + 1.0e-7*zhenryvalcor(ji,jsv) + zvar1(jsv) )
312  !
313  ELSEIF ( prs(ji) > 0. ) THEN
314  !
315  zextrc(ji,jsv) = zlandext(ji) / ( 1.0e-5 * zhenryvalcor(ji,jsv) + xsrealreactval(jsv) )
316  !
317  ELSE
318  !
319  zextrc(ji,jsv) = 9999.
320  !
321  ENDIF
322  !
323  ! IF ((HSV(JSV)=='O3').OR.(HSV(JSV)=='O_3')) &
324  ! ZEXTRC(:,JSV) = ZEXTRC_O3(:)*ZTCOR(:)
325  ! IF ((HSV(JSV)=='SO2').OR.(HSV(JSV)=='SO_2')) &
326  ! ZEXTRC(:,JSV) = ZTCOR(:)* 1./ (1/5000 +1./(3 *ZLANDEXT(:)))
327  !
328  ! Temperature correction
329  ! ----------------------
330  !
331  zextrc(ji,jsv) = zextrc(ji,jsv) + zcoef4(ji)
332  !
333  ENDDO
334  !
335 ENDDO
336 !
337 DO jsv = 1, SIZE(hsv,1)
338  !
339  DO ji = 1, SIZE(pta)
340  !
341  ! 3.2.6 Surface resistance on soil under veg
342  ! -------------------------------------
343  !
344  zsoilrc(ji,jsv) = 1. / ( zhenryvalcor(ji,jsv)*zinv1(ji) + xsrealreactval(jsv)*zinv2(ji) )
345  !
346  IF ( zstomrc(ji,jsv)>0. .AND. zincrc(ji)>0. .AND. zextrc(ji,jsv)>0. ) THEN
347  !
348  ! 3.2.7 Compute surface resistance on vegetation
349  ! -----------------------------------------
350  !
351  znatrc(ji,jsv) = 1./ &
352  ( 1./(zstomrc(ji,jsv)+zmesorc(ji,jsv)) + 1./(zincrc(ji)+zsoilrc(ji,jsv)) + 1./zextrc(ji,jsv) )
353  !
354  ELSE
355  znatrc(ji,jsv) = 1.e-4
356  ENDIF
357  !
358  ! 3.3 Surface resistance on NATURE with NO VEG (bare soil, rock, snow)
359  ! -----------------------------------------------------------------
360  !
361  ! 3.3.1 Surface resistance on clay
362  ! ---------------------------
363  !
364  zclayrc(ji,jsv) = ( 1.e5 * ztype1_clay(ji) * ztype2_clay ) / &
365  ( zhenryvalcor(ji,jsv)*ztype2_clay + ztype1_clay(ji)*1.e5*xsrealreactval(jsv) )
366  !
367  ! 3.3.2 Surface resistance on sand
368  ! ---------------------------
369  !
370  zsandrc(ji,jsv) = ( 1.e5 * ztype1_sand(ji) * ztype2_sand ) / &
371  ( zhenryvalcor(ji,jsv)*ztype2_sand + ztype1_sand(ji)*1.e5*xsrealreactval(jsv) )
372  !
373  ! 3.3.3 Compute surface resistance on bare soil
374  ! ---------------------------------------
375  !
376  zbarerc(ji,jsv) = 1./ ( psand(ji,1)/zsandrc(ji,jsv) + (1.-psand(ji,1))/zclayrc(ji,jsv) )
377  !
378  ! 3.3.4 Surface temperature correction
379  ! ------------------------------
380  !
381  zbarerc(ji,jsv) = zbarerc(ji,jsv) + ztcor(ji)
382  !
383  ! 3.3.5 Compute surface resistance on ROCK AREA
384  ! ---------------------------------------
385  !
386  zrockrc(ji,jsv) = ( 1.e5 * psoilrc_so2(ji) * psoilrc_o3(ji) ) / &
387  (zhenryvalcor(ji,jsv)*psoilrc_o3(ji) + psoilrc_so2(ji)*1.e5*xsrealreactval(jsv) )
388  !
389  ! 3.3.6 Surface temperature correction
390  ! ------------------------------
391  !
392  zrockrc(ji,jsv) = zrockrc(ji,jsv) + ztcor(ji)
393  !
394  ! 3.4 Surface resistance on snow
395  ! ----------------------------------
396  !
397  ! 3.4.1 Compute surface resistance on snow
398  ! ----------------------------------
399  !
400  zsnowrc(ji,jsv) = ( 1.e5 * ztype1_snow(ji) * ztype2_snow ) / &
401  ( zhenryvalcor(ji,jsv)*ztype2_snow + ztype1_snow(ji)*1.e5*xsrealreactval(jsv) )
402  !
403  ! 3.4.2 Surface temperature correction
404  ! ------------------------------
405  !
406  zsnowrc(ji,jsv) = zsnowrc(ji,jsv) + ztcor(ji)
407  !
408  ! 3.5 Surface resistance on snow (eternal or explicit)
409  ! --------------------------------------------
410  !
411  ! add rocks into bare soil resistance computation, when present
412  IF ( prock(ji)>0. ) THEN
413  zbarerc(ji,jsv) = ( pno(ji)+prock(ji) )/( pno(ji)/zbarerc(ji,jsv) + prock(ji)/zrockrc(ji,jsv) )
414  ENDIF
415  !
416  ! computes resistance due to soil and vegetation
417  znatrc(ji,jsv) = 1./ ( pveg(ji)/znatrc(ji,jsv) + (1.-pveg(ji))/zbarerc(ji,jsv) )
418  !
419  ENDDO
420  !
421 ENDDO
422 !
423 DO jsv = 1, SIZE(hsv,1)
424  !
425  DO ji = 1, SIZE(pta)
426  !
427  !---------------------------------------------------------------------
428  !
429  ! 4.0 Compute nature resistance
430  ! --------------------------
431  !
432  zres_vegtype(ji,jsv) = presa(ji) + znatrb(ji,jsv) + znatrc(ji,jsv)
433  zres_snowtype(ji,jsv) = presa(ji) + znatrb(ji,jsv) + zsnowrc(ji,jsv)
434  !
435  pdep(ji,jsv) = ( 1-ppsn(ji) )/zres_vegtype(ji,jsv) + ppsn(ji)/zres_snowtype(ji,jsv)
436  !
437  ENDDO
438  !
439 ENDDO
440 !
441 IF (lhook) CALL dr_hook('CH_DEP_ISBA',1,zhook_handle)
442 !
443 !---------------------------------------------------------------------
444 !
445 END SUBROUTINE ch_dep_isba
modd_isba_par
Definition: modd_isba_par.F90:6
modd_ch_surf
Definition: modd_ch_surf.F90:6
modd_data_cover_par
Definition: modd_data_cover_par.F90:6
ch_dep_isba
subroutine ch_dep_isba(PUSTAR, PHU, PPSN, PVEG, PLAI, PSAND, PCLAY, PRESA, PRS, PZ0, PTA, PPA, PTRAD, PNO, PROCK, HSV, PSOILRC_SO2, PSOILRC_O3, PDEP)
Definition: ch_dep_isba.F90:6
modd_surf_par
Definition: modd_surf_par.F90:6
modd_csts
Definition: modd_csts.F90:6
modd_ch_isba
Definition: modd_ch_isba.F90:6
Generated on Fri May 13 2016 09:31:00 for SURFEX by   doxygen 1.8.5