Index: LMDZ6/trunk/libf/phylmd/thermcell_down.F90
===================================================================
--- LMDZ6/trunk/libf/phylmd/thermcell_down.F90	(revision 4363)
+++ LMDZ6/trunk/libf/phylmd/thermcell_down.F90	(revision 4365)
@@ -1,8 +1,65 @@
-   SUBROUTINE thermcell_down(ngrid,nlay,po,pt,pu,pv,pplay,  &
-     &           pplev,deltaz,lmax,fup,eup,dup)
+   SUBROUTINE thermcell_updown_dq(ngrid,nlay,lmax,eup,dup,edn,ddn,masse,theta,dtheta,thetau,thetad)
 
 !--------------------------------------------------------------
-!thermcell_env: calcule les caracteristiques de l environnement
-!necessaires au calcul des proprietes dans le thermique
+! thermcell_updown_dq: calcul du transport d'un traceur en présence 
+! d'up/down drafts
+!--------------------------------------------------------------
+
+   ! Suite du travail :
+   ! Calculer les tendances d'un traceur (ici theta) en tenant compte
+   ! des up et down drafts et de la subsidence compensatoire.
+
+
+   IMPLICIT NONE
+
+! arguments
+
+   integer,intent(in) :: ngrid,nlay
+   real,intent(in), dimension(ngrid,nlay) :: eup,dup,edn,ddn,masse
+   real,intent(in), dimension(ngrid,nlay) :: theta
+   real,intent(out), dimension(ngrid,nlay) :: thetau,thetad,dtheta
+   integer, intent(out), dimension(ngrid) :: lmax
+
+   
+! Local
+
+   real, dimension(ngrid,nlay+1) :: fup,fdn
+
+   integer ig,ilay
+
+   fdn(:,:)=0.
+   thetad(:,:)=0.
+
+   ! lmax : indice tel que fu(kmax+1)=0
+   
+   ! Dans ce cas, pas besoin d'initialiser thetad(lmax) ( =theta(lmax) )
+
+   do ilay=nlay,1,-1
+      do ig=1,ngrid
+         if (ilay.le.lmax(ig)) then
+            fdn(ig,ilay)=fdn(ig,ilay+1)+edn(ig,ilay)-ddn(ig,ilay)
+            thetad(ig,ilay)=( fdn(ig,ilay+1)*thetad(ig,ilay+1) + edn(ig,ilay)*theta(ig,ilay) ) / (fdn(ig,ilay)+ddn(ig,ilay))
+         endif
+      enddo 
+   enddo
+
+
+
+!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+! Initialisations :
+!------------------
+
+
+!
+ RETURN
+   END
+
+!=========================================================================
+
+   SUBROUTINE thermcell_down(ngrid,nlay,po,pt,pu,pv,pplay,pplev,  &
+     &           lmax,fup,eup,dup,theta)
+
+!--------------------------------------------------------------
+!thermcell_down: calcul des propriétés du panache descendant.
 !--------------------------------------------------------------
 
@@ -14,5 +71,6 @@
 
    integer,intent(in) :: ngrid,nlay
-   real,intent(in), dimension(ngrid,nlay) :: po,pt,pu,pv,pplay,eup,dup,deltaz
+   real,intent(in), dimension(ngrid,nlay) :: po,pt,pu,pv,pplay,eup,dup
+   real,intent(in), dimension(ngrid,nlay) :: theta
    real,intent(in), dimension(ngrid,nlay+1) :: pplev,fup
    integer, intent(in), dimension(ngrid) :: lmax
@@ -22,5 +80,5 @@
 ! Local
 
-   real, dimension(ngrid,nlay) :: edn,ddn
+   real, dimension(ngrid,nlay) :: edn,ddn,thetad
    real, dimension(ngrid,nlay+1) :: fdn
 
@@ -32,17 +90,31 @@
    ddn(:,:)=0.
    fdn(:,:)=0.
+   thetad(:,:)=0.
 
-   ! lmax : indice de la derniere couche ou les thermiques sont actifs
+   ! lmax : indice tel que fu(kmax+1)=0
    
-   ! Convention de flux : fdn négatif vers le bas
-   ! Pas évident qu'on veuille conserver ce choix
+   ! Dans ce cas, pas besoin d'initialiser thetad(lmax) ( =theta(lmax) )
 
-   print*,'LMAX ',lmax(1)
-   do ilay=nlay+1,1,-1
+! FH MODIFS APRES REUNIONS POUR COMMISSIONS
+! quelques erreurs de declaration
+! probleme si lmax=1 ce qui a l'air d'être le cas en début de simu. Devrait être 0 ?
+! Remarques :
+! on pourrait écrire la formule de thetad
+!    www=fdn(ig,ilay+1)/ (fdn(ig,ilay)+ddn(ig,ilay))
+!    thetad(ig,ilay)= www * thetad(ig,ilay+1) + (1.-www) * theta(ig,ilay) 
+! Elle a l'avantage de bien montré la conservation, l'idée fondamentale dans le 
+!   transport qu'on ne fait que sommer des "sources" au travers d'un "propagateur"
+!   (Green)
+! Elle montre aussi beaucoup plus clairement pourquoi on n'a pas à se souccier (trop)
+!   de la possible nulité du dénominateur
+
+
+   do ilay=nlay,1,-1
       do ig=1,ngrid
-         if (ilay.le.lmax(ig)) then
-            edn(ig,ilay)=0.5*dup(ig,ilay)*deltaz(ig,ilay)
-            ddn(ig,ilay)=0.
-            fdn(ig,ilay)=-(-fdn(ig,ilay+1)+edn(ig,ilay)-ddn(ig,ilay))
+         if (ilay.le.lmax(ig).and.lmax(ig)>1) then
+            edn(ig,ilay)=0.5*dup(ig,ilay)
+            ddn(ig,ilay)=0.5*eup(ig,ilay)
+            fdn(ig,ilay)=fdn(ig,ilay+1)+edn(ig,ilay)-ddn(ig,ilay)
+            thetad(ig,ilay)=( fdn(ig,ilay+1)*thetad(ig,ilay+1) + edn(ig,ilay)*theta(ig,ilay) ) / (fdn(ig,ilay)+ddn(ig,ilay))
          endif
       enddo 
Index: LMDZ6/trunk/libf/phylmd/thermcell_main.F90
===================================================================
--- LMDZ6/trunk/libf/phylmd/thermcell_main.F90	(revision 4363)
+++ LMDZ6/trunk/libf/phylmd/thermcell_main.F90	(revision 4365)
@@ -459,6 +459,7 @@
 #ifdef DevThermcellDown
       print*,'WARNING !!! routine thermcell_down en cours de developpement'
-      CALL thermcell_down(ngrid,nlay,po,pt,pu,pv,pplay,  &
-     &           pplev,deltaz,lmax,fm,entr,detr)
+      CALL thermcell_down(ngrid,nlay,po,pt,pu,pv,pplay,pplev,  &
+     &           lmax,fm,entr,detr,zthl)
+
 #endif