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

Diff of /trunk/dyn3d/bilan_dyn.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

trunk/libf/dyn3d/bilan_dyn.f90 revision 55 by guez, Mon Dec 12 13:25:01 2011 UTC trunk/dyn3d/bilan_dyn.f revision 254 by guez, Mon Feb 5 10:39:38 2018 UTC
# Line 5  module bilan_dyn_m Line 5  module bilan_dyn_m
5  contains  contains
6    
7    SUBROUTINE bilan_dyn(ps, masse, pk, flux_u, flux_v, teta, phi, ucov, vcov, &    SUBROUTINE bilan_dyn(ps, masse, pk, flux_u, flux_v, teta, phi, ucov, vcov, &
8         trac, dt_app, dt_cum)         trac)
9    
10      ! From LMDZ4/libf/dyn3d/bilan_dyn.F, version 1.5 2005/03/16      ! From LMDZ4/libf/dyn3d/bilan_dyn.F, version 1.5 2005/03/16 10:12:17
     ! 10:12:17 fairhead  
11    
12      ! Sous-programme consacré à des diagnostics dynamiques de base.      ! Sous-programme consacr\'e \`a des diagnostics dynamiques de
13      ! De façon générale, les moyennes des scalaires Q sont pondérées      ! base.  De fa\c{}con g\'en\'erale, les moyennes des scalaires Q
14      ! par la masse. Les flux de masse sont, eux, simplement moyennés.      ! sont pond\'er\'ees par la masse. Les flux de masse sont, eux,
15        ! simplement moyenn\'es.
16    
17      USE histcom, ONLY: histbeg_totreg, histdef, histend, histvert      USE comconst, ONLY: cpp
18      USE calendar, ONLY: ymds2ju      USE comgeom, ONLY: constang_2d, cu_2d, cv_2d
19      USE histwrite_m, ONLY: histwrite      use covcont_m, only: covcont
20      USE dimens_m, ONLY: iim, jjm, llm      USE dimens_m, ONLY: iim, jjm, llm
21        use enercin_m, only: enercin
22        USE histwrite_m, ONLY: histwrite
23        use init_dynzon_m, only: ncum, fileid, znom, ntr, nq, nom
24        use massbar_m, only: massbar
25      USE paramet_m, ONLY: iip1, jjp1      USE paramet_m, ONLY: iip1, jjp1
26      USE comconst, ONLY: cpp  
27      USE comvert, ONLY: presnivs      real, intent(in):: ps(iip1, jjp1)
28      USE comgeom, ONLY: constang_2d, cu_2d, cv_2d, rlatv      real, intent(in):: masse(iip1, jjp1, llm), pk(iip1, jjp1, llm)
29      USE temps, ONLY: annee_ref, day_ref, itau_dyn      real, intent(in):: flux_u(iip1, jjp1, llm)
30      USE inigrads_m, ONLY: inigrads      real, intent(in):: flux_v(iip1, jjm, llm)
     USE nr_util, ONLY: pi  
   
     ! Arguments:  
   
     real, intent(in):: dt_app, dt_cum  
     real ps(iip1, jjp1)  
     real masse(iip1, jjp1, llm), pk(iip1, jjp1, llm)  
     real flux_u(iip1, jjp1, llm)  
     real flux_v(iip1, jjm, llm)  
31      real, intent(in):: teta(iip1, jjp1, llm)      real, intent(in):: teta(iip1, jjp1, llm)
32      real phi(iip1, jjp1, llm)      real, intent(in):: phi(iip1, jjp1, llm)
33      real ucov(iip1, jjp1, llm)      real, intent(in):: ucov(:, :, :) ! (iip1, jjp1, llm)
34      real vcov(iip1, jjm, llm)      real, intent(in):: vcov(iip1, jjm, llm)
35      real, intent(in):: trac(:, :, :) ! (iim + 1, jjm + 1, llm)      real, intent(in):: trac(:, :, :) ! (iim + 1, jjm + 1, llm)
36    
37      ! Local:      ! Local:
38    
39      integer:: icum  = 0      integer:: icum  = 0
     integer, save:: ncum  
     logical:: first = .true.  
     real zz, zqy, zfactv(jjm, llm)  
   
     integer, parameter:: nQ = 7  
     character(len=4), parameter:: nom(nQ) = (/'T   ', 'gz  ', 'K   ', 'ang ', &  
          'u   ', 'ovap', 'un  '/)  
     character(len=5), parameter:: unites(nQ) = (/'K    ', 'm2/s2', 'm2/s2', &  
          'ang  ', 'm/s  ', 'kg/kg', 'un   '/)  
   
     real:: time = 0.  
40      integer:: itau = 0      integer:: itau = 0
41      real ww      real qy, factv(jjm, llm)
42    
43      ! Variables dynamiques intermédiaires      ! Variables dynamiques interm\'ediaires
44      REAL vcont(iip1, jjm, llm), ucont(iip1, jjp1, llm)      REAL vcont(iip1, jjm, llm), ucont(iip1, jjp1, llm)
45      REAL ang(iip1, jjp1, llm), unat(iip1, jjp1, llm)      REAL ang(iip1, jjp1, llm), unat(iip1, jjp1, llm)
46      REAL massebx(iip1, jjp1, llm), masseby(iip1, jjm, llm)      REAL massebx(iip1, jjp1, llm), masseby(iip1, jjm, llm)
47      REAL w(iip1, jjp1, llm), ecin(iip1, jjp1, llm), convm(iip1, jjp1, llm)      REAL ecin(iip1, jjp1, llm)
48    
49      ! Champ contenant les scalaires advectés      ! Champ contenant les scalaires advect\'es
50      real Q(iip1, jjp1, llm, nQ)      real Q(iip1, jjp1, llm, nQ)
51    
52      ! Champs cumulés      ! Champs cumul\'es
53      real, save:: ps_cum(iip1, jjp1)      real, save:: ps_cum(iip1, jjp1)
54      real, save:: masse_cum(iip1, jjp1, llm)      real, save:: masse_cum(iip1, jjp1, llm)
55      real, save:: flux_u_cum(iip1, jjp1, llm)      real, save:: flux_u_cum(iip1, jjp1, llm)
# Line 73  contains Line 57  contains
57      real, save:: Q_cum(iip1, jjp1, llm, nQ)      real, save:: Q_cum(iip1, jjp1, llm, nQ)
58      real, save:: flux_uQ_cum(iip1, jjp1, llm, nQ)      real, save:: flux_uQ_cum(iip1, jjp1, llm, nQ)
59      real, save:: flux_vQ_cum(iip1, jjm, llm, nQ)      real, save:: flux_vQ_cum(iip1, jjm, llm, nQ)
     real dQ(iip1, jjp1, llm, nQ)  
60    
61      ! champs de tansport en moyenne zonale      ! champs de tansport en moyenne zonale
62      integer itr      integer itr
     integer, parameter:: ntr = 5  
   
     character(len=10), save:: znom(ntr, nQ)  
     character(len=26), save:: znoml(ntr, nQ)  
     character(len=12), save:: zunites(ntr, nQ)  
   
63      integer, parameter:: iave = 1, itot = 2, immc = 3, itrs = 4, istn = 5      integer, parameter:: iave = 1, itot = 2, immc = 3, itrs = 4, istn = 5
     character(len=3), parameter:: ctrs(ntr) = (/'   ', 'TOT', 'MMC', 'TRS', &  
          'STN'/)  
64    
65      real zvQ(jjm, llm, ntr, nQ), zvQtmp(jjm, llm)      real vq(jjm, llm, ntr, nQ), vqtmp(jjm, llm)
66      real zavQ(jjm, 2: ntr, nQ), psiQ(jjm, llm + 1, nQ)      real avq(jjm, 2: ntr, nQ), psiQ(jjm, llm + 1, nQ)
67      real zmasse(jjm, llm)      real zmasse(jjm, llm)
68        real v(jjm, llm), psi(jjm, llm + 1)
     real zv(jjm, llm), psi(jjm, llm + 1)  
   
69      integer i, j, l, iQ      integer i, j, l, iQ
70    
     ! Initialisation du fichier contenant les moyennes zonales.  
   
     integer, save:: fileid  
     integer thoriid, zvertiid  
   
     real zjulian  
     integer zan, dayref  
   
     real rlong(jjm), rlatg(jjm)  
   
71      !-----------------------------------------------------------------      !-----------------------------------------------------------------
72    
     !!print *, "Call sequence information: bilan_dyn"  
   
     ! Initialisation  
   
     time = time + dt_app  
     itau = itau + 1  
   
     first_call: if (first) then  
        ! initialisation des fichiers  
        first = .false.  
        ! ncum est la frequence de stokage en pas de temps  
        ncum = dt_cum / dt_app  
        if (abs(ncum * dt_app - dt_cum) > 1e-5 * dt_app) then  
           print *, 'Problème : le pas de cumul doit être multiple du pas'  
           print *, 'dt_app = ', dt_app  
           print *, 'dt_cum = ', dt_cum  
           stop 1  
        endif  
   
        call inigrads(i_f=4, x=(/0./), fx=180./pi, xmin=0., xmax=0., y=rlatv, &  
             ymin=-90., ymax=90., fy=180./pi, z=presnivs, fz=1., dt=dt_cum, &  
             file='dynzon', titlel='dyn_zon ')  
   
        ! Initialisation du fichier contenant les moyennes zonales  
   
        zan = annee_ref  
        dayref = day_ref  
        CALL ymds2ju(zan, 1, dayref, 0.0, zjulian)  
   
        rlong = 0.  
        rlatg = rlatv*180./pi  
   
        call histbeg_totreg('dynzon', rlong(:1), rlatg, 1, 1, 1, jjm, itau_dyn, &  
             zjulian, dt_cum, thoriid, fileid)  
   
        ! Appel à histvert pour la grille verticale  
   
        call histvert(fileid, 'presnivs', 'Niveaux sigma', 'mb', llm, presnivs, &  
             zvertiid)  
   
        ! Appels à histdef pour la définition des variables à sauvegarder  
        do iQ = 1, nQ  
           do itr = 1, ntr  
              if (itr == 1) then  
                 znom(itr, iQ) = nom(iQ)  
                 znoml(itr, iQ) = nom(iQ)  
                 zunites(itr, iQ) = unites(iQ)  
              else  
                 znom(itr, iQ) = ctrs(itr)//'v'//nom(iQ)  
                 znoml(itr, iQ) = 'transport : v * '//nom(iQ)//' '//ctrs(itr)  
                 zunites(itr, iQ) = 'm/s * '//unites(iQ)  
              endif  
           enddo  
        enddo  
   
        ! Déclarations des champs avec dimension verticale  
        do iQ = 1, nQ  
           do itr = 1, ntr  
              call histdef(fileid, znom(itr, iQ), znoml(itr, iQ), &  
                   zunites(itr, iQ), 1, jjm, thoriid, llm, 1, llm, zvertiid, &  
                   'ave(X)', dt_cum, dt_cum)  
           enddo  
           ! Declarations pour les fonctions de courant  
           call histdef(fileid, 'psi'//nom(iQ), 'stream fn. '//znoml(itot, iQ), &  
                zunites(itot, iQ), 1, jjm, thoriid, llm, 1, llm, zvertiid, &  
                'ave(X)', dt_cum, dt_cum)  
        enddo  
   
        ! Declarations pour les champs de transport d'air  
        call histdef(fileid, 'masse', 'masse', &  
             'kg', 1, jjm, thoriid, llm, 1, llm, zvertiid, &  
             'ave(X)', dt_cum, dt_cum)  
        call histdef(fileid, 'v', 'v', &  
             'm/s', 1, jjm, thoriid, llm, 1, llm, zvertiid, &  
             'ave(X)', dt_cum, dt_cum)  
        ! Declarations pour les fonctions de courant  
        call histdef(fileid, 'psi', 'stream fn. MMC ', 'mega t/s', &  
             1, jjm, thoriid, llm, 1, llm, zvertiid, &  
             'ave(X)', dt_cum, dt_cum)  
   
        ! Declaration des champs 1D de transport en latitude  
        do iQ = 1, nQ  
           do itr = 2, ntr  
              call histdef(fileid, 'a'//znom(itr, iQ), znoml(itr, iQ), &  
                   zunites(itr, iQ), 1, jjm, thoriid, 1, 1, 1, -99, &  
                   'ave(X)', dt_cum, dt_cum)  
           enddo  
        enddo  
   
        CALL histend(fileid)  
     endif first_call  
   
73      ! Calcul des champs dynamiques      ! Calcul des champs dynamiques
74    
75      ! Énergie cinétique      ! \'Energie cin\'etique
76      ucont = 0      ucont = 0
77      CALL covcont(llm, ucov, vcov, ucont, vcont)      CALL covcont(llm, ucov, vcov, ucont, vcont)
78      CALL enercin(vcov, ucov, vcont, ucont, ecin)      CALL enercin(vcov, ucov, vcont, ucont, ecin)
79    
80      ! moment cinétique      ! moment cin\'etique
81      do l = 1, llm      forall (l = 1: llm)
82         ang(:, :, l) = ucov(:, :, l) + constang_2d         ang(:, :, l) = ucov(:, :, l) + constang_2d
83         unat(:, :, l) = ucont(:, :, l)*cu_2d         unat(:, :, l) = ucont(:, :, l) * cu_2d
84      enddo      end forall
85    
86      Q(:, :, :, 1) = teta * pk / cpp      Q(:, :, :, 1) = teta * pk / cpp
87      Q(:, :, :, 2) = phi      Q(:, :, :, 2) = phi
# Line 232  contains Line 103  contains
103         flux_uQ_cum = 0.         flux_uQ_cum = 0.
104      endif      endif
105    
106        itau = itau + 1
107      icum = icum + 1      icum = icum + 1
108    
109      ! Accumulation des flux de masse horizontaux      ! Accumulation des flux de masse horizontaux
# Line 239  contains Line 111  contains
111      masse_cum = masse_cum + masse      masse_cum = masse_cum + masse
112      flux_u_cum = flux_u_cum + flux_u      flux_u_cum = flux_u_cum + flux_u
113      flux_v_cum = flux_v_cum + flux_v      flux_v_cum = flux_v_cum + flux_v
114      do iQ = 1, nQ      forall (iQ = 1: nQ) Q_cum(:, :, :, iQ) = Q_cum(:, :, :, iQ) &
115         Q_cum(:, :, :, iQ) = Q_cum(:, :, :, iQ) + Q(:, :, :, iQ)*masse           + Q(:, :, :, iQ) * masse
     enddo  
   
     ! FLUX ET TENDANCES  
116    
117      ! Flux longitudinal      ! Flux longitudinal
118      forall (iQ = 1: nQ, i = 1: iim) flux_uQ_cum(i, :, :, iQ) &      forall (iQ = 1: nQ, i = 1: iim) flux_uQ_cum(i, :, :, iQ) &
# Line 251  contains Line 120  contains
120           + flux_u(i, :, :) * 0.5 * (Q(i, :, :, iQ) + Q(i + 1, :, :, iQ))           + flux_u(i, :, :) * 0.5 * (Q(i, :, :, iQ) + Q(i + 1, :, :, iQ))
121      flux_uQ_cum(iip1, :, :, :) = flux_uQ_cum(1, :, :, :)      flux_uQ_cum(iip1, :, :, :) = flux_uQ_cum(1, :, :, :)
122    
123      ! Flux méridien      ! Flux m\'eridien
124      forall (iQ = 1: nQ, j = 1: jjm) flux_vQ_cum(:, j, :, iQ) &      forall (iQ = 1: nQ, j = 1: jjm) flux_vQ_cum(:, j, :, iQ) &
125           = flux_vQ_cum(:, j, :, iQ) &           = flux_vQ_cum(:, j, :, iQ) &
126           + flux_v(:, j, :) * 0.5 * (Q(:, j, :, iQ) + Q(:, j + 1, :, iQ))           + flux_v(:, j, :) * 0.5 * (Q(:, j, :, iQ) + Q(:, j + 1, :, iQ))
127    
     ! tendances  
   
     ! convergence horizontale  
     call convflu(flux_uQ_cum, flux_vQ_cum, llm*nQ, dQ)  
   
     ! calcul de la vitesse verticale  
     call convmas(flux_u_cum, flux_v_cum, convm)  
     CALL vitvert(convm, w)  
   
     do iQ = 1, nQ  
        do l = 1, llm-1  
           do j = 1, jjp1  
              do i = 1, iip1  
                 ww = -0.5*w(i, j, l + 1)*(Q(i, j, l, iQ) + Q(i, j, l + 1, iQ))  
                 dQ(i, j, l, iQ) = dQ(i, j, l, iQ)-ww  
                 dQ(i, j, l + 1, iQ) = dQ(i, j, l + 1, iQ) + ww  
              enddo  
           enddo  
        enddo  
     enddo  
   
     ! PAS DE TEMPS D'ECRITURE  
   
128      writing_step: if (icum == ncum) then      writing_step: if (icum == ncum) then
129         ! Normalisation         ! Normalisation
130         do iQ = 1, nQ         forall (iQ = 1: nQ) Q_cum(:, :, :, iQ) = Q_cum(:, :, :, iQ) / masse_cum
131            Q_cum(:, :, :, iQ) = Q_cum(:, :, :, iQ)/masse_cum         ps_cum = ps_cum / ncum
132         enddo         masse_cum = masse_cum / ncum
133         zz = 1. / real(ncum)         flux_u_cum = flux_u_cum / ncum
134         ps_cum = ps_cum*zz         flux_v_cum = flux_v_cum / ncum
135         masse_cum = masse_cum*zz         flux_uQ_cum = flux_uQ_cum / ncum
136         flux_u_cum = flux_u_cum*zz         flux_vQ_cum = flux_vQ_cum / ncum
        flux_v_cum = flux_v_cum*zz  
        flux_uQ_cum = flux_uQ_cum*zz  
        flux_vQ_cum = flux_vQ_cum*zz  
        dQ = dQ*zz  
   
        ! A retravailler eventuellement  
        ! division de dQ par la masse pour revenir aux bonnes grandeurs  
        do iQ = 1, nQ  
           dQ(:, :, :, iQ) = dQ(:, :, :, iQ)/masse_cum  
        enddo  
137    
138         ! Transport méridien         ! Transport m\'eridien
139    
140         ! cumul zonal des masses des mailles         ! Cumul zonal des masses des mailles
141    
142         zv = 0.         v = 0.
143         zmasse = 0.         zmasse = 0.
144         call massbar(masse_cum, massebx, masseby)         call massbar(masse_cum, massebx, masseby)
145         do l = 1, llm         do l = 1, llm
146            do j = 1, jjm            do j = 1, jjm
147               do i = 1, iim               do i = 1, iim
148                  zmasse(j, l) = zmasse(j, l) + masseby(i, j, l)                  zmasse(j, l) = zmasse(j, l) + masseby(i, j, l)
149                  zv(j, l) = zv(j, l) + flux_v_cum(i, j, l)                  v(j, l) = v(j, l) + flux_v_cum(i, j, l)
150               enddo               enddo
151               zfactv(j, l) = cv_2d(1, j)/zmasse(j, l)               factv(j, l) = cv_2d(1, j) / zmasse(j, l)
152            enddo            enddo
153         enddo         enddo
154    
155         ! Transport dans le plan latitude-altitude         ! Transport dans le plan latitude-altitude
156    
157         zvQ = 0.         vq = 0.
158         psiQ = 0.         psiQ = 0.
159         do iQ = 1, nQ         do iQ = 1, nQ
160            zvQtmp = 0.            vqtmp = 0.
161            do l = 1, llm            do l = 1, llm
162               do j = 1, jjm               do j = 1, jjm
163                  ! Calcul des moyennes zonales du transort total et de zvQtmp                  ! Calcul des moyennes zonales du transport total et de vqtmp
164                  do i = 1, iim                  do i = 1, iim
165                     zvQ(j, l, itot, iQ) = zvQ(j, l, itot, iQ) &                     vq(j, l, itot, iQ) = vq(j, l, itot, iQ) &
166                          + flux_vQ_cum(i, j, l, iQ)                          + flux_vQ_cum(i, j, l, iQ)
167                     zqy =  0.5 * (Q_cum(i, j, l, iQ) * masse_cum(i, j, l) &                     qy =  0.5 * (Q_cum(i, j, l, iQ) * masse_cum(i, j, l) &
168                          + Q_cum(i, j + 1, l, iQ) * masse_cum(i, j + 1, l))                          + Q_cum(i, j + 1, l, iQ) * masse_cum(i, j + 1, l))
169                     zvQtmp(j, l) = zvQtmp(j, l) + flux_v_cum(i, j, l) * zqy &                     vqtmp(j, l) = vqtmp(j, l) + flux_v_cum(i, j, l) * qy &
170                          / (0.5 * (masse_cum(i, j, l) + masse_cum(i, j + 1, l)))                          / (0.5 * (masse_cum(i, j, l) + masse_cum(i, j + 1, l)))
171                     zvQ(j, l, iave, iQ) = zvQ(j, l, iave, iQ) + zqy                     vq(j, l, iave, iQ) = vq(j, l, iave, iQ) + qy
172                  enddo                  enddo
173                  ! Decomposition                  ! Decomposition
174                  zvQ(j, l, iave, iQ) = zvQ(j, l, iave, iQ)/zmasse(j, l)                  vq(j, l, iave, iQ) = vq(j, l, iave, iQ) / zmasse(j, l)
175                  zvQ(j, l, itot, iQ) = zvQ(j, l, itot, iQ)*zfactv(j, l)                  vq(j, l, itot, iQ) = vq(j, l, itot, iQ) * factv(j, l)
176                  zvQtmp(j, l) = zvQtmp(j, l)*zfactv(j, l)                  vqtmp(j, l) = vqtmp(j, l) * factv(j, l)
177                  zvQ(j, l, immc, iQ) = zv(j, l)*zvQ(j, l, iave, iQ)*zfactv(j, l)                  vq(j, l, immc, iQ) = v(j, l) * vq(j, l, iave, iQ) * factv(j, l)
178                  zvQ(j, l, itrs, iQ) = zvQ(j, l, itot, iQ)-zvQtmp(j, l)                  vq(j, l, itrs, iQ) = vq(j, l, itot, iQ) - vqtmp(j, l)
179                  zvQ(j, l, istn, iQ) = zvQtmp(j, l)-zvQ(j, l, immc, iQ)                  vq(j, l, istn, iQ) = vqtmp(j, l) - vq(j, l, immc, iQ)
180               enddo               enddo
181            enddo            enddo
182            ! fonction de courant meridienne pour la quantite Q            ! Fonction de courant m\'eridienne pour la quantit\'e Q
183            do l = llm, 1, -1            do l = llm, 1, -1
184               do j = 1, jjm               do j = 1, jjm
185                  psiQ(j, l, iQ) = psiQ(j, l + 1, iQ) + zvQ(j, l, itot, iQ)                  psiQ(j, l, iQ) = psiQ(j, l + 1, iQ) + vq(j, l, itot, iQ)
186               enddo               enddo
187            enddo            enddo
188         enddo         enddo
189    
190         ! fonction de courant pour la circulation meridienne moyenne         ! Fonction de courant pour la circulation m\'eridienne moyenne
191         psi = 0.         psi = 0.
192         do l = llm, 1, -1         do l = llm, 1, -1
193            do j = 1, jjm            do j = 1, jjm
194               psi(j, l) = psi(j, l + 1) + zv(j, l)               psi(j, l) = psi(j, l + 1) + v(j, l)
195               zv(j, l) = zv(j, l)*zfactv(j, l)               v(j, l) = v(j, l) * factv(j, l)
196            enddo            enddo
197         enddo         enddo
198    
199         ! sorties proprement dites         ! Sorties proprement dites
200         do iQ = 1, nQ         do iQ = 1, nQ
201            do itr = 1, ntr            do itr = 1, ntr
202               call histwrite(fileid, znom(itr, iQ), itau, zvQ(:, :, itr, iQ))               call histwrite(fileid, znom(itr, iQ), itau, vq(:, :, itr, iQ))
203            enddo            enddo
204            call histwrite(fileid, 'psi'//nom(iQ), itau, psiQ(:, :llm, iQ))            call histwrite(fileid, 'psi' // nom(iQ), itau, psiQ(:, :llm, iQ))
205         enddo         enddo
206    
207         call histwrite(fileid, 'masse', itau, zmasse)         call histwrite(fileid, 'masse', itau, zmasse)
208         call histwrite(fileid, 'v', itau, zv)         call histwrite(fileid, 'v', itau, v)
209         psi = psi*1.e-9         psi = psi * 1e-9
210         call histwrite(fileid, 'psi', itau, psi(:, :llm))         call histwrite(fileid, 'psi', itau, psi(:, :llm))
211    
212         ! Intégrale verticale         ! Int\'egrale verticale
213    
214         forall (iQ = 1: nQ, itr = 2: ntr) zavQ(:, itr, iQ) &         forall (iQ = 1: nQ, itr = 2: ntr) avq(:, itr, iQ) &
215              = sum(zvQ(:, :, itr, iQ) * zmasse, dim=2) / cv_2d(1, :)              = sum(vq(:, :, itr, iQ) * zmasse, dim=2) / cv_2d(1, :)
216    
217         do iQ = 1, nQ         do iQ = 1, nQ
218            do itr = 2, ntr            do itr = 2, ntr
219               call histwrite(fileid, 'a'//znom(itr, iQ), itau, zavQ(:, itr, iQ))               call histwrite(fileid, 'a' // znom(itr, iQ), itau, avq(:, itr, iQ))
220            enddo            enddo
221         enddo         enddo
222    
        ! On doit pouvoir tracer systematiquement la fonction de courant.  
223         icum = 0         icum = 0
224      endif writing_step      endif writing_step
225    
Legend:
Removed from v.55  
changed lines
  Added in v.254
  ViewVC Help
Powered by ViewVC 1.1.21