source: LMDZ4/branches/LMDZ4_par_0/libf/phylmd/coefcdrag.F90 @ 2190

Last change on this file since 2190 was 524, checked in by lmdzadmin, 20 years ago

Initial revision

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 5.3 KB
Line 
1!
2! $Header$
3!
4!
5!
6!
7      SUBROUTINE coefcdrag (klon, knon, nsrf, zxli, &
8                            speed, t, q, zgeop, psol, &
9                            ts, qsurf, rugos, okri, ri1, &
10                            cdram, cdrah, cdran, zri1, pref)
11      IMPLICIT none
12!-------------------------------------------------------------------------
13! Objet : calcul des cdrags pour le moment (cdram) et les flux de chaleur
14!         sensible et latente (cdrah), du cdrag neutre (cdran),
15!         du nombre de Richardson entre la surface et le niveau de reference
16!         (zri1) et de la pression au niveau de reference (pref).   
17!
18! I. Musat, 01.07.2002
19!-------------------------------------------------------------------------
20!
21! klon----input-I- dimension de la grille physique (= nb_pts_latitude X nb_pts_longitude)
22! knon----input-I- nombre de points pour un type de surface
23! nsrf----input-I- indice pour le type de surface; voir indicesol.inc
24! zxli----input-L- TRUE si calcul des cdrags selon Laurent Li
25! speed---input-R- module du vent au 1er niveau du modele
26! t-------input-R- temperature de l'air au 1er niveau du modele
27! q-------input-R- humidite de l'air au 1er niveau du modele
28! zgeop---input-R- geopotentiel au 1er niveau du modele
29! psol----input-R- pression au sol
30! ts------input-R- temperature de l'air a la surface
31! qsurf---input-R- humidite de l'air a la surface
32! rugos---input-R- rugosite
33! okri----input-L- TRUE si on veut tester le nb. Richardson entre la sfce
34!                  et zref par rapport au Ri entre la sfce et la 1ere couche
35! ri1-----input-R- nb. Richardson entre la surface et la 1ere couche
36!
37! cdram--output-R- cdrag pour le moment
38! cdrah--output-R- cdrag pour les flux de chaleur latente et sensible
39! cdran--output-R- cdrag neutre
40! zri1---output-R- nb. Richardson entre la surface et la couche zgeop/RG
41! pref---output-R- pression au niveau zgeop/RG
42!
43      INTEGER, intent(in) :: klon, knon, nsrf
44      LOGICAL, intent(in) :: zxli
45      REAL, dimension(klon), intent(in) :: speed, t, q, zgeop, psol
46      REAL, dimension(klon), intent(in) :: ts, qsurf, rugos, ri1
47      LOGICAL, intent(in) :: okri   
48!
49      REAL, dimension(klon), intent(out) :: cdram, cdrah, cdran, zri1, pref
50!-------------------------------------------------------------------------
51!
52#include "YOMCST.inc"
53#include "YOETHF.inc"
54#include "indicesol.inc"
55! Quelques constantes :
56      REAL, parameter :: RKAR=0.40, CB=5.0, CC=5.0, CD=5.0
57!
58! Variables locales :
59      INTEGER :: i
60      REAL, dimension(klon) :: zdu2, zdphi, ztsolv, ztvd
61      REAL, dimension(klon) :: zscf, friv, frih, zucf, zcr
62      REAL, dimension(klon) :: zcfm1, zcfh1
63      REAL, dimension(klon) :: zcfm2, zcfh2
64      REAL, dimension(klon) :: trm0, trm1
65!-------------------------------------------------------------------------
66      REAL :: fsta, fins, x
67      fsta(x) = 1.0 / (1.0+10.0*x*(1+8.0*x))
68      fins(x) = SQRT(1.0-18.0*x)
69!-------------------------------------------------------------------------
70!
71      DO i = 1, knon
72!
73       zdphi(i) = zgeop(i)
74       zdu2(i) = speed(i)**2
75       pref(i) = exp(log(psol(i)) - zdphi(i)/(RD*t(i)* &
76                 (1.+ RETV * max(q(i),0.0))))
77       ztsolv(i) = ts(i)
78       ztvd(i) = t(i) * (psol(i)/pref(i))**RKAPPA
79       trm0(i) = 1. + RETV * max(qsurf(i),0.0)
80       trm1(i) = 1. + RETV * max(q(i),0.0)
81       ztsolv(i) = ztsolv(i) * trm0(i)
82       ztvd(i) = ztvd(i) * trm1(i)
83       zri1(i) = zdphi(i)*(ztvd(i)-ztsolv(i))/(zdu2(i)*ztvd(i))
84!
85! on teste zri1 par rapport au Richardson de la 1ere couche ri1
86!
87!IM +++
88       IF(1.EQ.0) THEN
89       IF (okri) THEN
90         IF (ri1(i).GE.0.0.AND.zri1(i).LT.0.0) THEN
91           zri1(i) = ri1(i)
92         ELSE IF(ri1(i).LT.0.0.AND.zri1(i).GE.0.0) THEN
93           zri1(i) = ri1(i)
94         ENDIF
95       ENDIF
96       ENDIF
97!IM ---
98!
99       cdran(i) = (RKAR/log(1.+zdphi(i)/(RG*rugos(i))))**2
100
101       IF (zri1(i) .ge. 0.) THEN
102!
103! situation stable : pour eviter les inconsistances dans les cas
104! tres stables on limite zri1 a 20. cf Hess et al. (1995)
105!
106         zri1(i) = min(20.,zri1(i))
107!
108         IF (.NOT.zxli) THEN
109           zscf(i) = SQRT(1.+CD*ABS(zri1(i)))
110           friv(i) = max(1. / (1.+2.*CB*zri1(i)/ zscf(i)), 0.1)
111           zcfm1(i) = cdran(i) * friv(i)
112           frih(i) = max(1./ (1.+3.*CB*zri1(i)*zscf(i)), 0.1 )
113           zcfh1(i) = cdran(i) * frih(i)
114           cdram(i) = zcfm1(i)
115           cdrah(i) = zcfh1(i)
116         ELSE
117           cdram(i) = cdran(i)* fsta(zri1(i))
118           cdrah(i) = cdran(i)* fsta(zri1(i))
119         ENDIF
120!
121       ELSE
122!
123! situation instable
124!
125         IF (.NOT.zxli) THEN
126           zucf(i) = 1./(1.+3.0*CB*CC*cdran(i)*SQRT(ABS(zri1(i)) &
127                 *(1.0+zdphi(i)/(RG*rugos(i)))))
128           zcfm2(i) = cdran(i)*max((1.-2.0*CB*zri1(i)*zucf(i)),0.1)
129           zcfh2(i) = cdran(i)*max((1.-3.0*CB*zri1(i)*zucf(i)),0.1)
130           cdram(i) = zcfm2(i)
131           cdrah(i) = zcfh2(i)
132         ELSE
133           cdram(i) = cdran(i)* fins(zri1(i))
134           cdrah(i) = cdran(i)* fins(zri1(i))
135         ENDIF
136!
137! cdrah sur l'ocean cf. Miller et al. (1992)
138!
139         zcr(i) = (0.0016/(cdran(i)*SQRT(zdu2(i))))*ABS(ztvd(i)-ztsolv(i)) &
140               **(1./3.)
141         IF (nsrf.EQ.is_oce) cdrah(i) = cdran(i)*(1.0+zcr(i)**1.25) &
142                  **(1./1.25)
143       ENDIF
144!
145      END DO
146      RETURN
147      END SUBROUTINE coefcdrag
Note: See TracBrowser for help on using the repository browser.