source: LMDZ4/trunk/libf/phylmd/clcdrag.F90 @ 1085

Last change on this file since 1085 was 1071, checked in by Laurent Fairhead, 16 years ago

Mauvaise integration des modifications faites sous cvs
LF

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 4.6 KB
RevLine 
[524]1!
[1067]2SUBROUTINE clcdrag(knon, nsrf, paprs, pplay,&
3     u1, v1, t1, q1, &
4     tsurf, qsurf, rugos, &
5     pcfm, pcfh)
6
7  USE dimphy
8  IMPLICIT NONE
[524]9! ================================================================= c
10!
11! Objet : calcul des cdrags pour le moment (pcfm) et
12!         les flux de chaleur sensible et latente (pcfh).   
13!
14! ================================================================= c
15!
16! knon----input-I- nombre de points pour un type de surface
[793]17! nsrf----input-I- indice pour le type de surface; voir indicesol.h
[1067]18! u1-------input-R- vent zonal au 1er niveau du modele
19! v1-------input-R- vent meridien au 1er niveau du modele
20! t1-------input-R- temperature de l'air au 1er niveau du modele
21! q1-------input-R- humidite de l'air au 1er niveau du modele
22! tsurf------input-R- temperature de l'air a la surface
[524]23! qsurf---input-R- humidite de l'air a la surface
24! rugos---input-R- rugosite
25!
26! pcfm---output-R- cdrag pour le moment
27! pcfh---output-R- cdrag pour les flux de chaleur latente et sensible
28!
[1071]29  INTEGER, INTENT(IN)                      :: knon, nsrf
30  REAL, DIMENSION(klon,klev+1), INTENT(IN) :: paprs
31  REAL, DIMENSION(klon,klev), INTENT(IN)   :: pplay
32  REAL, DIMENSION(klon), INTENT(IN)        :: u1, v1, t1, q1
33  REAL, DIMENSION(klon), INTENT(IN)        :: tsurf, qsurf
34  REAL, DIMENSION(klon), INTENT(IN)        :: rugos
35  REAL, DIMENSION(klon), INTENT(OUT)       :: pcfm, pcfh
36!
[524]37! ================================================================= c
38!
[1071]39  INCLUDE "YOMCST.h"
40  INCLUDE "YOETHF.h"
41  INCLUDE "indicesol.h"
42  INCLUDE "clesphys.h"
[524]43!
44! Quelques constantes et options:
45!!$PB      REAL, PARAMETER :: ckap=0.35, cb=5.0, cc=5.0, cd=5.0, cepdu2=(0.1)**2
[1067]46  REAL, PARAMETER :: ckap=0.40, cb=5.0, cc=5.0, cd=5.0, cepdu2=(0.1)**2
[524]47!
48! Variables locales :
[1067]49  INTEGER               :: i
50  REAL                  :: zdu2, ztsolv
51  REAL                  :: ztvd, zscf
52  REAL                  :: zucf, zcr
53  REAL                  :: friv, frih
54  REAL, DIMENSION(klon) :: zcfm1, zcfm2
55  REAL, DIMENSION(klon) :: zcfh1, zcfh2
56  REAL, DIMENSION(klon) :: zcdn
57  REAL, DIMENSION(klon) :: zri
58  REAL, DIMENSION(klon) :: zgeop1       ! geopotentiel au 1er niveau du modele
59  LOGICAL, PARAMETER    :: zxli=.FALSE. ! calcul des cdrags selon Laurent Li
[524]60!
61! Fonctions thermodynamiques et fonctions d'instabilite
[1067]62  REAL                  :: fsta, fins, x
63  fsta(x) = 1.0 / (1.0+10.0*x*(1+8.0*x))
64  fins(x) = SQRT(1.0-18.0*x)
65
[524]66! ================================================================= c
67!
[1067]68! Calculer le geopotentiel du premier couche de modele
69!
70  DO i = 1, knon
71     zgeop1(i) = RD * t1(i) / (0.5*(paprs(i,1)+pplay(i,1))) &
72          * (paprs(i,1)-pplay(i,1))
73  END DO
74! ================================================================= c
75!
[524]76! Calculer le frottement au sol (Cdrag)
77!
[1067]78  DO i = 1, knon
79     zdu2 = MAX(cepdu2,u1(i)**2+v1(i)**2)
80     ztsolv = tsurf(i) * (1.0+RETV*qsurf(i))
81     ztvd = (t1(i)+zgeop1(i)/RCPD/(1.+RVTMP2*q1(i))) &
82          *(1.+RETV*q1(i))
83     zri(i) = zgeop1(i)*(ztvd-ztsolv)/(zdu2*ztvd)
84     zcdn(i) = (ckap/LOG(1.+zgeop1(i)/(RG*rugos(i))))**2
85
[524]86!!$        IF (zri(i) .ge. 0.) THEN      ! situation stable
[1067]87     IF (zri(i) .GT. 0.) THEN      ! situation stable
88        zri(i) = MIN(20.,zri(i))
89        IF (.NOT.zxli) THEN
90           zscf = SQRT(1.+cd*ABS(zri(i)))
91           FRIV = AMAX1(1. / (1.+2.*CB*zri(i)/ZSCF), 0.1)
92           zcfm1(i) = zcdn(i) * FRIV
93           FRIH = AMAX1(1./ (1.+3.*CB*zri(i)*ZSCF), 0.1 )
[524]94!!$  PB          zcfh1(i) = zcdn(i) * FRIH
[1067]95           zcfh1(i) = 0.8 * zcdn(i) * FRIH
96           pcfm(i) = zcfm1(i)
97           pcfh(i) = zcfh1(i)
98        ELSE
99           pcfm(i) = zcdn(i)* fsta(zri(i))
100           pcfh(i) = zcdn(i)* fsta(zri(i))
101        ENDIF
102     ELSE                          ! situation instable
103        IF (.NOT.zxli) THEN
104           zucf = 1./(1.+3.0*cb*cc*zcdn(i)*SQRT(ABS(zri(i)) &
105                *(1.0+zgeop1(i)/(RG*rugos(i)))))
106           zcfm2(i) = zcdn(i)*amax1((1.-2.0*cb*zri(i)*zucf),0.1)
[524]107!!$PB            zcfh2(i) = zcdn(i)*amax1((1.-3.0*cb*zri(i)*zucf),0.1)
[1067]108           zcfh2(i) = 0.8 * zcdn(i)*amax1((1.-3.0*cb*zri(i)*zucf),0.1)
109           pcfm(i) = zcfm2(i)
110           pcfh(i) = zcfh2(i)
111        ELSE
112           pcfm(i) = zcdn(i)* fins(zri(i))
113           pcfh(i) = zcdn(i)* fins(zri(i))
[524]114        ENDIF
[1067]115        zcr = (0.0016/(zcdn(i)*SQRT(zdu2)))*ABS(ztvd-ztsolv)**(1./3.)
116        IF(nsrf.EQ.is_oce) pcfh(i) =0.8* zcdn(i)*(1.0+zcr**1.25)**(1./1.25)
117     ENDIF
118  END DO
119
120! ================================================================= c
121     
122  ! IM cf JLD : on seuille cdrag_m et cdrag_h
123  IF (nsrf == is_oce) THEN
124     DO i=1,knon
125        pcfm(i)=MIN(pcfm(i),cdmmax)
126        pcfh(i)=MIN(pcfh(i),cdhmax)
127     END DO
128  END IF
129
130END SUBROUTINE clcdrag
Note: See TracBrowser for help on using the repository browser.