/[lmdze]/trunk/phylmd/ozonecm.f90
ViewVC logotype

Annotation of /trunk/phylmd/ozonecm.f90

Parent Directory Parent Directory | Revision Log Revision Log


Revision 328 - (hide annotations)
Thu Jun 13 14:40:06 2019 UTC (4 years, 11 months ago) by guez
File size: 3159 byte(s)
Change all `.f` suffixes to `.f90`. (The opposite was done in revision
82.)  Because of change of philosopy in GNUmakefile: we already had a
rewritten rule for `.f`, so it does not make the makefile longer to
replace it by a rule for `.f90`. And it spares us options of
makedepf90 and of the compiler. Also we prepare the way for a simpler
`CMakeLists.txt`.

1 guez 22 module ozonecm_m
2 guez 3
3 guez 22 IMPLICIT NONE
4 guez 3
5 guez 22 contains
6 guez 3
7 guez 22 function ozonecm(rjour, paprs)
8    
9     ! From phylmd/ozonecm.F, version 1.3 2005/06/06 13:16:33
10    
11     ! The ozone climatology is based on an analytic formula which fits the
12     ! Krueger and Mintzner (1976) profile, as well as the variations with
13     ! altitude and latitude of the maximum ozone concentrations and the total
14     ! column ozone concentration of Keating and Young (1986). The analytic
15     ! formula have been established by J.-F. Royer (CRNM, Meteo France), who
16     ! also provided us the code.
17    
18     ! A. J. Krueger and R. A. Minzner, A Mid-Latitude Ozone Model for the
19 guez 315 ! 1976 U. S. Standard Atmosphere, J. Geophys. Res., 81, 4477 (1976).
20 guez 22
21 guez 322 ! Keating 1985 k1032
22 guez 22
23 guez 315 use nr_util, only: assert, twopi, deg_to_rad
24 guez 278
25 guez 265 use dimensions, only: llm
26 guez 22 use phyetat0_m, only: rlat
27    
28     REAL, INTENT (IN) :: rjour
29 guez 118
30 guez 266 REAL, INTENT (IN) :: paprs(:, :) ! (klon, llm + 1)
31 guez 118 ! pression pour chaque inter-couche, en Pa
32 guez 22
33 guez 266 REAL ozonecm(size(paprs, 1), llm) ! (klon, llm)
34 guez 22 ! "ozonecm(j, k)" is the column-density of ozone in cell "(j, k)", that is
35     ! between interface "k" and interface "k + 1", in kDU.
36    
37 guez 98 ! Local:
38 guez 22
39 guez 118 REAL tozon ! equivalent pressure of ozone above interface "k", in Pa
40 guez 22 INTEGER i, k
41    
42 guez 266 REAL field(llm + 1)
43 guez 118 ! "field(k)" is the column-density of ozone between interface
44 guez 22 ! "k" and the top of the atmosphere (interface "llm + 1"), in kDU.
45    
46 guez 118 real, PARAMETER:: ps = 101325.
47     REAL, parameter:: an = 360., zo3q3 = 4E-8
48     real, parameter:: zo3a3 = zo3q3 / ps / 2.
49     REAL, parameter:: dobson_unit = 2.1415E-5 ! in kg m-2
50     REAL gms, slat, slat2, sint, cost, ppm, a
51     REAL asec, bsec, aprim
52 guez 22
53     !----------------------------------------------------------
54    
55 guez 266 call assert(size(paprs, 2) == llm + 1, "ozonecm")
56 guez 22
57 guez 315 sint = sin(twopi * (rjour + 15.) / an)
58     cost = cos(twopi * (rjour + 15.) / an)
59 guez 118 field(llm + 1) = 0.
60    
61 guez 266 DO i = 1, size(paprs, 1)
62 guez 315 slat = sin(deg_to_rad * rlat(i))
63 guez 118 slat2 = slat * slat
64     gms = 0.0531 + sint * (- 0.001595 + 0.009443 * slat) + cost &
65     * (- 0.001344 - 0.00346 * slat) + slat2 * (0.056222 + slat2 &
66     * (- 0.037609 + 0.012248 * sint + 0.00521 * cost + 0.00889 &
67     * slat)) - zo3q3 * ps
68     ppm = 800. - 500. * slat2 - 150. * cost * slat
69     bsec = 2650. + 5000. * slat2
70     a = 4. * bsec**1.5 * ppm**1.5 * (1. + (bsec / ps)**1.5) &
71     / (bsec**1.5 + ppm**1.5)**2
72 guez 316 asec = max(0., (gms - max(0., (2.666666 * (1.74E-5 - 7.5E-6 * slat2 &
73     - 1.7E-6 * cost * slat) * ppm - a * gms) / (1. - a))) &
74     * (1. + (bsec / ps)**1.5))
75 guez 118 aprim = gms - asec / (1. + (bsec / ps)**1.5)
76    
77     DO k = 1, llm
78     tozon = aprim / (1. + 3. * (ppm / paprs(i, k))**2) &
79     + asec / (1. + (bsec / paprs(i, k))**1.5) &
80     + zo3a3 * paprs(i, k) * paprs(i, k)
81 guez 22 ! Convert from Pa to kDU:
82 guez 118 field(k) = tozon / 9.81 / dobson_unit / 1e3
83 guez 22 END DO
84 guez 118
85     forall (k = 1: llm) ozonecm(i, k) = field(k) - field(k + 1)
86 guez 22 END DO
87    
88 guez 118 ozonecm = max(ozonecm, 1e-12)
89 guez 22
90     END function ozonecm
91    
92     end module ozonecm_m

  ViewVC Help
Powered by ViewVC 1.1.21