SURFEX  V8_0
Surfex V8_0 release
 All Classes Files Functions Variables
ecume_flux.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 ecume_flux(PZ0SEA,PTA,PEXNA,PRHOA,PSST,PEXNS,PQA,PVMOD, &
7  pzref,puref,pps,pichce,oprecip,opwebb,opwg,&
8  pqsat,psfth,psftq,pustar,pcd,pcdn,pch,pce, &
9  pri,presa,prain,pz0hsea,opertflux,ppertflux)
10 !###############################################################################
11 !!
12 !!**** *ECUME_FLUX*
13 !!
14 !! PURPOSE
15 !! -------
16 ! Calculate the surface turbulent fluxes of heat, moisture, and momentum
17 ! over sea surface + corrections due to rainfall & Webb effect.
18 !!
19 !!** METHOD
20 !! ------
21 ! The estimation of the transfer coefficients relies on the iterative
22 ! computation of the scaling parameters U*/Teta*/q*. The convergence is
23 ! supposed to be reached in NITERFL iterations maximum.
24 ! The neutral transfer coefficients for momentum/temperature/humidity
25 ! are computed as a function of the 10m-height neutral wind speed using
26 ! the ECUME_v0 formulation based on the multi-campaign (POMME,FETCH,CATCH,
27 ! SEMAPHORE,EQUALANT) ALBATROS dataset. See MERSEA report for more
28 ! details on the ECUME formulation.
29 !!
30 !! EXTERNAL
31 !! --------
32 !!
33 !! IMPLICIT ARGUMENTS
34 !! ------------------
35 !!
36 !! REFERENCE
37 !! ---------
38 !! Fairall et al (1996), JGR, 3747-3764
39 !! Gosnell et al (1995), JGR, 437-442
40 !! Fairall et al (1996), JGR, 1295-1308
41 !!
42 !! AUTHOR
43 !! ------
44 !! C. Lebeaupin *Météo-France* (adapted from S. Belamari's code)
45 !!
46 !! MODIFICATIONS
47 !! -------------
48 !! Original 15/03/2005
49 !! Modified 01/2006 C. Lebeaupin (adapted from A. Pirani's code)
50 !! Modified 20/07/2009 S. Belamari
51 !! Modified 08/2009 B. Decharme: limitation of Ri
52 !! Modified 09/2012 B. Decharme: CD correction
53 !! Modified 09/2012 B. Decharme: limitation of Ri in surface_ri.F90
54 !! Modified 10/2012 P. Le Moigne: extra inputs for FLake use
55 !! Modified 06/2013 B. Decharme: bug in z0 (output) computation
56 !! Modified 06/2013 J.Escobar : for REAL4/8 add EPSILON management
57 !!!
58 !-------------------------------------------------------------------------------
59 
60 ! 0. DECLARATIONS
61 ! ------------
62 !
63 USE modd_csts, ONLY : xkarman, xg, xstefan, xrd, xrv, &
64  xlvtt, xcl, xcpd, xcpv, xrholw, &
65  xtt,xp00
66 USE modd_surf_par, ONLY : xundef, xsurf_epsilon
67 !
68 USE modd_reprod_oper, ONLY : ccharnock
69 !
70 USE modd_surf_par, ONLY : xundef
71 USE modd_snow_par, ONLY : xz0sn, xz0hsn
72 USE modd_surf_atm, ONLY : xvchrnk, xvz0cm
73 !
74 USE modd_water_par
75 !
76 USE modi_wind_threshold
77 USE modi_surface_ri
78 !
79 USE mode_thermos
80 !
81 !
82 USE yomhook ,ONLY : lhook, dr_hook
83 USE parkind1 ,ONLY : jprb
84 !
85 IMPLICIT NONE
86 
87 ! 0.1. Declarations of arguments
88 
89 REAL, DIMENSION(:), INTENT(IN) :: pta ! air temperature, atm.lev (K)
90 REAL, DIMENSION(:), INTENT(IN) :: pqa ! air spec. hum., atm.lev (kg/kg)
91 REAL, DIMENSION(:), INTENT(IN) :: prhoa ! air density, atm.lev (kg/m3)
92 REAL, DIMENSION(:), INTENT(IN) :: pvmod ! module of wind, atm.lev (m/s)
93 REAL, DIMENSION(:), INTENT(IN) :: pzref ! atm.level for temp./hum. (m)
94 REAL, DIMENSION(:), INTENT(IN) :: puref ! atm.level for wind (m)
95 REAL, DIMENSION(:), INTENT(IN) :: psst ! Sea Surface Temperature (K)
96 REAL, DIMENSION(:), INTENT(IN) :: pps ! air pressure at sea surf. (Pa)
97 REAL, DIMENSION(:), INTENT(IN) :: prain ! precipitation rate (kg/s/m2)
98 REAL, DIMENSION(:), INTENT(IN) :: pexna ! Exner function at atm. level
99 REAL, DIMENSION(:), INTENT(IN) :: pexns ! Exner function at sea surface
100 REAL, DIMENSION(:), INTENT(IN) :: ppertflux ! stochastic flux perturbation pattern
101 
102 REAL, INTENT(IN) :: pichce !
103 LOGICAL, INTENT(IN) :: oprecip !
104 LOGICAL, INTENT(IN) :: opwebb !
105 LOGICAL, INTENT(IN) :: opwg !
106 LOGICAL, INTENT(IN) :: opertflux !
107 
108 REAL, DIMENSION(:), INTENT(INOUT) :: pz0sea ! roughness length over the ocean
109 
110 ! surface fluxes : latent heat, sensible heat, friction fluxes
111 REAL, DIMENSION(:), INTENT(OUT) :: psfth ! heat flux (W/m2)
112 REAL, DIMENSION(:), INTENT(OUT) :: psftq ! water flux (kg/m2/s)
113 REAL, DIMENSION(:), INTENT(OUT) :: pustar ! friction velocity (m/s)
114 
115 ! diagnostics
116 REAL, DIMENSION(:), INTENT(OUT) :: pqsat ! sea surface spec. hum. (kg/kg)
117 REAL, DIMENSION(:), INTENT(OUT) :: pcd ! transfer coef. for momentum
118 REAL, DIMENSION(:), INTENT(OUT) :: pch ! transfer coef. for temperature
119 REAL, DIMENSION(:), INTENT(OUT) :: pce ! transfer coef. for humidity
120 REAL, DIMENSION(:), INTENT(OUT) :: pcdn ! neutral coef. for momentum
121 REAL, DIMENSION(:), INTENT(OUT) :: presa ! aerodynamical resistance
122 REAL, DIMENSION(:), INTENT(OUT) :: pri ! Richardson number
123 REAL, DIMENSION(:), INTENT(OUT) :: pz0hsea ! heat roughness length
124 
125 ! 0.2. Declarations of local variables
126 
127 REAL, DIMENSION(SIZE(PTA)) :: ztau ! momentum flux (N/m2)
128 REAL, DIMENSION(SIZE(PTA)) :: zhf ! sensible heat flux (W/m2)
129 REAL, DIMENSION(SIZE(PTA)) :: zef ! latent heat flux (W/m2)
130 REAL, DIMENSION(SIZE(PTA)) :: ztaur ! momentum flx due to rain (N/m2)
131 REAL, DIMENSION(SIZE(PTA)) :: zrf ! sensible flx due to rain (W/m2)
132 REAL, DIMENSION(SIZE(PTA)) :: zefwebb ! Webb corr. on latent flx (W/m2)
133 
134 REAL, DIMENSION(SIZE(PTA)) :: zvmod ! wind intensity at atm.lev (m/s)
135 REAL, DIMENSION(SIZE(PTA)) :: zqsata ! sat.spec.hum., atm.lev (kg/kg)
136 REAL, DIMENSION(SIZE(PTA)) :: zpa ! air pressure at atm. level (Pa)
137 REAL, DIMENSION(SIZE(PTA)) :: zusr ! velocity scaling param. (m/s)
138  ! =friction velocity
139 REAL, DIMENSION(SIZE(PTA)) :: ztsr ! temperature scaling param. (K)
140 REAL, DIMENSION(SIZE(PTA)) :: zqsr ! humidity scaling param. (kg/kg)
141 REAL, DIMENSION(SIZE(PTA)) :: zwg ! gustiness factor (m/s)
142 
143 REAL, DIMENSION(SIZE(PTA)) :: zustar2 ! square of friction velocity
144 REAL, DIMENSION(SIZE(PTA)) :: zac ! aerodynamical conductance
145 REAL, DIMENSION(SIZE(PTA)) :: zdircoszw ! orography slope cosine
146  ! (=1 on water!)
147 
148 REAL, DIMENSION(SIZE(PTA)) :: zlv,zlr ! vap.heat, sea/atm level (J/kg)
149 REAL, DIMENSION(SIZE(PTA)) :: zdu,zdth,zdq,zduwg
150  ! vert. gradients (real atm.)
151 REAL, DIMENSION(SIZE(PTA)) :: zdeltau10n,zdeltat10n,zdeltaq10n
152  ! vert. gradients (10-m, neutral)
153 REAL, DIMENSION(SIZE(PTA)) :: zchn,zcen ! neutral coef. for T,Q
154 REAL, DIMENSION(SIZE(PTA)) :: zd0
155 REAL, DIMENSION(SIZE(PTA)) :: zcharn !Charnock number
156 REAL, DIMENSION(SIZE(PTA)) :: ztvsr,zbf ! constants to compute gustiness factor
157 !
158 REAL :: zetv,zrdsrv ! thermodynamic constants
159 REAL :: zlmou,zlmot ! Obukhovs stability param. z/l for U, T/Q
160 REAL :: zpsi_u,zpsi_t ! PSI funct. for U, T/Q
161 REAL :: zlmomin,zlmomax ! min/max value of Obukhovs stability parameter z/l
162 REAL :: zbetagust ! gustiness factor
163 REAL :: zzbl !atm. boundary layer depth (m)
164 !
165 REAL :: zchic,zchik,zeps,zloghs10,zlogts10,zpi,zpis2,zpsic,zpsik, &
166  zsqr3,zzdq,zzdth
167 !
168 REAL :: zalfac,zcplw,zdqsdt,zdtmp,zdwat,zp00,ztac,zww
169  ! to compute rainfall impact & Webb correction
170 !
171 INTEGER :: niterfl ! maximum number of iterations (5 or 6)
172 INTEGER :: jlon, jj
173 REAL(KIND=JPRB) :: zhook_handle
174 !
175 !-------------------------------------------------------------------------------
176 IF (lhook) CALL dr_hook('ECUME_FLUX',0,zhook_handle)
177 !
178 niterfl = 10
179 !
180 !-------------------------------------------------------------------------------
181 !
182 ! 1. AUXILIARY CONSTANTS & ARRAY INITIALISATION BY UNDEFINED VALUES.
183 ! --------------------------------------------------------------------
184 !
185 zlmomin = -200.0
186 zlmomax = 0.25
187 zp00 = 1013.25e+02
188 zpis2 = 2.0*atan(1.0)
189 zpi = 2.0*zpis2
190 zsqr3 = sqrt(3.0)
191 zeps = 1.e-8
192 zetv = xrv/xrd-1.0
193 zrdsrv = xrd/xrv
194 !
195 zbetagust=1.2 ! value based on TOGA-COARE experiment
196 zzbl =600. ! Set a default value for boundary layer depth
197 !
198 zdircoszw(:)=1.
199 !
200 pcd(:) = xundef
201 pch(:) = xundef
202 pce(:) = xundef
203 pcdn(:) = xundef
204 zusr(:) = xundef
205 ztsr(:) = xundef
206 zqsr(:) = xundef
207 ztau(:) = xundef
208 zhf(:) = xundef
209 zef(:) = xundef
210 !
211 psfth(:) = xundef
212 psftq(:) = xundef
213 pustar(:) = xundef
214 presa(:) = xundef
215 pri(:) = xundef
216 !
217 zwg(:) = 0.0
218 ztaur(:) = 0.0
219 zrf(:) = 0.0
220 zefwebb(:) = 0.0
221 !
222 !-------------------------------------------------------------------------------
223 !
224 ! 2. INITIALISATIONS BEFORE ITERATIVE LOOP.
225 ! -------------------------------------------
226 !
227 zvmod(:) = wind_threshold(pvmod(:),puref(:)) !set a minimum value to wind
228 !
229 ! 2.1. Specific humidity at saturation
230 !
231 pqsat(:) = qsat_seawater(psst(:),pps(:)) !at sea surface
232 zpa(:) = xp00*(pexna(:)**(xcpd/xrd))
233 zqsata(:) = qsat(pta(:),zpa(:)) !at atm. level
234 !
235 ! 2.2. Gradients at the air-sea interface
236 !
237 zdu(:) = zvmod(:) !one assumes u is measured / sea surface current
238 zdth(:) = pta(:)/pexna(:)-psst(:)/pexns(:)
239 zdq(:) = pqa(:)-pqsat(:)
240 !
241 ! 2.3. Initial guess
242 !
243 zd0(:) = 1.2+6.3e-03*max(zdu(:)-10.0,0.0)
244 !
245 IF (opwg) THEN !initial guess for gustiness factor
246  zwg(:) = 0.5
247  zduwg(:) = sqrt(zdu(:)**2+zwg(:)**2)
248 ELSE
249  zduwg(:) = zdu(:)
250 ENDIF
251 !
252 zdeltau10n(:) = zduwg(:)
253 zdeltat10n(:) = zdth(:)*zd0(:)
254 zdeltaq10n(:) = zdq(:)
255 !
256 ! 2.4. Latent heat of vaporisation
257 !
258 zlv(:) = xlvtt+(xcpv-xcl)*(psst(:)-xtt) !at sea surface
259 zlr(:) = xlvtt+(xcpv-xcl)*(pta(:)-xtt) !at atm.level
260 !
261 ! 2.5. Charnock number
262 !
263 IF(ccharnock=='OLD')THEN
264  zcharn(:) = xvchrnk
265 ELSE
266 ! vary between 0.011 et 0.018 according to Chris Fairall's data as in coare3.0
267  zcharn(:) = max(0.011,min(0.018,0.011+0.007*(zduwg(:)-10.)/8.))
268 ENDIF
269 !
270 !-------------------------------------------------------------------------------
271 !
272 ! 3. ITERATIVE LOOP TO COMPUTE U*, T*, Q*.
273 ! ------------------------------------------
274 !
275 DO jj=1,niterfl
276  DO jlon=1,SIZE(pta)
277 !
278 ! 3.1. Neutral coefficient for wind speed cdn (ECUME_v0 formulation)
279 !
280  IF (zdeltau10n(jlon) <= 16.8) THEN
281  pcdn(jlon) = 1.3013e-03 &
282  + (-1.2719e-04 * zdeltau10n(jlon) ) &
283  + (+1.3067e-05 * zdeltau10n(jlon)**2) &
284  + (-2.2261e-07 * zdeltau10n(jlon)**3)
285  ELSEIF (zdeltau10n(jlon) <= 50.0) THEN
286  pcdn(jlon) = 1.3633e-03 &
287  + (-1.3056e-04 * zdeltau10n(jlon) ) &
288  + (+1.6212e-05 * zdeltau10n(jlon)**2) &
289  + (-4.8208e-07 * zdeltau10n(jlon)**3) &
290  + (+4.2684e-09 * zdeltau10n(jlon)**4)
291  ELSE
292  pcdn(jlon) = 1.7828e-03
293  ENDIF
294 !
295 ! 3.2. Neutral coefficient for temperature chn (ECUME_v0 formulation)
296 !
297  IF (zdeltau10n(jlon) <= 33.0) THEN
298  zchn(jlon) = 1.2536e-03 &
299  + (-1.2455e-04 * zdeltau10n(jlon) ) &
300  + (+1.6038e-05 * zdeltau10n(jlon)**2) &
301  + (-4.3701e-07 * zdeltau10n(jlon)**3) &
302  + (+3.4517e-09 * zdeltau10n(jlon)**4) &
303  + (+3.5763e-12 * zdeltau10n(jlon)**5)
304  ELSE
305  zchn(jlon) = 3.1374e-03
306  ENDIF
307 !
308 ! 3.3. Neutral coefficient for humidity cen (ECUME_v0 formulation)
309 !
310  IF (zdeltau10n(jlon) <= 29.0) THEN
311  zcen(jlon) = 1.2687e-03 &
312  + (-1.1384e-04 * zdeltau10n(jlon) ) &
313  + (+1.1467e-05 * zdeltau10n(jlon)**2) &
314  + (-3.9144e-07 * zdeltau10n(jlon)**3) &
315  + (+5.0864e-09 * zdeltau10n(jlon)**4)
316  ELSEIF (zdeltau10n(jlon) <= 33.0) THEN
317  zcen(jlon) = -1.3526e-03 &
318  + (+1.8229e-04 * zdeltau10n(jlon) ) &
319  + (-2.6995e-06 * zdeltau10n(jlon)**2)
320  ELSE
321  zcen(jlon) = 1.7232e-03
322  ENDIF
323  zcen(jlon) = zcen(jlon)*(1.0-pichce)+zchn(jlon)*pichce
324 !
325 ! 3.4. Scaling parameters and roughness lenght
326 !
327  zusr(jlon) = sqrt(pcdn(jlon))*zdeltau10n(jlon)
328  ztsr(jlon) = zchn(jlon)/sqrt(pcdn(jlon))*zdeltat10n(jlon)
329  zqsr(jlon) = zcen(jlon)/sqrt(pcdn(jlon))*zdeltaq10n(jlon)
330 !
331 ! 3.5. Gustiness factor ZWG following Mondon & Redelsperger (1998)
332 !
333  IF(opwg) THEN
334  ztvsr(jlon)=ztsr(jlon)*(1.0+zetv*pqa(jlon))+zetv*pta(jlon)*zqsr(jlon)
335  zbf(jlon)=max(0.0,-xg/pta(jlon)*zusr(jlon)*ztvsr(jlon))
336  zwg(jlon)=zbetagust*(zbf(jlon)*zzbl)**(1./3.)
337  ENDIF
338 !
339 ! 3.6. Obukhovs stability param. z/l following Liu et al. (JAS, 1979)
340 !
341 ! For U
342  zlmou = puref(jlon)*xg*xkarman*(ztsr(jlon)/(pta(jlon)) &
343  +zetv*zqsr(jlon)/(1.0+zetv*pqa(jlon)))/max(zusr(jlon),zeps)**2
344 ! For T/Q
345  zlmot = zlmou*pzref(jlon)/puref(jlon)
346  zlmou = max(min(zlmou,zlmomax),zlmomin)
347  zlmot = max(min(zlmot,zlmomax),zlmomin)
348 !
349 ! 3.7. Stability function psi (see Liu et al, 1979 ; Dyer and Hicks, 1970)
350 ! Modified to include convective form following Fairall (unpublished)
351 !
352 ! For U
353  IF (zlmou == 0.0) THEN
354  zpsi_u = 0.0
355  ELSEIF (zlmou > 0.0) THEN
356  zpsi_u = -7.0*zlmou
357  ELSE
358  zchik = (1.0-16.0*zlmou)**0.25
359  zpsik = 2.0*log((1.0+zchik)/2.0) &
360  +log((1.0+zchik**2)/2.0) &
361  -2.0*atan(zchik)+zpis2
362  zchic = (1.0-12.87*zlmou)**(1.0/3.0) !for very unstable conditions
363  zpsic = 1.5*log((zchic**2+zchic+1.0)/3.0) &
364  -zsqr3*atan((2.0*zchic+1.0)/zsqr3) &
365  +zpi/zsqr3
366  zpsi_u = zpsic+(zpsik-zpsic)/(1.0+zlmou**2)
367  !match Kansas & free-conv. forms
368  ENDIF
369 ! For T/Q
370  IF (zlmot == 0.0) THEN
371  zpsi_t = 0.0
372  ELSEIF (zlmot > 0.0) THEN
373  zpsi_t = -7.0*zlmot
374  ELSE
375  zchik = (1.0-16.0*zlmot)**0.25
376  zpsik = 2.0*log((1.0+zchik**2)/2.0)
377  zchic = (1.0-12.87*zlmot)**(1.0/3.0) !for very unstable conditions
378  zpsic = 1.5*log((zchic**2+zchic+1.0)/3.0) &
379  -zsqr3*atan((2.0*zchic+1.0)/zsqr3) &
380  +zpi/zsqr3
381  zpsi_t = zpsic+(zpsik-zpsic)/(1.0+zlmot**2)
382  !match Kansas & free-conv. forms
383  ENDIF
384 !
385 ! 3.8. Update ZDELTAU10N, ZDELTAT10N and ZDELTAQ10N
386 !
387  zloghs10 = log(puref(jlon)/10.0)
388  zlogts10 = log(pzref(jlon)/10.0)
389  zduwg(jlon) = sqrt(zdu(jlon)**2+zwg(jlon)**2)
390  zdeltau10n(jlon) = zduwg(jlon)-zusr(jlon)*(zloghs10-zpsi_u)/xkarman
391  zdeltat10n(jlon) = zdth(jlon)-ztsr(jlon)*(zlogts10-zpsi_t)/xkarman
392  zdeltaq10n(jlon) = zdq(jlon)-zqsr(jlon)*(zlogts10-zpsi_t)/xkarman
393 
394  ENDDO
395 ENDDO
396 !
397 !-------------------------------------------------------------------------------
398 !
399 ! 4. COMPUTATION OF EXCHANGE COEFFICIENTS AND TURBULENT FLUXES.
400 ! ---------------------------------------------------------------
401 !
402 DO jlon=1,SIZE(pta)
403 !
404 ! 4.1. Exchange coefficients PCD, PCH, PCE
405 !
406  zzdth = 0.5* &
407  ((1.0+sign(1.0,zdth(jlon)))*max(zdth(jlon),zeps) &
408  +(1.0-sign(1.0,zdth(jlon)))*min(zdth(jlon),-zeps))
409  zzdq = 0.5* &
410  ((1.0+sign(1.0,zdq(jlon)))*max(zdq(jlon),zeps) &
411  +(1.0-sign(1.0,zdq(jlon)))*min(zdq(jlon),-zeps))
412  pcd(jlon) = (zusr(jlon)/zduwg(jlon))**2
413  pch(jlon) = zusr(jlon)*ztsr(jlon)/(zduwg(jlon)*zzdth)
414  pce(jlon) = zusr(jlon)*zqsr(jlon)/(zduwg(jlon)*zzdq)
415 !
416 ! 4.2. Surface turbulent fluxes
417 ! (ATM CONV.: ZTAU<<0 ; ZHF,ZEF<0 if atm looses heat)
418 !
419  ztau(jlon) = -prhoa(jlon)*pcd(jlon)*zduwg(jlon)**2
420  zhf(jlon) = -prhoa(jlon)*xcpd*pch(jlon)*zduwg(jlon)*zdth(jlon)
421  zef(jlon) = -prhoa(jlon)*zlv(jlon)*pce(jlon)*zduwg(jlon)*zdq(jlon)
422 !
423 ! 4.3. Stochastic perturbation of turbulent fluxes
424 
425  IF( opertflux )THEN
426  ztau(jlon) = ztau(jlon)* ( 1. + ppertflux(jlon) / 2. )
427  zhf(jlon) = zhf(jlon)* ( 1. + ppertflux(jlon) / 2. )
428  zef(jlon) = zef(jlon)* ( 1. + ppertflux(jlon) / 2. )
429  ENDIF
430 !
431 ENDDO
432 !
433 !-------------------------------------------------------------------------------
434 !
435 ! 5. COMPUTATION OF FLUX CORRECTIONS DUE TO RAINFALL.
436 ! (ATM conv: ZRF<0 if atm. looses heat, ZTAUR<<0)
437 ! -----------------------------------------------------
438 
439 IF(oprecip) THEN
440  DO jlon=1,SIZE(pta)
441 !
442 ! 5.1. Momentum flux due to rainfall (ZTAUR, N/m2)
443 !
444 ! See pp3752 in FBR96.
445  ztaur(jlon) = -prain(jlon)*zduwg(jlon)
446 !
447 ! 5.2. Sensible heat flux due to rainfall (ZRF, W/m2)
448 !
449 ! See Eq.12 in GoF95, with ZCPLW as specific heat of water (J/kg/K), ZDWAT as
450 ! water vapor diffusivity (Eq.13-3 of Pruppacher and Klett, 1978), ZDTMP as
451 ! heat diffusivity, ZDQSDT from Clausius-Clapeyron relation and ZALFAC as
452 ! wet-bulb factor (Eq.11 in GoF95).
453  ztac = pta(jlon)-xtt
454  zcplw = 4224.8482+ztac*(-4.707+ztac*(0.08499 &
455  +ztac*(1.2826e-03+ztac*(4.7884e-05 &
456  -2.0027e-06*ztac))))
457  zdwat = 2.11e-05*(zp00/zpa(jlon)) &
458  *(pta(jlon)/xtt)**1.94
459  zdtmp = (1.0+3.309e-03*ztac-1.44e-06*ztac**2) &
460  *0.02411/(prhoa(jlon)*xcpd)
461  zdqsdt = zqsata(jlon)*zlr(jlon)/(xrd*pta(jlon)**2)
462  zalfac = 1.0/(1.0+zdqsdt*(zlr(jlon)*zdwat)/(zdtmp*xcpd))
463  zrf(jlon) = zcplw*prain(jlon)*zalfac*((psst(jlon)-pta(jlon)) &
464  +(pqsat(jlon)-pqa(jlon))*(zlr(jlon)*zdwat)/(zdtmp*xcpd))
465 
466  ENDDO
467 ENDIF
468 !
469 !-------------------------------------------------------------------------------
470 !
471 ! 6. COMPUTATION OF WEBB CORRECTION TO LATENT HEAT FLUX (ZEFWEBB, W/m2).
472 ! ------------------------------------------------------------------------
473 !
474 ! See Eq.21 and Eq.22 in FBR96.
475 IF (opwebb) THEN
476  DO jlon=1,SIZE(pta)
477  zww = -(1.0+zetv)*(pce(jlon)*zduwg(jlon)*zdq(jlon)) &
478  -(1.0+(1.0+zetv)*pqa(jlon))* &
479  (pch(jlon)*zduwg(jlon)*zdth(jlon))/(pta(jlon))
480  zefwebb(jlon) = prhoa(jlon)*zww*zlv(jlon)*pqa(jlon)
481  ENDDO
482 ENDIF
483 !
484 !-------------------------------------------------------------------------------
485 !
486 ! 7. FINAL STEP : TOTAL SURFACE FLUXES AND DERIVED DIAGNOSTICS.
487 ! ---------------------------------------------------------------
488 !
489 ! 7.1. Richardson number
490 !
491  CALL surface_ri(psst,pqsat,pexns,pexna,pta,zqsata, &
492  pzref,puref,zdircoszw,pvmod,pri )
493 !
494 ! 7.2. Friction velocity which contains correction du to rain
495 !
496 zustar2(:)=-(ztau(:)+ztaur(:))/prhoa(:)
497 !
498 IF(oprecip) THEN
499  pcd(:)=zustar2(:)/(zduwg(:)**2)
500 ENDIF
501 !
502 pustar(:)=sqrt(zustar2(:))
503 !
504 ! 7.3. Aerodynamical conductance and resistance
505 !
506 zac(:)=pch(:)*zvmod(:)
507 presa(:)=1./ max(zac(:),xsurf_epsilon)
508 !
509 ! 7.4. Total surface fluxes
510 !
511 psfth(:)=zhf(:)+zrf(:)
512 psftq(:)=(zef(:)+zefwebb(:))/zlv(:)
513 !
514 ! 7.5. Z0 and Z0H over water
515 !
516 pz0sea(:) = zcharn(:) * zustar2(:) / xg + xvz0cm * pcd(:) / pcdn(:)
517 !
518 pz0hsea(:)=pz0sea(:)
519 !
520 IF (lhook) CALL dr_hook('ECUME_FLUX',1,zhook_handle)
521 !-------------------------------------------------------------------------------
522 
523 END SUBROUTINE ecume_flux
mode_thermos::qsat_seawater
Definition: mode_thermos.F90:62
wind_threshold
real function, dimension(size(pwind)) wind_threshold(PWIND, PUREF)
Definition: wind_threshold.F90:6
surface_ri
subroutine surface_ri(PTG, PQS, PEXNS, PEXNA, PTA, PQA, PZREF, PUREF, PDIRCOSZW, PVMOD, PRI)
Definition: surface_ri.F90:6
mode_thermos
Definition: mode_thermos.F90:6
modd_surf_atm
Definition: modd_surf_atm.F90:6
ecume_flux
subroutine ecume_flux(PZ0SEA, PTA, PEXNA, PRHOA, PSST, PEXNS, PQA, PVMOD, PZREF, PUREF, PPS, PICHCE, OPRECIP, OPWEBB, OPWG, PQSAT, PSFTH, PSFTQ, PUSTAR, PCD, PCDN, PCH, PCE, PRI, PRESA, PRAIN, PZ0HSEA, OPERTFLUX, PPERTFLUX)
Definition: ecume_flux.F90:6
modd_water_par
Definition: modd_water_par.F90:6
mode_thermos::qsat
Definition: mode_thermos.F90:57
modd_surf_par
Definition: modd_surf_par.F90:6
modd_snow_par
Definition: modd_snow_par.F90:6
modd_reprod_oper
Definition: modd_reprod_oper.F90:6
modd_csts
Definition: modd_csts.F90:6
Generated on Fri May 13 2016 09:31:03 for SURFEX by   doxygen 1.8.5