Index: trunk/LMDZ.COMMON/libf/dyn3d/caladvtrac.F
===================================================================
--- trunk/LMDZ.COMMON/libf/dyn3d/caladvtrac.F	(revision 2108)
+++ trunk/LMDZ.COMMON/libf/dyn3d/caladvtrac.F	(revision 2126)
@@ -9,7 +9,8 @@
 c
       USE infotrac, ONLY : nqtot
-      USE control_mod, ONLY : iapp_tracvl,planet_type
+      USE control_mod, ONLY : iapp_tracvl,planet_type,
+     &                        force_conserv_tracer
       USE comconst_mod, ONLY: dtvr
- 
+      USE planetary_operations, ONLY: planetary_tracer_amount_from_mass
       IMPLICIT NONE
 c
@@ -47,5 +48,14 @@
       INTEGER ij,l, iq, iapptrac
       REAL finmasse(ip1jmp1,llm), dtvrtrac
+      REAL :: totaltracer_old(nqtot),totaltracer_new(nqtot)
+      REAL :: ratio
 
+! Ehouarn : try to fix tracer conservation issues:
+      if (force_conserv_tracer) then
+        do iq=1,nqtot
+          call planetary_tracer_amount_from_mass(masse,q(:,:,iq),
+     &                                totaltracer_old(iq))
+        enddo
+      endif
 cc
 c
@@ -107,5 +117,16 @@
 c
         endif ! of if (planet_type.eq."earth")
-      ELSE
+        
+        ! Ehouarn : try to fix tracer conservation after tracer advection
+        if (force_conserv_tracer) then
+          do iq=1,nqtot
+            call planetary_tracer_amount_from_mass(masse,q(:,:,iq),
+     &                                  totaltracer_new(iq))
+            ratio=totaltracer_old(iq)/totaltracer_new(iq)
+            q(:,:,iq)=q(:,:,iq)*ratio
+          enddo
+        endif !of if (force_conserv_tracer)
+        
+      ELSE ! i.e. iapptrac.NE.iapp_tracvl
         if (planet_type.eq."earth") then
 ! Earth-specific treatment for the first 2 tracers (water)
Index: trunk/LMDZ.COMMON/libf/dyn3d/conf_gcm.F90
===================================================================
--- trunk/LMDZ.COMMON/libf/dyn3d/conf_gcm.F90	(revision 2108)
+++ trunk/LMDZ.COMMON/libf/dyn3d/conf_gcm.F90	(revision 2126)
@@ -17,5 +17,5 @@
                          output_grads_dyn, periodav, planet_type, &
                          raz_date, resetvarc, starttime, timestart, &
-                         ecritstart
+                         ecritstart,force_conserv_tracer
   USE infotrac, ONLY : type_trac
   use assert_m, only: assert
@@ -233,4 +233,8 @@
   iapp_tracvl = iperiod
   CALL getin('iapp_tracvl',iapp_tracvl)
+
+! Enforce tracer conservation in integrd & caladvtrac ?
+  force_conserv_tracer = .false. 
+  CALL getin('force_conserv_tracer',force_conserv_tracer)
 
 !Config  Key  = iconser
@@ -949,4 +953,6 @@
  write(lunout,*)' iphysiq = ', iphysiq
  write(lunout,*)' iflag_trac = ', iflag_trac
+ write(lunout,*)' iapp_tracvl = ', iapp_tracvl
+ write(lunout,*)' force_conserv_tracer = ', force_conserv_tracer
  write(lunout,*)' clon = ', clon
  write(lunout,*)' clat = ', clat
Index: trunk/LMDZ.COMMON/libf/dyn3d/integrd.F
===================================================================
--- trunk/LMDZ.COMMON/libf/dyn3d/integrd.F	(revision 2108)
+++ trunk/LMDZ.COMMON/libf/dyn3d/integrd.F	(revision 2126)
@@ -7,5 +7,5 @@
      &  )
 
-      use control_mod, only : planet_type
+      use control_mod, only : planet_type,force_conserv_tracer
       USE comvert_mod, ONLY: ap,bp
       USE comconst_mod, ONLY: pi
@@ -66,4 +66,5 @@
       REAL vscr( ip1jm ),uscr( ip1jmp1 ),hscr( ip1jmp1 ),pscr(ip1jmp1)
       REAL massescr( ip1jmp1,llm )
+      REAL :: massratio(ip1jmp1,llm)
 !      REAL finvmasse(ip1jmp1,llm)
       REAL p(ip1jmp1,llmp1)
@@ -242,4 +243,15 @@
 
       endif ! of if (planet_type.eq."earth")
+
+      if (force_conserv_tracer) then
+        ! Ehouarn: try to keep total amont of tracers fixed 
+        ! by acounting for mass change in each cell
+        massratio(1:ip1jmp1,1:llm)=massescr(1:ip1jmp1,1:llm)
+     &                             /masse(1:ip1jmp1,1:llm)
+        do iq=1,nq
+        q(1:ip1jmp1,1:llm,iq)=q(1:ip1jmp1,1:llm,iq)
+     &                        *massratio(1:ip1jmp1,1:llm)
+        enddo
+      endif ! of if (force_conserv_tracer)
 c
 c
Index: trunk/LMDZ.COMMON/libf/dyn3d/planetary_operations.F90
===================================================================
--- trunk/LMDZ.COMMON/libf/dyn3d/planetary_operations.F90	(revision 2126)
+++ trunk/LMDZ.COMMON/libf/dyn3d/planetary_operations.F90	(revision 2126)
@@ -0,0 +1,133 @@
+!
+! $Id: $
+!
+module planetary_operations
+! module which contains functions to obtain total values over the
+! entire globe (trivial in serial mode, but not in parallel)
+
+implicit none
+
+contains
+
+  subroutine planetary_atm_mass_from_ps(ps,value)
+  implicit none
+  include "dimensions.h"
+  include "paramet.h"
+  include "comgeom.h"
+  real,intent(in) :: ps(ip1jmp1)
+  real,intent(out) :: value
+  integer :: i,j,ij
+
+  ! compute total atmospheric mass over the whole planet
+  ! North Pole
+  value=ps(1)*airesurg(1)*iim
+  do j=2,jjm
+    do i=1,iim
+      ij=i+(j-1)*iip1
+      value=value+ps(ij)*airesurg(ij)
+    enddo
+  enddo
+  ! South pole
+  value=value+ps(ip1jmp1-iim)*airesurg(ip1jmp1-iim)*iim
+  
+  end subroutine planetary_atm_mass_from_ps
+
+  subroutine planetary_atm_mass_from_mass(mass,value)
+  implicit none
+  include "dimensions.h"
+  include "paramet.h"
+  real,intent(in) :: mass(ip1jmp1,llm) ! air mass (kg)
+  real,intent(out) :: value
+  integer :: i,j,ij,l
+  real :: column(ip1jmp1) ! columns amount of air (kg)
+
+  ! compute total atmospheric mass over the whole planet
+  ! 1. build column amout of tracer (kg)
+  column(:)=0
+  do ij=1,ip1jmp1
+    do l=1,llm
+      column(ij)=column(ij)+mass(ij,l)
+    enddo
+  enddo
+  
+  ! 2. Compute total tracer mass over the planet
+  !North pole
+  value=column(1)*iim
+  do j=2,jjm
+    do i=1,iim
+      ij=i+(j-1)*iip1
+      value=value+column(ij)
+    enddo
+  enddo
+  ! South pole
+  value=value+column(ip1jmp1-iim)*iim
+  
+  end subroutine planetary_atm_mass_from_mass
+  
+  subroutine planetary_tracer_amount_from_p(p,q,amount)
+  use comconst_mod, ONLY: g
+  implicit none
+  include "dimensions.h"
+  include "paramet.h"
+  include "comgeom.h"
+  real,intent(in) :: p(ip1jmp1,llmp1) ! pressure at vertical mesh interfaces
+  real,intent(in) :: q(ip1jmp1,llm) ! 3D tracer (kg/kg_air)
+  real,intent(out) :: amount
+  integer :: i,j,ij,l
+  real :: column(ip1jmp1) ! columns amount of tracer (kg.m-2)
+  
+  ! 1. build column amout of tracer (kg.m-2)
+  column(:)=0
+  do ij=1,ip1jmp1
+    do l=1,llm
+      column(ij)=column(ij)+q(ij,l)*(p(ij,l)-p(ij,l+1))/g
+    enddo
+  enddo
+  
+  ! 2. Compute total tracer mass over the planet
+  !North pole
+  amount=column(1)*aire(1)*iim
+  do j=2,jjm
+    do i=1,iim
+      ij=i+(j-1)*iip1
+      amount=amount+column(ij)*aire(ij)
+    enddo
+  enddo
+  ! South pole
+  amount=amount+column(ip1jmp1-iim)*aire(ip1jmp1-iim)*iim
+
+  end subroutine planetary_tracer_amount_from_p
+
+  subroutine planetary_tracer_amount_from_mass(mass,q,amount)
+  implicit none
+  include "dimensions.h"
+  include "paramet.h"
+  real,intent(in) :: mass(ip1jmp1,llm) ! air mass (kg)
+  real,intent(in) :: q(ip1jmp1,llm) ! 3D tracer (kg/kg_air)
+  real,intent(out) :: amount
+  integer :: i,j,ij,l
+  real :: column(ip1jmp1) ! columns amount of tracer (kg)
+  
+  ! 1. build column amout of tracer (kg)
+  column(:)=0
+  do ij=1,ip1jmp1
+    do l=1,llm
+      column(ij)=column(ij)+q(ij,l)*mass(ij,l)
+    enddo
+  enddo
+  
+  ! 2. Compute total tracer mass over the planet
+  !North pole
+  amount=column(1)*iim
+  do j=2,jjm
+    do i=1,iim
+      ij=i+(j-1)*iip1
+      amount=amount+column(ij)
+    enddo
+  enddo
+  ! South pole
+  amount=amount+column(ip1jmp1-iim)*iim
+
+  end subroutine planetary_tracer_amount_from_mass
+
+end module planetary_operations