/[lmdze]/trunk/libf/dyn3d/leapfrog.f90
ViewVC logotype

Diff of /trunk/libf/dyn3d/leapfrog.f90

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

revision 29 by guez, Tue Mar 30 10:44:42 2010 UTC revision 30 by guez, Thu Apr 1 09:07:28 2010 UTC
# Line 8  contains Line 8  contains
8    
9      ! From dyn3d/leapfrog.F, version 1.6, 2005/04/13 08:58:34      ! From dyn3d/leapfrog.F, version 1.6, 2005/04/13 08:58:34
10      ! Authors: P. Le Van, L. Fairhead, F. Hourdin      ! Authors: P. Le Van, L. Fairhead, F. Hourdin
11        ! schema matsuno + leapfrog
12    
13      USE calfis_m, ONLY: calfis      USE calfis_m, ONLY: calfis
14      USE com_io_dyn, ONLY: histaveid      USE com_io_dyn, ONLY: histaveid
# Line 77  contains Line 78  contains
78      REAL time ! time of day, as a fraction of day length      REAL time ! time of day, as a fraction of day length
79      real finvmaold(ip1jmp1, llm)      real finvmaold(ip1jmp1, llm)
80      LOGICAL:: lafin=.false.      LOGICAL:: lafin=.false.
81      INTEGER l      INTEGER i, j, l
82    
83      REAL rdayvrai, rdaym_ini      REAL rdayvrai, rdaym_ini
84    
# Line 96  contains Line 97  contains
97      print *, "Call sequence information: leapfrog"      print *, "Call sequence information: leapfrog"
98    
99      itaufin = nday * day_step      itaufin = nday * day_step
100        ! "day_step" is a multiple of "iperiod", therefore "itaufin" is one too
101    
102      itau = 0      itau = 0
103      iday = day_ini      iday = day_ini
104      time = time_0      time = time_0
# Line 105  contains Line 108  contains
108      CALL exner_hyb(ps, p3d, pks, pk, pkf)      CALL exner_hyb(ps, p3d, pks, pk, pkf)
109    
110      ! Début de l'integration temporelle :      ! Début de l'integration temporelle :
111      outer_loop:do      outer_loop:do i = 1, itaufin / iperiod
112           ! {itau is a multiple of iperiod}
113    
114           ! 1. Matsuno forward:
115    
116         if (ok_guide .and. (itaufin - itau - 1) * dtvr > 21600.) &         if (ok_guide .and. (itaufin - itau - 1) * dtvr > 21600.) &
117              call guide(itau, ucov, vcov, teta, q, masse, ps)              call guide(itau, ucov, vcov, teta, q, masse, ps)
118         vcovm1 = vcov         vcovm1 = vcov
# Line 113  contains Line 120  contains
120         tetam1 = teta         tetam1 = teta
121         massem1 = masse         massem1 = masse
122         psm1 = ps         psm1 = ps
        forward = .TRUE.  
        leapf = .FALSE.  
        dt = dtvr  
123         finvmaold = masse         finvmaold = masse
124         CALL filtreg(finvmaold, jjm + 1, llm, - 2, 2, .TRUE., 1)         CALL filtreg(finvmaold, jjm + 1, llm, - 2, 2, .TRUE., 1)
125    
126         do         ! Calcul des tendances dynamiques:
127           CALL geopot(ip1jmp1, teta, pk, pks, phis, phi)
128           CALL caldyn(itau, ucov, vcov, teta, ps, masse, pk, pkf, phis, phi, &
129                MOD(itau, iconser) == 0, du, dv, dteta, dp, w, pbaru, pbarv, &
130                time + iday - day_ini)
131    
132           ! Calcul des tendances advection des traceurs (dont l'humidité)
133           CALL caladvtrac(q, pbaru, pbarv, p3d, masse, dq, teta, pk)
134           ! Stokage du flux de masse pour traceurs offline:
135           IF (offline) CALL fluxstokenc(pbaru, pbarv, masse, teta, phi, phis, &
136                dtvr, itau)
137    
138           ! integrations dynamique et traceurs:
139           CALL integrd(2, vcovm1, ucovm1, tetam1, psm1, massem1, dv, du, dteta, &
140                dq, dp, vcov, ucov, teta, q, ps, masse, phis, finvmaold, .false., &
141                dtvr)
142    
143           CALL pression(ip1jmp1, ap, bp, ps, p3d)
144           CALL exner_hyb(ps, p3d, pks, pk, pkf)
145    
146           ! 2. Matsuno backward:
147    
148           itau = itau + 1
149           iday = day_ini + itau / day_step
150           time = REAL(itau - (iday - day_ini) * day_step) / day_step + time_0
151           IF (time > 1.) THEN
152              time = time - 1.
153              iday = iday + 1
154           ENDIF
155    
156           ! Calcul des tendances dynamiques:
157           CALL geopot(ip1jmp1, teta, pk, pks, phis, phi)
158           CALL caldyn(itau, ucov, vcov, teta, ps, masse, pk, pkf, phis, phi, &
159                .false., du, dv, dteta, dp, w, pbaru, pbarv, time + iday - day_ini)
160    
161           ! integrations dynamique et traceurs:
162           CALL integrd(2, vcovm1, ucovm1, tetam1, psm1, massem1, dv, du, dteta, &
163                dq, dp, vcov, ucov, teta, q, ps, masse, phis, finvmaold, .false., &
164                dtvr)
165    
166           CALL pression(ip1jmp1, ap, bp, ps, p3d)
167           CALL exner_hyb(ps, p3d, pks, pk, pkf)
168    
169           ! 3. Leapfrog:
170    
171           do j = 1, iperiod - 1
172            ! Calcul des tendances dynamiques:            ! Calcul des tendances dynamiques:
173            CALL geopot(ip1jmp1, teta, pk, pks, phis, phi)            CALL geopot(ip1jmp1, teta, pk, pks, phis, phi)
174            CALL caldyn(itau, ucov, vcov, teta, ps, masse, pk, pkf, phis, phi, &            CALL caldyn(itau, ucov, vcov, teta, ps, masse, pk, pkf, phis, phi, &
175                 MOD(itau, iconser) == 0, du, dv, dteta, dp, w, pbaru, pbarv, &                 .false., du, dv, dteta, dp, w, pbaru, pbarv, &
176                 time + iday - day_ini)                 time + iday - day_ini)
177    
178            IF (forward .OR. leapf) THEN            ! Calcul des tendances advection des traceurs (dont l'humidité)
179               ! Calcul des tendances advection des traceurs (dont l'humidité)            CALL caladvtrac(q, pbaru, pbarv, p3d, masse, dq, teta, pk)
180               CALL caladvtrac(q, pbaru, pbarv, p3d, masse, dq, teta, pk)            ! Stokage du flux de masse pour traceurs off-line:
181               IF (offline) THEN            IF (offline) CALL fluxstokenc(pbaru, pbarv, masse, teta, phi, phis, &
182                  ! Stokage du flux de masse pour traceurs off-line                 dtvr, itau)
                 CALL fluxstokenc(pbaru, pbarv, masse, teta, phi, phis, dtvr, &  
                      itau)  
              ENDIF  
           ENDIF  
183    
184            ! integrations dynamique et traceurs:            ! integrations dynamique et traceurs:
185            CALL integrd(2, vcovm1, ucovm1, tetam1, psm1, massem1, dv, du, &            CALL integrd(2, vcovm1, ucovm1, tetam1, psm1, massem1, dv, du, &
186                 dteta, dq, dp, vcov, ucov, teta, q, ps, masse, phis, &                 dteta, dq, dp, vcov, ucov, teta, q, ps, masse, phis, &
187                 finvmaold, leapf, dt)                 finvmaold, .true., 2 * dtvr)
188    
189            IF (MOD(itau + 1, iphysiq) == 0 .AND. iflag_phys /= 0) THEN            IF (MOD(itau + 1, iphysiq) == 0 .AND. iflag_phys /= 0) THEN
190               ! calcul des tendances physiques:               ! calcul des tendances physiques:
# Line 196  contains Line 241  contains
241                    / apols                    / apols
242            END IF            END IF
243    
244            ! fin de l'intégration dynamique et physique pour le pas "itau"            itau = itau + 1
245            ! préparation du pas d'intégration suivant            iday = day_ini + itau / day_step
246              time = REAL(itau - (iday - day_ini) * day_step) / day_step + time_0
247            ! schema matsuno + leapfrog            IF (time > 1.) THEN
248            IF (forward .OR. leapf) THEN               time = time - 1.
249               itau = itau + 1               iday = iday + 1
              iday = day_ini + itau / day_step  
              time = REAL(itau - (iday - day_ini) * day_step) / day_step &  
                   + time_0  
              IF (time > 1.) THEN  
                 time = time - 1.  
                 iday = iday + 1  
              ENDIF  
250            ENDIF            ENDIF
251    
252            IF (itau == itaufin + 1) exit outer_loop            IF (MOD(itau, iperiod) == 0) THEN
   
           IF (MOD(itau, iperiod) == 0 .OR. itau == itaufin) THEN  
253               ! ecriture du fichier histoire moyenne:               ! ecriture du fichier histoire moyenne:
254               CALL writedynav(histaveid, nqmx, itau, vcov, &               CALL writedynav(histaveid, nqmx, itau, vcov, &
255                    ucov, teta, pk, phi, q, masse, ps, phis)                    ucov, teta, pk, phi, q, masse, ps, phis)
256               call bilan_dyn(2, dtvr * iperiod, dtvr * day_step * periodav, &               call bilan_dyn(2, dtvr * iperiod, dtvr * day_step * periodav, &
257                    ps, masse, pk, pbaru, pbarv, teta, phi, ucov, vcov, q)                    ps, masse, pk, pbaru, pbarv, teta, phi, ucov, vcov, q)
258            ENDIF            ENDIF
   
           IF (itau == itaufin) CALL dynredem1("restart.nc", vcov, ucov, teta, &  
                q, masse, ps, itau=itau_dyn+itaufin)  
   
           ! gestion de l'integration temporelle:  
           IF (MOD(itau, iperiod) == 0) exit  
           IF (MOD(itau - 1, iperiod) == 0) THEN  
              IF (forward) THEN  
                 ! fin du pas forward et debut du pas backward  
                 forward = .FALSE.  
                 leapf = .FALSE.  
              ELSE  
                 ! fin du pas backward et debut du premier pas leapfrog  
                 leapf = .TRUE.  
                 dt = 2. * dtvr  
              END IF  
           ELSE  
              ! pas leapfrog  
              leapf = .TRUE.  
              dt = 2. * dtvr  
           END IF  
259         end do         end do
260      end do outer_loop      end do outer_loop
261    
262        ! {itau == itaufin}
263        CALL dynredem1("restart.nc", vcov, ucov, teta, q, masse, ps, &
264             itau=itau_dyn+itaufin)
265    
266        vcovm1 = vcov
267        ucovm1 = ucov
268        tetam1 = teta
269        massem1 = masse
270        psm1 = ps
271        finvmaold = masse
272        CALL filtreg(finvmaold, jjm + 1, llm, - 2, 2, .TRUE., 1)
273    
274        ! Calcul des tendances dynamiques:
275        CALL geopot(ip1jmp1, teta, pk, pks, phis, phi)
276        CALL caldyn(itaufin, ucov, vcov, teta, ps, masse, pk, pkf, phis, phi, &
277             MOD(itaufin, iconser) == 0, du, dv, dteta, dp, w, pbaru, pbarv, &
278             time + iday - day_ini)
279    
280        ! Calcul des tendances advection des traceurs (dont l'humidité)
281        CALL caladvtrac(q, pbaru, pbarv, p3d, masse, dq, teta, pk)
282        ! Stokage du flux de masse pour traceurs off-line:
283        IF (offline) CALL fluxstokenc(pbaru, pbarv, masse, teta, phi, phis, dtvr, &
284             itaufin)
285    
286        ! integrations dynamique et traceurs:
287        CALL integrd(2, vcovm1, ucovm1, tetam1, psm1, massem1, dv, du, dteta, dq, &
288             dp, vcov, ucov, teta, q, ps, masse, phis, finvmaold, .false., dtvr)
289        
290        CALL pression(ip1jmp1, ap, bp, ps, p3d)
291        CALL exner_hyb(ps, p3d, pks, pk, pkf)
292    
293    END SUBROUTINE leapfrog    END SUBROUTINE leapfrog
294    
295  end module leapfrog_m  end module leapfrog_m

Legend:
Removed from v.29  
changed lines
  Added in v.30

  ViewVC Help
Powered by ViewVC 1.1.21