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

Contents of /trunk/dyn3d/bilan_dyn.f

Parent Directory Parent Directory | Revision Log Revision Log

Revision 91 - (show annotations)
Wed Mar 26 17:18:58 2014 UTC (10 years, 1 month ago) by guez
File size: 7219 byte(s) 
Removed unused variables lock_startdate and time_stamp of module
calendar.

Noticed that physiq does not change the surface pressure. So removed
arguments ps and dpfi of subroutine addfi. dpfi was always 0. The
computation of ps in addfi included some averaging at the poles. In
principle, this does not change ps but in practice it does because of
finite numerical precision. So the results of the simulation are
changed. Removed arguments ps and dpfi of calfis. Removed argument
d_ps of physiq.

du at the poles is not computed by dudv1, so declare only the
corresponding latitudes in dudv1. caldyn passes only a section of the
array dudyn as argument.

Removed variable niadv of module iniadvtrac_m.

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