/[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 26 by guez, Tue Mar 9 15:27:15 2010 UTC revision 28 by guez, Fri Mar 26 18:33:04 2010 UTC
# Line 6  contains Line 6  contains
6    
7    SUBROUTINE leapfrog(ucov, vcov, teta, ps, masse, phis, q, time_0)    SUBROUTINE leapfrog(ucov, vcov, teta, ps, masse, phis, q, time_0)
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      ! Auteurs: P. Le Van, L. Fairhead, F. Hourdin      ! Authors: P. Le Van, L. Fairhead, F. Hourdin
11    
12      USE calfis_m, ONLY: calfis      USE calfis_m, ONLY: calfis
13      USE com_io_dyn, ONLY: histaveid      USE com_io_dyn, ONLY: histaveid
14      USE comconst, ONLY: daysec, dtphys, dtvr      USE comconst, ONLY: daysec, dtphys, dtvr
15      USE comgeom, ONLY: aire, apoln, apols      USE comgeom, ONLY: aire, apoln, apols
16      USE comvert, ONLY: ap, bp      USE comvert, ONLY: ap, bp
17      USE conf_gcm_m, ONLY: day_step, iconser, iperiod, iphysiq, &      USE conf_gcm_m, ONLY: day_step, iconser, iperiod, iphysiq, nday, offline, &
18           nday, offline, periodav           periodav
19      USE dimens_m, ONLY: iim, llm, nqmx      USE dimens_m, ONLY: iim, llm, nqmx
20      USE dynetat0_m, ONLY: day_ini      USE dynetat0_m, ONLY: day_ini
21        use dynredem1_m, only: dynredem1
22      USE exner_hyb_m, ONLY: exner_hyb      USE exner_hyb_m, ONLY: exner_hyb
23        use filtreg_m, only: filtreg
24      USE guide_m, ONLY: guide      USE guide_m, ONLY: guide
25      use inidissip_m, only: idissip      use inidissip_m, only: idissip
26      USE logic, ONLY: iflag_phys, ok_guide      USE logic, ONLY: iflag_phys, ok_guide
27      USE paramet_m, ONLY: iip1, ip1jm, ip1jmp1, jjp1      USE paramet_m, ONLY: iip1, ip1jm, ip1jmp1, jjp1
28      USE pression_m, ONLY: pression      USE pression_m, ONLY: pression
29      USE pressure_var, ONLY: p3d      USE pressure_var, ONLY: p3d
30      USE temps, ONLY: dt, itaufin      USE temps, ONLY: itau_dyn
31    
32      ! Variables dynamiques:      ! Variables dynamiques:
33      REAL vcov(ip1jm, llm), ucov(ip1jmp1, llm) ! vents covariants      REAL vcov(ip1jm, llm), ucov(ip1jmp1, llm) ! vents covariants
# Line 72  contains Line 74  contains
74      REAL tppn(iim), tpps(iim), tpn, tps      REAL tppn(iim), tpps(iim), tpn, tps
75    
76      INTEGER itau ! index of the time step of the dynamics, starts at 0      INTEGER itau ! index of the time step of the dynamics, starts at 0
77        INTEGER itaufin
78      INTEGER iday ! jour julien      INTEGER iday ! jour julien
79      REAL time ! time of day, as a fraction of day length      REAL time ! time of day, as a fraction of day length
80      real finvmaold(ip1jmp1, llm)      real finvmaold(ip1jmp1, llm)
# Line 87  contains Line 90  contains
90      ! cree par la dissipation      ! cree par la dissipation
91      REAL dtetaecdt(ip1jmp1, llm)      REAL dtetaecdt(ip1jmp1, llm)
92      REAL vcont(ip1jm, llm), ucont(ip1jmp1, llm)      REAL vcont(ip1jm, llm), ucont(ip1jmp1, llm)
93      CHARACTER*15 ztit      logical forward, leapf
94      INTEGER:: ip_ebil_dyn = 0 ! PRINT level for energy conserv. diag.      REAL dt
   
     logical:: dissip_conservative = .true.  
     logical forward, leapf, apphys, conser, apdiss  
95    
96      !---------------------------------------------------      !---------------------------------------------------
97    
# Line 106  contains Line 106  contains
106      CALL pression(ip1jmp1, ap, bp, ps, p3d)      CALL pression(ip1jmp1, ap, bp, ps, p3d)
107      CALL exner_hyb(ps, p3d, pks, pk, pkf)      CALL exner_hyb(ps, p3d, pks, pk, pkf)
108    
109      ! Debut de l'integration temporelle:      ! Début de l'integration temporelle :
110      outer_loop:do      outer_loop:do
111         if (ok_guide .and. (itaufin - itau - 1) * dtvr > 21600.) &         if (ok_guide .and. (itaufin - itau - 1) * dtvr > 21600.) &
112              call guide(itau, ucov, vcov, teta, q, masse, ps)              call guide(itau, ucov, vcov, teta, q, masse, ps)
# Line 122  contains Line 122  contains
122         CALL filtreg(finvmaold, jjp1, llm, - 2, 2, .TRUE., 1)         CALL filtreg(finvmaold, jjp1, llm, - 2, 2, .TRUE., 1)
123    
124         do         do
125            ! gestion des appels de la physique et des dissipations:            ! Calcul des tendances dynamiques:
           apphys = MOD(itau + 1, iphysiq) == 0 .AND. iflag_phys /= 0  
           conser = MOD(itau, iconser) == 0  
           apdiss = MOD(itau + 1, idissip) == 0  
   
           ! calcul des tendances dynamiques:  
126            CALL geopot(ip1jmp1, teta, pk, pks, phis, phi)            CALL geopot(ip1jmp1, teta, pk, pks, phis, phi)
127            CALL caldyn(itau, ucov, vcov, teta, ps, masse, pk, pkf, phis, phi, &            CALL caldyn(itau, ucov, vcov, teta, ps, masse, pk, pkf, phis, phi, &
128                 conser, du, dv, dteta, dp, w, pbaru, pbarv, &                 MOD(itau, iconser) == 0, du, dv, dteta, dp, w, pbaru, pbarv, &
129                 time + iday - day_ini)                 time + iday - day_ini)
130    
131            IF (forward .OR. leapf) THEN            IF (forward .OR. leapf) THEN
132               ! calcul des tendances advection des traceurs (dont l'humidite)               ! Calcul des tendances advection des traceurs (dont l'humidité)
133               CALL caladvtrac(q, pbaru, pbarv, p3d, masse, dq, teta, pk)               CALL caladvtrac(q, pbaru, pbarv, p3d, masse, dq, teta, pk)
134               IF (offline) THEN               IF (offline) THEN
135                  ! Stokage du flux de masse pour traceurs off-line                  ! Stokage du flux de masse pour traceurs off-line
# Line 146  contains Line 141  contains
141            ! integrations dynamique et traceurs:            ! integrations dynamique et traceurs:
142            CALL integrd(2, vcovm1, ucovm1, tetam1, psm1, massem1, dv, du, &            CALL integrd(2, vcovm1, ucovm1, tetam1, psm1, massem1, dv, du, &
143                 dteta, dq, dp, vcov, ucov, teta, q, ps, masse, phis, &                 dteta, dq, dp, vcov, ucov, teta, q, ps, masse, phis, &
144                 finvmaold, leapf)                 finvmaold, leapf, dt)
145    
146            IF (apphys) THEN            IF (MOD(itau + 1, iphysiq) == 0 .AND. iflag_phys /= 0) THEN
147               ! calcul des tendances physiques:               ! calcul des tendances physiques:
148               IF (itau + 1 == itaufin) lafin = .TRUE.               IF (itau + 1 == itaufin) lafin = .TRUE.
149    
# Line 158  contains Line 153  contains
153               rdaym_ini = itau * dtvr / daysec               rdaym_ini = itau * dtvr / daysec
154               rdayvrai = rdaym_ini + day_ini               rdayvrai = rdaym_ini + day_ini
155    
              ! Diagnostique de conservation de l'énergie : initialisation  
              IF (ip_ebil_dyn >= 1) THEN  
                 ztit='bil dyn'  
                 CALL diagedyn(ztit, 2, 1, 1, dtphys, ucov, vcov, ps, p3d, pk, &  
                      teta, q(:, :, 1), q(:, :, 2))  
              ENDIF  
   
156               CALL calfis(nqmx, lafin, rdayvrai, time, ucov, vcov, teta, q, &               CALL calfis(nqmx, lafin, rdayvrai, time, ucov, vcov, teta, q, &
157                    masse, ps, pk, phis, phi, du, dv, dteta, dq, w, &                    masse, ps, pk, phis, phi, du, dv, dteta, dq, w, &
158                    dufi, dvfi, dtetafi, dqfi, dpfi)                    dufi, dvfi, dtetafi, dqfi, dpfi)
# Line 173  contains Line 161  contains
161               CALL addfi(nqmx, dtphys, &               CALL addfi(nqmx, dtphys, &
162                    ucov, vcov, teta, q, ps, &                    ucov, vcov, teta, q, ps, &
163                    dufi, dvfi, dtetafi, dqfi, dpfi)                    dufi, dvfi, dtetafi, dqfi, dpfi)
   
              ! Diagnostique de conservation de l'énergie : difference  
              IF (ip_ebil_dyn >= 1) THEN  
                 ztit = 'bil phys'  
                 CALL diagedyn(ztit, 2, 1, 1, dtphys, ucov, vcov, ps, p3d, pk, &  
                      teta, q(:, :, 1), q(:, :, 2))  
              ENDIF  
164            ENDIF            ENDIF
165    
166            CALL pression(ip1jmp1, ap, bp, ps, p3d)            CALL pression(ip1jmp1, ap, bp, ps, p3d)
167            CALL exner_hyb(ps, p3d, pks, pk, pkf)            CALL exner_hyb(ps, p3d, pks, pk, pkf)
168    
169            IF (apdiss) THEN            IF (MOD(itau + 1, idissip) == 0) THEN
170               ! dissipation horizontale et verticale des petites echelles:               ! dissipation horizontale et verticale des petites echelles:
171    
172               ! calcul de l'energie cinetique avant dissipation               ! calcul de l'energie cinetique avant dissipation
# Line 197  contains Line 178  contains
178               ucov=ucov + dudis               ucov=ucov + dudis
179               vcov=vcov + dvdis               vcov=vcov + dvdis
180    
181               if (dissip_conservative) then               ! On rajoute la tendance due à la transformation Ec -> E
182                  ! On rajoute la tendance due a la transform. Ec -> E               ! thermique créée lors de la dissipation
183                  ! therm. cree lors de la dissipation               call covcont(llm, ucov, vcov, ucont, vcont)
184                  call covcont(llm, ucov, vcov, ucont, vcont)               call enercin(vcov, ucov, vcont, ucont, ecin)
185                  call enercin(vcov, ucov, vcont, ucont, ecin)               dtetaecdt= (ecin0 - ecin) / pk
186                  dtetaecdt= (ecin0 - ecin) / pk               dtetadis=dtetadis + dtetaecdt
                 dtetadis=dtetadis + dtetaecdt  
              endif  
187               teta=teta + dtetadis               teta=teta + dtetadis
188    
189               ! Calcul de la valeur moyenne, unique de h aux poles .....               ! Calcul de la valeur moyenne unique de h aux pôles
190               DO l = 1, llm               DO l = 1, llm
191                  DO ij = 1, iim                  DO ij = 1, iim
192                     tppn(ij) = aire(ij) * teta(ij, l)                     tppn(ij) = aire(ij) * teta(ij, l)
# Line 261  contains Line 240  contains
240            ENDIF            ENDIF
241    
242            IF (itau == itaufin) THEN            IF (itau == itaufin) THEN
243               CALL dynredem1("restart.nc", vcov, ucov, teta, q, masse, ps)               CALL dynredem1("restart.nc", vcov, ucov, teta, q, masse, ps, &
244                      itau=itau_dyn+itaufin)
245            ENDIF            ENDIF
246    
247            ! gestion de l'integration temporelle:            ! gestion de l'integration temporelle:
Legend:
Removed from v.26  
changed lines
  Added in v.28
  ViewVC Help
Powered by ViewVC 1.1.21