/[lmdze]/trunk/Sources/phylmd/Orography/gwprofil.f
ViewVC logotype

Contents of /trunk/Sources/phylmd/Orography/gwprofil.f

Parent Directory Parent Directory | Revision Log Revision Log

Revision 178 - (show annotations)
Fri Mar 11 18:47:26 2016 UTC (8 years, 2 months ago) by guez
File size: 5215 byte(s) 
Moved variables date0, deltat, datasz_max, ncvar_ids, point, buff_pos,
buffer, regular from module histcom_var to modules where they are
defined.

Removed procedure ioipslmpp, useless for a sequential program.

Added argument datasz_max to histwrite_real (to avoid circular
dependency with histwrite).

Removed useless variables and computations everywhere.

Changed real litteral constants from default kind to double precision
in lwb, lwu, lwvn, sw1s, swtt, swtt1, swu.

Removed unused arguments: paer of sw, sw1s, sw2s, swclr; pcldsw of
sw1s, sw2s; pdsig, prayl of swr; co2_ppm of clmain, clqh; tsol of
transp_lay; nsrf of screenp; kcrit and kknu of gwstress; pstd of
orosetup.

Added output of relative humidity.
1 module gwprofil_m
2
3 IMPLICIT NONE
4
5 contains
6
7 SUBROUTINE gwprofil(nlon, nlev, ktest, kkcrith, kcrit, paphm1, prho, pstab, &
8 pvph, pri, ptau, pdmod, psig, pvar)
9
10 ! Method. The stress profile for gravity waves is computed as
11 ! follows: it is constant (no gwd) at the levels between the
12 ! ground and the top of the blocked layer (kkenvh). It decreases
13 ! linearly with height from the top of the blocked layer to 3 *
14 ! varor (kknu), to simulate lee waves or nonlinear gravity wave
15 ! breaking. Above it is constant, except when the wave encounters
16 ! a critical level (kcrit) or when it breaks.
17
18 ! Reference. See ECMWF research department documentation of the
19 ! "I.F.S."
20
21 ! Modifications. Passage of the new gwdrag TO I.F.S. (F. LOTT,
22 ! 22/11/93)
23
24 USE dimphy, ONLY : klev, klon
25 USE yoegwd, ONLY : gkdrag, grahilo, grcrit, gssec, gtsec, nstra
26
27 INTEGER, intent(in):: nlon, nlev
28 INTEGER, intent(in):: ktest(nlon), kkcrith(nlon), kcrit(nlon)
29 REAL, intent(in):: paphm1(nlon, nlev+1), prho(nlon, nlev+1)
30 REAL, intent(in):: pstab(nlon, nlev+1)
31 real, intent(in):: pvph(nlon, nlev+1), pri(nlon, nlev+1)
32 real ptau(nlon, nlev+1)
33 REAL, intent(in):: pdmod(nlon)
34 REAL, INTENT (IN) :: psig(nlon)
35 REAL, INTENT (IN) :: pvar(nlon)
36
37 ! Local:
38 INTEGER jl, jk
39 REAL zsqr, zalfa, zriw, zdel, zb, zalpha, zdz2n
40 REAL zdelp, zdelpt
41 REAL zdz2(klon, klev), znorm(klon), zoro(klon)
42 REAL ztau(klon, klev+1)
43
44 !-----------------------------------------------------------------------
45
46 ! 1. INITIALIZATION
47
48 ! COMPUTATIONAL CONSTANTS.
49
50 DO jl = 1, klon
51 IF (ktest(jl)==1) THEN
52 zoro(jl) = psig(jl)*pdmod(jl)/4./max(pvar(jl), 1.0)
53 ztau(jl, klev+1) = ptau(jl, klev+1)
54 END IF
55 end DO
56
57 DO jk = klev, 2, -1
58 ! 4.1 CONSTANT WAVE STRESS UNTIL TOP OF THE
59 ! BLOCKING LAYER.
60 DO jl = 1, klon
61 IF (ktest(jl)==1) THEN
62 IF (jk>kkcrith(jl)) THEN
63 ptau(jl, jk) = ztau(jl, klev+1)
64 ELSE
65 ptau(jl, jk) = grahilo*ztau(jl, klev+1)
66 END IF
67 END IF
68 end DO
69
70 ! 4.2 WAVE DISPLACEMENT AT NEXT LEVEL.
71 DO jl = 1, klon
72 IF (ktest(jl)==1) THEN
73 IF (jk<kkcrith(jl)) THEN
74 znorm(jl) = gkdrag * prho(jl, jk) * sqrt(pstab(jl, jk)) &
75 * pvph(jl, jk)* zoro(jl)
76 zdz2(jl, jk) = ptau(jl, jk+1)/max(znorm(jl), gssec)
77 END IF
78 END IF
79 end DO
80
81 ! 4.3 WAVE RICHARDSON NUMBER, NEW WAVE DISPLACEMENT
82 ! AND STRESS: BREAKING EVALUATION AND CRITICAL
83 ! LEVEL
84 DO jl = 1, klon
85 IF (ktest(jl)==1) THEN
86 IF (jk<kkcrith(jl)) THEN
87 IF ((ptau(jl, jk+1)<gtsec) .OR. (jk<=kcrit(jl))) THEN
88 ptau(jl, jk) = 0.0
89 ELSE
90 zsqr = sqrt(pri(jl, jk))
91 zalfa = sqrt(pstab(jl, jk)*zdz2(jl, jk))/pvph(jl, jk)
92 zriw = pri(jl, jk)*(1.-zalfa)/(1+zalfa*zsqr)**2
93 IF (zriw<grcrit) THEN
94 zdel = 4./zsqr/grcrit + 1./grcrit**2 + 4./grcrit
95 zb = 1./grcrit + 2./zsqr
96 zalpha = 0.5*(-zb+sqrt(zdel))
97 zdz2n = (pvph(jl, jk)*zalpha)**2/pstab(jl, jk)
98 ptau(jl, jk) = znorm(jl)*zdz2n
99 ELSE
100 ptau(jl, jk) = znorm(jl)*zdz2(jl, jk)
101 END IF
102 ptau(jl, jk) = min(ptau(jl, jk), ptau(jl, jk+1))
103 END IF
104 END IF
105 END IF
106 end DO
107 end DO
108
109 ! REORGANISATION OF THE STRESS PROFILE AT LOW LEVEL
110
111 DO jl = 1, klon
112 IF (ktest(jl)==1) THEN
113 ztau(jl, kkcrith(jl)) = ptau(jl, kkcrith(jl))
114 ztau(jl, nstra) = ptau(jl, nstra)
115 END IF
116 end DO
117
118 DO jk = 1, klev
119 DO jl = 1, klon
120 IF (ktest(jl)==1) THEN
121 IF (jk>kkcrith(jl)) THEN
122 zdelp = paphm1(jl, jk) - paphm1(jl, klev+1)
123 zdelpt = paphm1(jl, kkcrith(jl)) - paphm1(jl, klev+1)
124 ptau(jl, jk) = ztau(jl, klev+1) &
125 + (ztau(jl, kkcrith(jl)) - ztau(jl, klev+1))*zdelp/zdelpt
126 END IF
127 END IF
128 end DO
129
130 ! REORGANISATION IN THE STRATOSPHERE
131 DO jl = 1, klon
132 IF (ktest(jl)==1) THEN
133 IF (jk < nstra) THEN
134 zdelp = paphm1(jl, nstra)
135 zdelpt = paphm1(jl, jk)
136 ptau(jl, jk) = ztau(jl, nstra) * zdelpt / zdelp
137 END IF
138 END IF
139 end DO
140
141 ! REORGANISATION IN THE TROPOSPHERE
142 DO jl = 1, klon
143 IF (ktest(jl)==1) THEN
144 IF (jk<kkcrith(jl) .AND. jk > nstra) THEN
145 zdelp = paphm1(jl, jk) - paphm1(jl, kkcrith(jl))
146 zdelpt = paphm1(jl, nstra) - paphm1(jl, kkcrith(jl))
147 ptau(jl, jk) = ztau(jl, kkcrith(jl)) &
148 + (ztau(jl, nstra) - ztau(jl, kkcrith(jl)))*zdelp/zdelpt
149 END IF
150 END IF
151 end DO
152 end DO
153
154 END SUBROUTINE gwprofil
155
156 end module gwprofil_m
  ViewVC Help
Powered by ViewVC 1.1.21