source: dynamico_lmdz/simple_physics/phyparam/param/vdif_cd.F @ 4185

Last change on this file since 4185 was 4184, checked in by dubos, 6 years ago

simple_physics : planete.h => MODULE planet.F90

File size: 2.6 KB
Line 
1      SUBROUTINE vdif_cd( ngrid,pz0,pg,pz,pu,pv,pts,ph,pcdv,pcdh)
2      USE planet, ONLY : karman
3      IMPLICIT NONE
4c=======================================================================
5c
6c   Subject: computation of the surface drag coefficient using the
7c   -------  approch developed by Loui for ECMWF.
8c
9c   Author: Frederic Hourdin  15 /10 /93
10c   -------
11c
12c   Arguments:
13c   ----------
14c
15c   inputs:
16c   ------
17c     ngrid            size of the horizontal grid
18c     pg               gravity (m s -2)
19c     pz(ngrid)        height of the first atmospheric layer
20c     pu(ngrid)        u component of the wind in that layer
21c     pv(ngrid)        v component of the wind in that layer
22c     pts(ngrid)       surfacte temperature
23c     ph(ngrid)        potential temperature T*(p/ps)^kappa
24c
25c   outputs:
26c   --------
27c     pcdv(ngrid)      Cd for the wind
28c     pcdh(ngrid)      Cd for potential temperature
29c
30c=======================================================================
31c
32c-----------------------------------------------------------------------
33c   Declarations:
34c   -------------
35
36c   Arguments:
37c   ----------
38
39      INTEGER ngrid,nlay
40      REAL pz0(ngrid)
41      REAL pg,pz(ngrid)
42      REAL pu(ngrid),pv(ngrid)
43      REAL pts(ngrid),ph(ngrid)
44      REAL pcdv(ngrid),pcdh(ngrid)
45
46c   Local:
47c   ------
48
49      INTEGER ig
50
51      REAL zu2,z1,zri,zcd0,zz
52
53      REAL b,c,d,c2b,c3bc,c3b,z0,umin2
54      LOGICAL firstcal
55      DATA b,c,d,umin2/5.,5.,5.,1.e-12/
56      DATA firstcal/.true./
57      SAVE b,c,d,c2b,c3bc,c3b,firstcal,umin2
58
59c-----------------------------------------------------------------------
60c   couche de surface:
61c   ------------------
62
63c     DO ig=1,ngrid
64c        zu2=pu(ig)*pu(ig)+pv(ig)*pv(ig)+umin2
65c        pcdv(ig)=pz0(ig)*(1.+sqrt(zu2))
66c        pcdh(ig)=pcdv(ig)
67c     ENDDO
68c     RETURN
69
70      IF (firstcal) THEN
71         c2b=2.*b
72         c3bc=3.*b*c
73         c3b=3.*b
74         firstcal=.false.
75      ENDIF
76
77c!!!! WARNING, verifier la formule originale de Louis!
78      DO ig=1,ngrid
79         zu2=pu(ig)*pu(ig)+pv(ig)*pv(ig)+umin2
80         zri=pg*pz(ig)*(ph(ig)-pts(ig))/(ph(ig)*zu2)
81         z1=1.+pz(ig)/pz0(ig)
82         zcd0=karman/log(z1)
83         zcd0=zcd0*zcd0*sqrt(zu2)
84         IF(zri.LT.0.) THEN
85            z1=b*zri/(1.+c3bc*zcd0*sqrt(-z1*zri))
86            pcdv(ig)=zcd0*(1.-2.*z1)
87            pcdh(ig)=zcd0*(1.-3.*z1)
88         ELSE
89            zz=sqrt(1.+d*zri)
90            pcdv(ig)=zcd0/(1.+c2b*zri/zz)
91            pcdh(ig)=zcd0/(1.+c3b*zri*zz)
92         ENDIF
93      ENDDO
94
95c-----------------------------------------------------------------------
96
97      RETURN
98      END
Note: See TracBrowser for help on using the repository browser.