source: LMDZ4/branches/LMDZ4-dev-20091210/libf/phylmd/clcdrag.F90 @ 4759

Last change on this file since 4759 was 1204, checked in by Laurent Fairhead, 15 years ago

Nettoyage du contrôle des coefficients de frottement en surface séparément
sur terre et sur océan (disparition de f_cdrag_stable)
PasB

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