/[lmdze]/trunk/Sources/phylmd/Radlwsw/swclr.f
ViewVC logotype

Annotation of /trunk/Sources/phylmd/Radlwsw/swclr.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 157 - (hide annotations)
Mon Jul 20 16:01:49 2015 UTC (8 years, 10 months ago) by guez
File size: 8497 byte(s)
Just encapsulated SUBROUTINE vlsplt in a module and cleaned it.

In procedure vlx, local variables dxqu and adxqu only need indices
iip2:ip1jm. Otherwise, just cleaned vlx.

Procedures dynredem0 and dynredem1 no longer have argument fichnom,
they just operate on a file named "restart.nc". The programming
guideline here is that gcm should not be more complex than it needs by
itself, other programs (ce0l etc.) just have to adapt to gcm. So ce0l
now creates files "restart.nc" and "restartphy.nc".

In order to facilitate decentralizing the writing of "restartphy.nc",
created a procedure phyredem0 out of phyredem. phyredem0 creates the
NetCDF header of "restartphy.nc" while phyredem writes the NetCDF
variables. As the global attribute itau_phy needs to be filled in
phyredem0, at the beginnig of the run, we must compute its value
instead of just using itap. So we have a dummy argument lmt_pas of
phyredem0. Also, the ncid of "startphy.nc" is upgraded from local
variable of phyetat0 to dummy argument. phyetat0 no longer closes
"startphy.nc".

Following the same decentralizing objective, the ncid of "restart.nc"
is upgraded from local variable of dynredem0 to module variable of
dynredem0_m. "restart.nc" is not closed at the end of dynredem0 nor
opened at the beginning of dynredem1.

In procedure etat0, instead of creating many vectors of size klon
which will be filled with zeroes, just create one array null_array.

In procedure phytrac, instead of writing trs(: 1) to a text file,
write it to "restartphy.nc" (following LMDZ). This is better because
now trs(: 1) is next to its coordinates. We can write to
"restartphy.nc" from phytrac directly, and not add trs(: 1) to the
long list of variables in physiq, thanks to the decentralizing of
"restartphy.nc".

In procedure phyetat0, we no longer write to standard output the
minimum and maximum values of read arrays. It is ok to check input and
abort on invalid values but just printing statistics on input seems too
much useless computation and out of place clutter.

1 guez 81 SUBROUTINE swclr(knu, paer, flag_aer, tauae, pizae, cgae, palbp, pdsig, &
2     prayl, psec, pcgaz, ppizaz, pray1, pray2, prefz, prj, prk, prmu0, ptauaz, &
3     ptra1, ptra2)
4     USE dimens_m
5     USE dimphy
6     USE raddim
7     USE radepsi
8     USE radopt
9     IMPLICIT NONE
10    
11     ! ------------------------------------------------------------------
12     ! PURPOSE.
13     ! --------
14     ! COMPUTES THE REFLECTIVITY AND TRANSMISSIVITY IN CASE OF
15     ! CLEAR-SKY COLUMN
16    
17     ! REFERENCE.
18     ! ----------
19    
20     ! SEE RADIATION'S PART OF THE ECMWF RESEARCH DEPARTMENT
21     ! DOCUMENTATION, AND FOUQUART AND BONNEL (1980)
22    
23     ! AUTHOR.
24     ! -------
25     ! JEAN-JACQUES MORCRETTE *ECMWF*
26    
27     ! MODIFICATIONS.
28     ! --------------
29     ! ORIGINAL : 94-11-15
30     ! ------------------------------------------------------------------
31     ! * ARGUMENTS:
32    
33     INTEGER knu
34     ! -OB
35     DOUBLE PRECISION flag_aer
36     DOUBLE PRECISION tauae(kdlon, kflev, 2)
37     DOUBLE PRECISION pizae(kdlon, kflev, 2)
38     DOUBLE PRECISION cgae(kdlon, kflev, 2)
39     DOUBLE PRECISION paer(kdlon, kflev, 5)
40     DOUBLE PRECISION palbp(kdlon, 2)
41     DOUBLE PRECISION pdsig(kdlon, kflev)
42     DOUBLE PRECISION prayl(kdlon)
43     DOUBLE PRECISION psec(kdlon)
44    
45     DOUBLE PRECISION pcgaz(kdlon, kflev)
46     DOUBLE PRECISION ppizaz(kdlon, kflev)
47     DOUBLE PRECISION pray1(kdlon, kflev+1)
48     DOUBLE PRECISION pray2(kdlon, kflev+1)
49     DOUBLE PRECISION prefz(kdlon, 2, kflev+1)
50     DOUBLE PRECISION prj(kdlon, 6, kflev+1)
51     DOUBLE PRECISION prk(kdlon, 6, kflev+1)
52     DOUBLE PRECISION prmu0(kdlon, kflev+1)
53     DOUBLE PRECISION ptauaz(kdlon, kflev)
54     DOUBLE PRECISION ptra1(kdlon, kflev+1)
55     DOUBLE PRECISION ptra2(kdlon, kflev+1)
56    
57     ! * LOCAL VARIABLES:
58    
59     DOUBLE PRECISION zc0i(kdlon, kflev+1)
60     DOUBLE PRECISION zcle0(kdlon, kflev)
61     DOUBLE PRECISION zclear(kdlon)
62     DOUBLE PRECISION zr21(kdlon)
63     DOUBLE PRECISION zr23(kdlon)
64     DOUBLE PRECISION zss0(kdlon)
65     DOUBLE PRECISION zscat(kdlon)
66     DOUBLE PRECISION ztr(kdlon, 2, kflev+1)
67    
68 guez 157 INTEGER jl, jk, ja, jkl, jklp1, jaj, jkm1
69 guez 81 DOUBLE PRECISION ztray, zgar, zratio, zff, zfacoa, zcorae
70     DOUBLE PRECISION zmue, zgap, zww, zto, zden, zmu1, zden1
71     DOUBLE PRECISION zbmu0, zbmu1, zre11
72    
73     ! ------------------------------------------------------------------
74    
75     ! * 1. OPTICAL PARAMETERS FOR AEROSOLS AND RAYLEIGH
76     ! --------------------------------------------
77    
78    
79     DO jk = 1, kflev + 1
80     DO ja = 1, 6
81     DO jl = 1, kdlon
82     prj(jl, ja, jk) = 0.
83     prk(jl, ja, jk) = 0.
84 guez 24 END DO
85 guez 81 END DO
86     END DO
87    
88     DO jk = 1, kflev
89     DO jl = 1, kdlon
90     ptauaz(jl, jk) = flag_aer*tauae(jl, jk, knu)
91     ppizaz(jl, jk) = flag_aer*pizae(jl, jk, knu)
92     pcgaz(jl, jk) = flag_aer*cgae(jl, jk, knu)
93     END DO
94    
95     IF (flag_aer>0) THEN
96     ! -OB
97     DO jl = 1, kdlon
98     ! PCGAZ(JL,JK)=PCGAZ(JL,JK)/PPIZAZ(JL,JK)
99     ! PPIZAZ(JL,JK)=PPIZAZ(JL,JK)/PTAUAZ(JL,JK)
100     ztray = prayl(jl)*pdsig(jl, jk)
101     zratio = ztray/(ztray+ptauaz(jl,jk))
102     zgar = pcgaz(jl, jk)
103     zff = zgar*zgar
104     ptauaz(jl, jk) = ztray + ptauaz(jl, jk)*(1.-ppizaz(jl,jk)*zff)
105     pcgaz(jl, jk) = zgar*(1.-zratio)/(1.+zgar)
106     ppizaz(jl, jk) = zratio + (1.-zratio)*ppizaz(jl, jk)*(1.-zff)/(1.- &
107     ppizaz(jl,jk)*zff)
108     END DO
109     ELSE
110     DO jl = 1, kdlon
111     ztray = prayl(jl)*pdsig(jl, jk)
112     ptauaz(jl, jk) = ztray
113     pcgaz(jl, jk) = 0.
114     ppizaz(jl, jk) = 1. - repsct
115     END DO
116     END IF ! check flag_aer
117     END DO
118    
119     ! ------------------------------------------------------------------
120    
121     ! * 2. TOTAL EFFECTIVE CLOUDINESS ABOVE A GIVEN LEVEL
122     ! ----------------------------------------------
123    
124    
125     DO jl = 1, kdlon
126     zr23(jl) = 0.
127     zc0i(jl, kflev+1) = 0.
128     zclear(jl) = 1.
129     zscat(jl) = 0.
130     END DO
131    
132     jk = 1
133     jkl = kflev + 1 - jk
134     jklp1 = jkl + 1
135     DO jl = 1, kdlon
136     zfacoa = 1. - ppizaz(jl, jkl)*pcgaz(jl, jkl)*pcgaz(jl, jkl)
137     zcorae = zfacoa*ptauaz(jl, jkl)*psec(jl)
138     zr21(jl) = exp(-zcorae)
139     zss0(jl) = 1. - zr21(jl)
140     zcle0(jl, jkl) = zss0(jl)
141    
142     IF (novlp==1) THEN
143     ! * maximum-random
144     zclear(jl) = zclear(jl)*(1.0-max(zss0(jl),zscat(jl)))/ &
145     (1.0-min(zscat(jl),1.-zepsec))
146     zc0i(jl, jkl) = 1.0 - zclear(jl)
147     zscat(jl) = zss0(jl)
148     ELSE IF (novlp==2) THEN
149     ! * maximum
150     zscat(jl) = max(zss0(jl), zscat(jl))
151     zc0i(jl, jkl) = zscat(jl)
152     ELSE IF (novlp==3) THEN
153     ! * random
154     zclear(jl) = zclear(jl)*(1.0-zss0(jl))
155     zscat(jl) = 1.0 - zclear(jl)
156     zc0i(jl, jkl) = zscat(jl)
157     END IF
158     END DO
159    
160     DO jk = 2, kflev
161     jkl = kflev + 1 - jk
162     jklp1 = jkl + 1
163     DO jl = 1, kdlon
164     zfacoa = 1. - ppizaz(jl, jkl)*pcgaz(jl, jkl)*pcgaz(jl, jkl)
165     zcorae = zfacoa*ptauaz(jl, jkl)*psec(jl)
166     zr21(jl) = exp(-zcorae)
167     zss0(jl) = 1. - zr21(jl)
168     zcle0(jl, jkl) = zss0(jl)
169    
170     IF (novlp==1) THEN
171     ! * maximum-random
172     zclear(jl) = zclear(jl)*(1.0-max(zss0(jl),zscat(jl)))/ &
173     (1.0-min(zscat(jl),1.-zepsec))
174     zc0i(jl, jkl) = 1.0 - zclear(jl)
175     zscat(jl) = zss0(jl)
176     ELSE IF (novlp==2) THEN
177     ! * maximum
178     zscat(jl) = max(zss0(jl), zscat(jl))
179     zc0i(jl, jkl) = zscat(jl)
180     ELSE IF (novlp==3) THEN
181     ! * random
182     zclear(jl) = zclear(jl)*(1.0-zss0(jl))
183     zscat(jl) = 1.0 - zclear(jl)
184     zc0i(jl, jkl) = zscat(jl)
185 guez 24 END IF
186 guez 81 END DO
187     END DO
188    
189     ! ------------------------------------------------------------------
190    
191     ! * 3. REFLECTIVITY/TRANSMISSIVITY FOR PURE SCATTERING
192     ! -----------------------------------------------
193    
194    
195     DO jl = 1, kdlon
196     pray1(jl, kflev+1) = 0.
197     pray2(jl, kflev+1) = 0.
198     prefz(jl, 2, 1) = palbp(jl, knu)
199     prefz(jl, 1, 1) = palbp(jl, knu)
200     ptra1(jl, kflev+1) = 1.
201     ptra2(jl, kflev+1) = 1.
202     END DO
203    
204     DO jk = 2, kflev + 1
205     jkm1 = jk - 1
206     DO jl = 1, kdlon
207    
208    
209     ! ------------------------------------------------------------------
210    
211     ! * 3.1 EQUIVALENT ZENITH ANGLE
212     ! -----------------------
213    
214    
215     zmue = (1.-zc0i(jl,jk))*psec(jl) + zc0i(jl, jk)*1.66
216     prmu0(jl, jk) = 1./zmue
217    
218    
219     ! ------------------------------------------------------------------
220    
221     ! * 3.2 REFLECT./TRANSMISSIVITY DUE TO RAYLEIGH AND AEROSOLS
222     ! ----------------------------------------------------
223    
224    
225     zgap = pcgaz(jl, jkm1)
226     zbmu0 = 0.5 - 0.75*zgap/zmue
227     zww = ppizaz(jl, jkm1)
228     zto = ptauaz(jl, jkm1)
229     zden = 1. + (1.-zww+zbmu0*zww)*zto*zmue + (1-zww)*(1.-zww+2.*zbmu0*zww) &
230     *zto*zto*zmue*zmue
231     pray1(jl, jkm1) = zbmu0*zww*zto*zmue/zden
232     ptra1(jl, jkm1) = 1./zden
233    
234     zmu1 = 0.5
235     zbmu1 = 0.5 - 0.75*zgap*zmu1
236     zden1 = 1. + (1.-zww+zbmu1*zww)*zto/zmu1 + (1-zww)*(1.-zww+2.*zbmu1*zww &
237     )*zto*zto/zmu1/zmu1
238     pray2(jl, jkm1) = zbmu1*zww*zto/zmu1/zden1
239     ptra2(jl, jkm1) = 1./zden1
240    
241    
242    
243     prefz(jl, 1, jk) = (pray1(jl,jkm1)+prefz(jl,1,jkm1)*ptra1(jl,jkm1)* &
244     ptra2(jl,jkm1)/(1.-pray2(jl,jkm1)*prefz(jl,1,jkm1)))
245    
246     ztr(jl, 1, jkm1) = (ptra1(jl,jkm1)/(1.-pray2(jl,jkm1)*prefz(jl,1, &
247     jkm1)))
248    
249     prefz(jl, 2, jk) = (pray1(jl,jkm1)+prefz(jl,2,jkm1)*ptra1(jl,jkm1)* &
250     ptra2(jl,jkm1))
251    
252     ztr(jl, 2, jkm1) = ptra1(jl, jkm1)
253    
254     END DO
255     END DO
256     DO jl = 1, kdlon
257     zmue = (1.-zc0i(jl,1))*psec(jl) + zc0i(jl, 1)*1.66
258     prmu0(jl, 1) = 1./zmue
259     END DO
260    
261    
262     ! ------------------------------------------------------------------
263    
264     ! * 3.5 REFLECT./TRANSMISSIVITY BETWEEN SURFACE AND LEVEL
265     ! -------------------------------------------------
266    
267    
268     IF (knu==1) THEN
269     jaj = 2
270     DO jl = 1, kdlon
271     prj(jl, jaj, kflev+1) = 1.
272     prk(jl, jaj, kflev+1) = prefz(jl, 1, kflev+1)
273     END DO
274    
275     DO jk = 1, kflev
276     jkl = kflev + 1 - jk
277     jklp1 = jkl + 1
278     DO jl = 1, kdlon
279     zre11 = prj(jl, jaj, jklp1)*ztr(jl, 1, jkl)
280     prj(jl, jaj, jkl) = zre11
281     prk(jl, jaj, jkl) = zre11*prefz(jl, 1, jkl)
282     END DO
283     END DO
284    
285     ELSE
286    
287     DO jaj = 1, 2
288     DO jl = 1, kdlon
289     prj(jl, jaj, kflev+1) = 1.
290     prk(jl, jaj, kflev+1) = prefz(jl, jaj, kflev+1)
291     END DO
292    
293     DO jk = 1, kflev
294     jkl = kflev + 1 - jk
295     jklp1 = jkl + 1
296     DO jl = 1, kdlon
297     zre11 = prj(jl, jaj, jklp1)*ztr(jl, jaj, jkl)
298     prj(jl, jaj, jkl) = zre11
299     prk(jl, jaj, jkl) = zre11*prefz(jl, jaj, jkl)
300     END DO
301     END DO
302     END DO
303    
304     END IF
305    
306     ! ------------------------------------------------------------------
307    
308     RETURN
309     END SUBROUTINE swclr

  ViewVC Help
Powered by ViewVC 1.1.21