/[lmdze]/trunk/phylmd/Interface_surf/cdrag.f
ViewVC logotype

Contents of /trunk/phylmd/Interface_surf/cdrag.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 291 - (show annotations)
Wed Jul 25 14:15:44 2018 UTC (5 years, 10 months ago) by guez
File size: 4073 byte(s)
Use named constant f_ri_cd_min in procedure cdrag (following LMDZ).
Remove some intermediary variables.

1 module cdrag_m
2
3 IMPLICIT NONE
4
5 contains
6
7 SUBROUTINE cdrag(nsrf, speed, t, q, zgeop, psol, ts, qsurf, rugos, cdragm, &
8 cdragh, pref)
9
10 ! From LMDZ4/libf/phylmd/clcdrag.F90 and
11 ! LMDZ4/libf/phylmd/coefcdrag.F90, version 1.1.1.1, 2004/05/19
12 ! 12:53:07
13
14 ! Objet : calcul des drag coefficients au sol pour le moment et
15 ! les flux de chaleurs sensible et latente et calcul de la
16 ! pression au niveau de reference.
17
18 ! Ionela MUSAT, July, 1st, 2002
19
20 ! Louis, J. F., Tiedtke, M. and Geleyn, J. F., 1982: `A short
21 ! history of the operational PBL parametrization at
22 ! ECMWF'. Workshop on boundary layer parametrization, November
23 ! 1981, ECMWF, Reading, England. Page: 19. Equations in Table 1.
24
25 use nr_util, only: assert_eq
26
27 use clesphys, only: f_cdrag_oce, f_cdrag_ter
28 use indicesol, only: is_oce
29 use SUPHEC_M, only: rcpd, rd, retv, rg
30 USE yoethf_m, ONLY: rvtmp2
31
32 INTEGER, intent(in):: nsrf ! indice pour le type de surface
33
34 REAL, intent(in):: speed(:) ! (knon)
35 ! norm of the wind at the first model level
36
37 REAL, intent(in):: t(:) ! (knon)
38 ! temperature de l'air au 1er niveau du modele
39
40 REAL, intent(in):: q(:) ! (knon) ! humidite de l'air au 1er niveau du modele
41
42 REAL, intent(in):: zgeop(:) ! (knon)
43 ! g\'eopotentiel au 1er niveau du mod\`ele
44
45 REAL, intent(in) :: psol(:) ! (knon) pression au sol
46 REAL, intent(in):: ts(:) ! (knon) temperature de l'air a la surface
47 REAL, intent(in):: qsurf(:) ! (knon) humidite de l'air a la surface
48 REAL, intent(in):: rugos(:) ! (knon) rugosit\'e
49 REAL, intent(out):: cdragm(:) ! (knon) drag coefficient pour le moment
50
51 REAL, intent(out):: cdragh(:) ! (knon)
52 ! drag coefficient pour les flux de chaleur latente et sensible
53
54 REAL, intent(out), optional:: pref(:) ! (knon) pression au niveau zgeop / RG
55
56 ! Local:
57
58 REAL, PARAMETER:: ckap = 0.40, cb = 5., cc = 5., cd = 5., cepdu2 = 0.1**2
59 real, parameter:: f_ri_cd_min = 0.1
60 INTEGER i, knon
61 REAL zdu2, ztsolv, ztvd, zscf, zucf
62 real zcdn ! drag coefficient neutre
63
64 REAL zri
65 ! nombre de Richardson entre la surface et le niveau de reference
66 ! zgeop / RG
67
68 !-------------------------------------------------------------------------
69
70 knon = assert_eq([size(speed), size(t), size(q), size(zgeop), size(ts), &
71 size(qsurf), size(rugos), size(cdragm), size(cdragh)], &
72 "cdrag knon")
73
74 DO i = 1, knon
75 zdu2 = max(cepdu2, speed(i)**2)
76 ztsolv = ts(i) * (1. + RETV * max(qsurf(i), 0.))
77 ztvd = (t(i) + zgeop(i) / RCPD / (1. + RVTMP2 * q(i))) &
78 * (1. + RETV * q(i))
79 zri = zgeop(i) * (ztvd - ztsolv) / (zdu2 * ztvd)
80 zcdn = (ckap / log(1. + zgeop(i) / (RG * rugos(i))))**2
81
82 IF (zri > 0.) THEN
83 ! Situation stable. Pour \'eviter les incoh\'erences dans
84 ! les cas tr\`es stables, on limite zri \`a 20. Cf Hess et
85 ! al. (1995).
86 zri = min(20., zri)
87 zscf = SQRT(1. + cd * ABS(zri))
88 cdragm(i) = zcdn * max(1. / (1. + 2. * CB * zri / zscf), f_ri_cd_min)
89 cdragh(i) = merge(f_cdrag_oce, f_cdrag_ter, nsrf == is_oce) * zcdn &
90 * max(1. / (1. + 3. * CB * zri * zscf), f_ri_cd_min)
91 ELSE
92 ! situation instable
93 zucf = 1. / (1. + 3. * cb * cc * zcdn &
94 * SQRT(ABS(zri) * (1. + zgeop(i) / (RG * rugos(i)))))
95 cdragm(i) = zcdn * max((1. - 2. * cb * zri * zucf), f_ri_cd_min)
96
97 IF (nsrf == is_oce) then
98 ! Cf. Miller et al. (1992).
99 cdragh(i) = f_cdrag_oce * zcdn * (1. + ((0.0016 &
100 / (zcdn * SQRT(zdu2))) * ABS(ztvd - ztsolv)**(1. &
101 / 3.))**1.25)**(1. / 1.25)
102 else
103 cdragh(i) = f_cdrag_ter * zcdn &
104 * max((1. - 3. * cb * zri * zucf), f_ri_cd_min)
105 end IF
106 ENDIF
107 END DO
108
109 if (present(pref)) &
110 pref = exp(log(psol) - zgeop / (RD * t * (1. + RETV * max(q, 0.))))
111
112 END SUBROUTINE cdrag
113
114 end module cdrag_m

  ViewVC Help
Powered by ViewVC 1.1.21