/[lmdze]/trunk/dyn3d/bilan_dyn.f
ViewVC logotype

Contents of /trunk/dyn3d/bilan_dyn.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 161 - (show annotations)
Fri Jul 24 14:27:59 2015 UTC (8 years, 10 months ago) by guez
Original Path: trunk/Sources/dyn3d/bilan_dyn.f
File size: 7301 byte(s)
rlon[uv] and rlat[uv] are already in start.nc.

Just encapsulated covcont in a module.

finvmaold was not used in leapfrog. Downgraded it from dummy argument
to local variable of SUBROUTINE integrd.

Simplified handling of mass in integrd: down from five 3-dimensional
arrays (masse, massem1, finvmaold, massescr and finvmasse) to three
(masse, massem1, finvmaold).

1 module bilan_dyn_m
2
3 IMPLICIT NONE
4
5 contains
6
7 SUBROUTINE bilan_dyn(ps, masse, pk, flux_u, flux_v, teta, phi, ucov, vcov, &
8 trac)
9
10 ! From LMDZ4/libf/dyn3d/bilan_dyn.F, version 1.5 2005/03/16 10:12:17
11
12 ! Sous-programme consacr\'e \`a des diagnostics dynamiques de
13 ! base. De fa\c{}con g\'en\'erale, les moyennes des scalaires Q
14 ! sont pond\'er\'ees par la masse. Les flux de masse sont, eux,
15 ! simplement moyenn\'es.
16
17 USE comconst, ONLY: cpp
18 USE comgeom, ONLY: constang_2d, cu_2d, cv_2d
19 use covcont_m, only: covcont
20 USE dimens_m, ONLY: iim, jjm, llm
21 USE histwrite_m, ONLY: histwrite
22 use init_dynzon_m, only: ncum, fileid, znom, ntr, nq, nom
23 use massbar_m, only: massbar
24 USE paramet_m, ONLY: iip1, jjp1
25
26 real, intent(in):: ps(iip1, jjp1)
27 real, intent(in):: masse(iip1, jjp1, llm), pk(iip1, jjp1, llm)
28 real, intent(in):: flux_u(iip1, jjp1, llm)
29 real, intent(in):: flux_v(iip1, jjm, llm)
30 real, intent(in):: teta(iip1, jjp1, llm)
31 real, intent(in):: phi(iip1, jjp1, llm)
32 real, intent(in):: ucov(:, :, :) ! (iip1, jjp1, llm)
33 real, intent(in):: vcov(iip1, jjm, llm)
34 real, intent(in):: trac(:, :, :) ! (iim + 1, jjm + 1, llm)
35
36 ! Local:
37
38 integer:: icum = 0
39 integer:: itau = 0
40 real qy, factv(jjm, llm)
41
42 ! Variables dynamiques interm\'ediaires
43 REAL vcont(iip1, jjm, llm), ucont(iip1, jjp1, llm)
44 REAL ang(iip1, jjp1, llm), unat(iip1, jjp1, llm)
45 REAL massebx(iip1, jjp1, llm), masseby(iip1, jjm, llm)
46 REAL ecin(iip1, jjp1, llm)
47
48 ! Champ contenant les scalaires advect\'es
49 real Q(iip1, jjp1, llm, nQ)
50
51 ! Champs cumul\'es
52 real, save:: ps_cum(iip1, jjp1)
53 real, save:: masse_cum(iip1, jjp1, llm)
54 real, save:: flux_u_cum(iip1, jjp1, llm)
55 real, save:: flux_v_cum(iip1, jjm, llm)
56 real, save:: Q_cum(iip1, jjp1, llm, nQ)
57 real, save:: flux_uQ_cum(iip1, jjp1, llm, nQ)
58 real, save:: flux_vQ_cum(iip1, jjm, llm, nQ)
59
60 ! champs de tansport en moyenne zonale
61 integer itr
62 integer, parameter:: iave = 1, itot = 2, immc = 3, itrs = 4, istn = 5
63
64 real vq(jjm, llm, ntr, nQ), vqtmp(jjm, llm)
65 real avq(jjm, 2: ntr, nQ), psiQ(jjm, llm + 1, nQ)
66 real zmasse(jjm, llm)
67 real v(jjm, llm), psi(jjm, llm + 1)
68 integer i, j, l, iQ
69
70 !-----------------------------------------------------------------
71
72 ! Calcul des champs dynamiques
73
74 ! \'Energie cin\'etique
75 ucont = 0
76 CALL covcont(llm, ucov, vcov, ucont, vcont)
77 CALL enercin(vcov, ucov, vcont, ucont, ecin)
78
79 ! moment cin\'etique
80 forall (l = 1: llm)
81 ang(:, :, l) = ucov(:, :, l) + constang_2d
82 unat(:, :, l) = ucont(:, :, l) * cu_2d
83 end forall
84
85 Q(:, :, :, 1) = teta * pk / cpp
86 Q(:, :, :, 2) = phi
87 Q(:, :, :, 3) = ecin
88 Q(:, :, :, 4) = ang
89 Q(:, :, :, 5) = unat
90 Q(:, :, :, 6) = trac
91 Q(:, :, :, 7) = 1.
92
93 ! Cumul
94
95 if (icum == 0) then
96 ps_cum = 0.
97 masse_cum = 0.
98 flux_u_cum = 0.
99 flux_v_cum = 0.
100 Q_cum = 0.
101 flux_vQ_cum = 0.
102 flux_uQ_cum = 0.
103 endif
104
105 itau = itau + 1
106 icum = icum + 1
107
108 ! Accumulation des flux de masse horizontaux
109 ps_cum = ps_cum + ps
110 masse_cum = masse_cum + masse
111 flux_u_cum = flux_u_cum + flux_u
112 flux_v_cum = flux_v_cum + flux_v
113 forall (iQ = 1: nQ) Q_cum(:, :, :, iQ) = Q_cum(:, :, :, iQ) &
114 + Q(:, :, :, iQ) * masse
115
116 ! Flux longitudinal
117 forall (iQ = 1: nQ, i = 1: iim) flux_uQ_cum(i, :, :, iQ) &
118 = flux_uQ_cum(i, :, :, iQ) &
119 + flux_u(i, :, :) * 0.5 * (Q(i, :, :, iQ) + Q(i + 1, :, :, iQ))
120 flux_uQ_cum(iip1, :, :, :) = flux_uQ_cum(1, :, :, :)
121
122 ! Flux m\'eridien
123 forall (iQ = 1: nQ, j = 1: jjm) flux_vQ_cum(:, j, :, iQ) &
124 = flux_vQ_cum(:, j, :, iQ) &
125 + flux_v(:, j, :) * 0.5 * (Q(:, j, :, iQ) + Q(:, j + 1, :, iQ))
126
127 writing_step: if (icum == ncum) then
128 ! Normalisation
129 forall (iQ = 1: nQ) Q_cum(:, :, :, iQ) = Q_cum(:, :, :, iQ) / masse_cum
130 ps_cum = ps_cum / ncum
131 masse_cum = masse_cum / ncum
132 flux_u_cum = flux_u_cum / ncum
133 flux_v_cum = flux_v_cum / ncum
134 flux_uQ_cum = flux_uQ_cum / ncum
135 flux_vQ_cum = flux_vQ_cum / ncum
136
137 ! Transport m\'eridien
138
139 ! Cumul zonal des masses des mailles
140
141 v = 0.
142 zmasse = 0.
143 call massbar(masse_cum, massebx, masseby)
144 do l = 1, llm
145 do j = 1, jjm
146 do i = 1, iim
147 zmasse(j, l) = zmasse(j, l) + masseby(i, j, l)
148 v(j, l) = v(j, l) + flux_v_cum(i, j, l)
149 enddo
150 factv(j, l) = cv_2d(1, j) / zmasse(j, l)
151 enddo
152 enddo
153
154 ! Transport dans le plan latitude-altitude
155
156 vq = 0.
157 psiQ = 0.
158 do iQ = 1, nQ
159 vqtmp = 0.
160 do l = 1, llm
161 do j = 1, jjm
162 ! Calcul des moyennes zonales du transport total et de vqtmp
163 do i = 1, iim
164 vq(j, l, itot, iQ) = vq(j, l, itot, iQ) &
165 + flux_vQ_cum(i, j, l, iQ)
166 qy = 0.5 * (Q_cum(i, j, l, iQ) * masse_cum(i, j, l) &
167 + Q_cum(i, j + 1, l, iQ) * masse_cum(i, j + 1, l))
168 vqtmp(j, l) = vqtmp(j, l) + flux_v_cum(i, j, l) * qy &
169 / (0.5 * (masse_cum(i, j, l) + masse_cum(i, j + 1, l)))
170 vq(j, l, iave, iQ) = vq(j, l, iave, iQ) + qy
171 enddo
172 ! Decomposition
173 vq(j, l, iave, iQ) = vq(j, l, iave, iQ) / zmasse(j, l)
174 vq(j, l, itot, iQ) = vq(j, l, itot, iQ) * factv(j, l)
175 vqtmp(j, l) = vqtmp(j, l) * factv(j, l)
176 vq(j, l, immc, iQ) = v(j, l) * vq(j, l, iave, iQ) * factv(j, l)
177 vq(j, l, itrs, iQ) = vq(j, l, itot, iQ) - vqtmp(j, l)
178 vq(j, l, istn, iQ) = vqtmp(j, l) - vq(j, l, immc, iQ)
179 enddo
180 enddo
181 ! Fonction de courant m\'eridienne pour la quantit\'e Q
182 do l = llm, 1, -1
183 do j = 1, jjm
184 psiQ(j, l, iQ) = psiQ(j, l + 1, iQ) + vq(j, l, itot, iQ)
185 enddo
186 enddo
187 enddo
188
189 ! Fonction de courant pour la circulation m\'eridienne moyenne
190 psi = 0.
191 do l = llm, 1, -1
192 do j = 1, jjm
193 psi(j, l) = psi(j, l + 1) + v(j, l)
194 v(j, l) = v(j, l) * factv(j, l)
195 enddo
196 enddo
197
198 ! Sorties proprement dites
199 do iQ = 1, nQ
200 do itr = 1, ntr
201 call histwrite(fileid, znom(itr, iQ), itau, vq(:, :, itr, iQ))
202 enddo
203 call histwrite(fileid, 'psi' // nom(iQ), itau, psiQ(:, :llm, iQ))
204 enddo
205
206 call histwrite(fileid, 'masse', itau, zmasse)
207 call histwrite(fileid, 'v', itau, v)
208 psi = psi * 1e-9
209 call histwrite(fileid, 'psi', itau, psi(:, :llm))
210
211 ! Int\'egrale verticale
212
213 forall (iQ = 1: nQ, itr = 2: ntr) avq(:, itr, iQ) &
214 = sum(vq(:, :, itr, iQ) * zmasse, dim=2) / cv_2d(1, :)
215
216 do iQ = 1, nQ
217 do itr = 2, ntr
218 call histwrite(fileid, 'a' // znom(itr, iQ), itau, avq(:, itr, iQ))
219 enddo
220 enddo
221
222 icum = 0
223 endif writing_step
224
225 end SUBROUTINE bilan_dyn
226
227 end module bilan_dyn_m

  ViewVC Help
Powered by ViewVC 1.1.21