Index: trunk/libf/dyn3dpar/addfi_p.F
===================================================================
--- trunk/libf/dyn3dpar/addfi_p.F	(revision 101)
+++ trunk/libf/dyn3dpar/addfi_p.F	(revision 108)
@@ -28,5 +28,5 @@
 c      pdufi(ip1jmp1,llm)     |
 c      pdvfi(ip1jm,llm)       |   respective
-c      pdhfi(ip1jmp1)         |      tendencies
+c      pdhfi(ip1jmp1)         |      tendencies  (unit/s)
 c      pdtsfi(ip1jmp1)        |
 c
Index: trunk/libf/dyn3dpar/calfis_p.F
===================================================================
--- trunk/libf/dyn3dpar/calfis_p.F	(revision 101)
+++ trunk/libf/dyn3dpar/calfis_p.F	(revision 108)
@@ -952,9 +952,4 @@
       
 ! ADAPTATION GCM POUR CP(T)
-!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
-      DO l=1,llm
-        ztfi(1:klon,l)=ztfi(1:klon,l)+zdtfi(1:klon,l)*dtphys
-      ENDDO
-!$OMP END DO
       call t2tpot_p(ngridmx,llm,ztfi,zteta,zpk)
 
Index: trunk/libf/dyn3dpar/comconst.h
===================================================================
--- trunk/libf/dyn3dpar/comconst.h	(revision 101)
+++ trunk/libf/dyn3dpar/comconst.h	(revision 108)
@@ -6,9 +6,9 @@
 
       COMMON/comconsti/im,jm,lllm,imp1,jmp1,lllmm1,lllmp1,lcl,          &
-     &                 iflag_top_bound
+     &                 iflag_top_bound,mode_top_bound
       COMMON/comconstr/dtvr,daysec,                                     &
      & pi,dtphys,dtdiss,rad,r,cpp,kappa,cotot,unsim,g,omeg              &
-     &                   ,dissip_factz,dissip_deltaz,dissip_zref        &
-     &                   ,tau_top_bound,                                &
+     & ,dissip_fac_mid,dissip_fac_up,dissip_deltaz,dissip_hdelta        &
+     & ,dissip_pupstart  ,tau_top_bound,                                &
      & daylen,year_day,molmass, ihf
       COMMON/cpdetvenus/nu_venus,t0_venus
@@ -28,6 +28,7 @@
       REAL g ! (m/s2) gravity
       REAL omeg ! (rad/s) rotation rate of the planet
-      REAL dissip_factz,dissip_deltaz,dissip_zref
-      INTEGER iflag_top_bound
+      REAL dissip_fac_mid,dissip_fac_up,dissip_deltaz,dissip_hdelta
+      REAL dissip_pupstart
+      INTEGER iflag_top_bound,mode_top_bound
       REAL tau_top_bound
       REAL daylen ! length of solar day, in 'standard' day length
Index: trunk/libf/dyn3dpar/conf_gcm.F
===================================================================
--- trunk/libf/dyn3dpar/conf_gcm.F	(revision 101)
+++ trunk/libf/dyn3dpar/conf_gcm.F	(revision 108)
@@ -314,17 +314,28 @@
 ! Parametres controlant la variation sur la verticale des constantes de
 ! dissipation.
-! Pour le moment actifs uniquement dans la version a 39 niveaux
-! avec ok_strato=y
-
-       dissip_factz=4.
-       dissip_deltaz=10.
-       dissip_zref=30.
-       CALL getin('dissip_factz',dissip_factz )
+! Actifs uniquement avec ok_strato=y
+
+       dissip_fac_mid=2.
+       dissip_fac_up=10.
+       dissip_deltaz=10.! Intervalle (km) pour le changement mid / up
+       dissip_hdelta=5. ! scale height (km) dans la zone de la transition(m)
+       dissip_pupstart=1.e3  ! pression (Pa) au bas la transition mid / up
+       CALL getin('dissip_fac_mid',dissip_fac_mid )
+       CALL getin('dissip_fac_up',dissip_fac_up )
        CALL getin('dissip_deltaz',dissip_deltaz )
-       CALL getin('dissip_zref',dissip_zref )
-
+       CALL getin('dissip_hdelta',dissip_hdelta )
+       CALL getin('dissip_pupstart',dissip_pupstart )
+
+! Parametres controlant la couche eponge
+! Actifs uniquement avec ok_strato=y
+!   mode = 0 : pas de sponge
+!   mode = 1 : u et v -> 0
+!   mode = 2 : u et v -> moyenne zonale
+!   mode = 3 : u, v et h -> moyenne zonale
        iflag_top_bound=1
+       mode_top_bound=3
        tau_top_bound=1.e-5
        CALL getin('iflag_top_bound',iflag_top_bound)
+       CALL getin('mode_top_bound',mode_top_bound)
        CALL getin('tau_top_bound',tau_top_bound)
 
Index: trunk/libf/dyn3dpar/leapfrog_p.F
===================================================================
--- trunk/libf/dyn3dpar/leapfrog_p.F	(revision 101)
+++ trunk/libf/dyn3dpar/leapfrog_p.F	(revision 108)
@@ -113,7 +113,9 @@
       REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: dqfi
 
-       !! Aymeric -- cp(T) comme dans leapfrog.F, SAVE OK ???
-      REAL,SAVE :: duspg(ip1jmp1,llm) ! for bilan_dyn
-
+c   tendances top_bound (sponge layer)
+      REAL,SAVE :: dvtop(ip1jm,llm),dutop(ip1jmp1,llm)
+      REAL,SAVE :: dtetatop(ip1jmp1,llm)
+      REAL,SAVE :: dptop(ip1jmp1)
+      REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: dqtop
 
 c   variables pour le fichier histoire
@@ -242,4 +244,5 @@
             ALLOCATE(dqfi(ip1jmp1,llm,nqtot))
             ALLOCATE(dqfi_tmp(iip1,llm,nqtot))
+            ALLOCATE(dqtop(ip1jmp1,llm,nqtot))
          END IF
 c$OMP END MASTER      
@@ -251,5 +254,20 @@
 
 c$OMP MASTER
+c INITIALISATIONS
+        dudis(:,:)   =0.
+        dvdis(:,:)   =0.
+        dtetadis(:,:)=0.
+        dutop(:,:)   =0.
+        dvtop(:,:)   =0.
+        dtetatop(:,:)=0.
+        dqtop(:,:,:) =0.
+        dptop(:)     =0.
+        dufi(:,:)   =0.
+        dvfi(:,:)   =0.
+        dtetafi(:,:)=0.
+        dqfi(:,:,:) =0.
+        dpfi(:)     =0.
       dq(:,:,:)=0.
+
       CALL pression ( ip1jmp1, ap, bp, ps, p       )
       CALL exner_hyb( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf )
@@ -914,12 +932,18 @@
 c      ajout des tendances physiques:
 c      ------------------------------
-         IF (ok_strato) THEN
-           CALL top_bound_p( vcov,ucov,teta,masse,dufi,dvfi,dtetafi)
-         ENDIF
-       
           CALL addfi_p( dtphys, leapf, forward   ,
      $                  ucov, vcov, teta , q   ,ps ,
      $                 dufi, dvfi, dtetafi , dqfi ,dpfi  )
 
+c      Couche superieure :
+c      -------------------
+         IF (ok_strato) THEN
+           CALL top_bound_p( vcov,ucov,teta,masse,dutop,dvtop,dtetatop)
+           CALL addfi_p( dtphys, leapf, forward   ,
+     $                  ucov, vcov, teta , q   ,ps ,
+     $                 dutop, dvtop, dtetatop , dqtop ,dptop  )
+
+         ENDIF
+       
 c$OMP BARRIER
 c$OMP MASTER
@@ -1444,5 +1468,5 @@
                  CALL bilan_dyn_p(dtvr*iperiod,dtvr*day_step*periodav,
      &                 ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,
-     &                 du,dudis,duspg,dufi)
+     &                 du,dudis,dutop,dufi)
 c$OMP END MASTER
               ENDIF !ok_dynzon
@@ -1670,5 +1694,5 @@
                  CALL bilan_dyn_p(dtvr*iperiod,dtvr*day_step*periodav,
      &                 ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,
-     &                 du,dudis,duspg,dufi)
+     &                 du,dudis,dutop,dufi)
 
 c$OMP END MASTER
Index: trunk/libf/dyn3dpar/top_bound_p.F
===================================================================
--- trunk/libf/dyn3dpar/top_bound_p.F	(revision 101)
+++ trunk/libf/dyn3dpar/top_bound_p.F	(revision 108)
@@ -64,5 +64,5 @@
              rdamp(llm-3)=tau_top_bound/8.
          else if (iflag_top_bound == 2) then
-! couce eponge dans toutes les couches de pression plus faible que
+! couche eponge dans toutes les couches de pression plus faible que
 ! 100 fois la pression de la derniere couche
              rdamp(:)=tau_top_bound
@@ -70,5 +70,15 @@
          endif
          first=.false.
-         print*,'TOP_BOUND rdamp=',rdamp
+         print*,'TOP_BOUND mode',mode_top_bound
+         print*,'Coeffs pour la couche eponge a l equateur'
+         print*,'p (Pa)  z(km)  tau (s)   dt*rdamp'
+         do l=1,llm
+           if (rdamp(l).ne.0.) then
+            zkm        = phi(iip1/2,jjp1/2,l)/(1000*g)
+          print*,presnivs(l),zkm,
+     .          1./rdamp(l),
+     .          dt*rdamp(l)
+           endif
+         enddo
 c$OMP END MASTER
 c$OMP BARRIER
@@ -77,4 +87,12 @@
 
       CALL massbar_p(masse,massebx,masseby)
+
+c   mode = 0 : pas de sponge
+c   mode = 1 : u et v -> 0
+c   mode = 2 : u et v -> moyenne zonale
+c   mode = 3 : u, v et h -> moyenne zonale
+
+C POUR V
+
 C  CALCUL DES CHAMPS EN MOYENNE ZONALE:
 
@@ -84,5 +102,6 @@
 
 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)      
-      do l=1,llm
+      if (mode_top_bound.ge.2) then
+       do l=1,llm
         do j=jjb,jje
           zm=0.
@@ -97,16 +116,31 @@
           vzon(j,l)=vzon(j,l)/zm
         enddo
-      enddo
+       enddo
+      else
+       do l=1,llm
+        do j=jjb,jje
+          vzon(j,l)=0.
+        enddo
+       enddo
+      endif
 c$OMP END DO NOWAIT   
 
-c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)      
-      do l=1,llm
+C   AMORTISSEMENTS LINEAIRES:
+
+c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)      
+      if (mode_top_bound.ge.1) then
+       do l=1,llm
         do j=jjb,jje
           do i=1,iip1
-            dv(i,j,l)=dv(i,j,l)-rdamp(l)*(vcov(i,j,l)-vzon(j,l))
-          enddo
-        enddo
-      enddo
-c$OMP END DO NOWAIT
+            dv(i,j,l)= -rdamp(l)*(vcov(i,j,l)-vzon(j,l))
+          enddo
+        enddo
+       enddo
+      endif
+c$OMP END DO NOWAIT
+
+C POUR U ET H
+
+C  CALCUL DES CHAMPS EN MOYENNE ZONALE:
 
       jjb=jj_begin
@@ -116,5 +150,6 @@
 
 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)      
-      do l=1,llm
+      if (mode_top_bound.ge.2) then
+       do l=1,llm
         do j=jjb,jje
           uzon(j,l)=0.
@@ -126,9 +161,17 @@
           uzon(j,l)=uzon(j,l)/zm
         enddo
-      enddo
+       enddo
+      else
+       do l=1,llm
+        do j=jjb,jje
+          uzon(j,l)=0.
+        enddo
+       enddo
+      endif
 c$OMP END DO NOWAIT
 
 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)   
-      do l=1,llm
+      if (mode_top_bound.ge.3) then
+       do l=1,llm
         do j=jjb,jje
           zm=0.
@@ -140,5 +183,6 @@
           tzon(j,l)=tzon(j,l)/zm
         enddo
-      enddo
+       enddo
+      endif
 c$OMP END DO NOWAIT
 
@@ -146,13 +190,25 @@
 
 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)      
-      do l=1,llm
+      if (mode_top_bound.ge.1) then
+       do l=1,llm
         do j=jjb,jje
           do i=1,iip1
-            du(i,j,l)=du(i,j,l)
-     s               -rdamp(l)*(ucov(i,j,l)-cu(i,j)*uzon(j,l))
-            dh(i,j,l)=dh(i,j,l)-rdamp(l)*(teta(i,j,l)-tzon(j,l))
-          enddo
-       enddo
-      enddo
+            du(i,j,l)= -rdamp(l)*(ucov(i,j,l)-cu(i,j)*uzon(j,l))
+          enddo
+        enddo
+       enddo
+      endif
+c$OMP END DO NOWAIT
+      
+c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)      
+      if (mode_top_bound.ge.3) then
+       do l=1,llm
+        do j=jjb,jje
+          do i=1,iip1
+            dh(i,j,l)= -rdamp(l)*(teta(i,j,l)-tzon(j,l))
+          enddo
+        enddo
+       enddo
+      endif
 c$OMP END DO NOWAIT
       
