Index: trunk/LMDZ.VENUS/libf/phyvenus/physiq.F
===================================================================
--- trunk/LMDZ.VENUS/libf/phyvenus/physiq.F	(revision 1310)
+++ trunk/LMDZ.VENUS/libf/phyvenus/physiq.F	(revision 1442)
@@ -64,5 +64,5 @@
       USE write_field_phy
       USE iophy
-      use cpdet_mod, only: cpdet, t2tpot
+      USE cpdet_mod, only: cpdet, t2tpot
       USE chemparam_mod
       USE conc
@@ -240,4 +240,6 @@
       EXTERNAL conduction
       EXTERNAL molvis
+      EXTERNAL moldiff_red
+
 c
 c Variables locales
@@ -260,5 +262,5 @@
       REAL zdtime, zlongi
 c
-      INTEGER i, k, iq, ig, j, ll
+      INTEGER i, k, iq, ig, j, ll, ilon, ilat, ilev
 c
       REAL zphi(klon,klev)
@@ -301,4 +303,7 @@
       real d_v_molvis(klon,klev)     ! (m/s) /s
 
+c Tendencies due to molecular diffusion
+      real d_q_moldif(klon,klev,nqmax)
+
 c
 c Variables liees a l'ecriture de la bande histoire physique
@@ -322,4 +327,10 @@
       REAL :: tr_seri(klon,klev,nqmax)
       REAL :: d_tr(klon,klev,nqmax)
+
+c Champ de modification de la temperature par rapport a VIRAII
+      REAL delta_temp(klon,klev)
+c      SAVE delta_temp
+      REAL mat_dtemp(33,50)
+      SAVE mat_dtemp
 
 c Variables tendance sedimentation
@@ -493,5 +504,5 @@
 C TRACEURS
 C source dans couche limite
-         source = 0.0 ! pas de source, pour l'instant
+         source(:,:) = 0.0 ! pas de source, pour l'instant
 c---------
 
@@ -510,7 +521,5 @@
             rnew(ig,j)=R
             cpnew(ig,j)=cpdet(tmoy(j))
-c            print*, ' physique  l503'
-c            print*,  j, cpdet(tmoy(j))
-             mmean(ig,j)=RMD
+            mmean(ig,j)=RMD
             akknew(ig,j)=1.e-4
           enddo
@@ -522,6 +531,6 @@
          call compo_hedin83_init2
       ENDIF
-      if (callnlte) call nlte_setup
-      if(callnirco2.and.(nircorr.eq.1)) call nir_leedat         
+      if (callnlte.and.nltemodel.eq.2) call nlte_setup
+      if (callnirco2.and.nircorr.eq.1) call nir_leedat         
 c---------
       
@@ -627,5 +636,5 @@
       endif
            
-      if ((nlon .GT. 1) .AND. ok_chem) then
+      if ((nlon .GT. 1) .AND. (ok_chem.OR.ok_cloud)) then
 c !!! DONC 3D !!!
         CALL chemparam_ini() 
@@ -636,7 +645,37 @@
         CALL cloud_ini(nlon,nlev)
       endif
+
+c======================================================================
+c      Lecture du fichier DeltaT
+c======================================================================
+
+c     ATTENTION tout ce qui suit est pour un 48*32*50
+
+      if (ok_deltatemp) then
+
+      print*,'lecture de VenusDeltaT.txt '
+      open(99, form = 'formatted', status = 'old', file = 
+     &     'VenusDeltaT.dat')
+      print*,'Ouverture de VenusDeltaT.txt '
+
+      DO ilev = 1, klev
+         read(99,'(33(1x,e13.6))') (mat_dtemp(ilat,ilev),ilat=1,33)
+         print*,'lecture de VenusDeltaT.txt ligne:',ilev
+      ENDDO
+      
+      close(99)
+      print*,'FIN lecture de VenusDeltaT.txt ok.'
+
+      DO k = 1, klev
+      DO i = 1, klon      
+         ilat=(rlatd(i)/5.625) + 17.
+         delta_temp(i,k)=mat_dtemp(INT(ilat),k)
+      ENDDO
+      ENDDO
+      
+      endif
         
       ENDIF ! debut
-c====================================================================
+c======================================================================
 c======================================================================
 
@@ -830,18 +869,45 @@
 ! Case 3: Full chemistry and/or clouds
 
-         call phytrac_chimie(                  
+         if (ok_deltatemp) then
+!           PRINT*,'Def de delta_temp'
+           DO k = 1, klev
+           DO i = 1, klon      
+             ilat=(rlatd(i)/5.625) + 17.
+!             PRINT*,INT(ilat),rlatd(i),mat_dtemp(INT(ilat),k)
+             delta_temp(i,k)=mat_dtemp(INT(ilat),k)
+           ENDDO
+           ENDDO
+      
+         endif
+
+         if (ok_deltatemp) then
+! Utilisation du champ de temperature modifie         
+           call phytrac_chimie(                              
      I             debut,
      I             gmtime,
      I             nqmax,
-     I             nlon,
+     I             klon,
      I             rlatd,
      I             rlond,
      I             nlev,
      I             dtime,
-     I             t_seri,pplay,
-     O             tr_seri,
-     O             NBRTOT,
-     O             WH2SO4,
-     O             rho_droplet)
+     I             t_seri+delta_temp,
+     I             pplay,
+     O             tr_seri)
+         else
+         
+           call phytrac_chimie(                              
+     I             debut,
+     I             gmtime,
+     I             nqmax,
+     I             klon,
+     I             rlatd,
+     I             rlond,
+     I             nlev,
+     I             dtime,
+     I             t_seri,
+     I             pplay,
+     O             tr_seri)
+         endif
 
 c        CALL WriteField_phy('Pression',pplay,nlev)
@@ -856,28 +922,37 @@
          if (ok_sedim) then 
 
+         if (ok_deltatemp) then
+! Utilisation du champ de temperature modifie  
            CALL new_cloud_sedim(
-     I               klon,
-     I               nlev,
-     I               dtime,
-     I               pplay,
-     I               paprs,
-     I               t_seri, 
-     I               WH2SO4,
-     I               tr_seri,
-     I               nqmax,
-     I               NBRTOT,
-     I               rho_droplet,
-     O               Fsedim,
-     O               d_tr_sed,
-     O               d_tr_ssed)
-
-          DO k = 1, klev
-           DO i = 1, klon
-     
-c        WRITE(88,"(11(e15.8,','))") pplay(5,25),
-c     &  t_seri(5,25),tr_seri(5,25,i_h2oliq),
-c     &  tr_seri(5,25,i_h2o),tr_seri(5,25,i_h2so4liq),
-c     &  tr_seri(5,25,i_h2so4),NBRTOT(5,25),WH2SO4(5,25),
-c     &  Fsedim(5,25),d_tr_sed(5,25,1),d_tr_sed(5,25,2)
+     I                 klon,
+     I		     nlev,
+     I		     dtime,
+     I                 pplay,
+     I		     paprs,
+     I		     t_seri+delta_temp,
+     I                 tr_seri,
+     O		     d_tr_sed,
+     O		     d_tr_ssed,
+     I		     nqmax,
+     O                 Fsedim)
+         else
+         
+           CALL new_cloud_sedim(
+     I                 klon,
+     I		     nlev,
+     I		     dtime,
+     I                 pplay,
+     I		     paprs,
+     I		     t_seri,
+     I                 tr_seri,
+     O		     d_tr_sed,
+     O		     d_tr_ssed,
+     I		     nqmax,
+     O                 Fsedim)
+
+         endif
+
+        DO k = 1, klev
+         DO i = 1, klon
 
 c--------------------
@@ -888,5 +963,5 @@
         PRINT*,'d_tr_sed Nan?',d_tr_sed(i,k,:),'Temp',t_seri(i,k)
         PRINT*,'lat-lon',i,'level',k,'dtime',dtime
-        PRINT*,'NBRTOT',NBRTOT(i,k),'F_sed',Fsedim(i,k)
+        PRINT*,'F_sed',Fsedim(i,k)
         PRINT*,'==============================================='
                 d_tr_sed(i,k,:)=0.
@@ -901,6 +976,6 @@
         Fsedim(i,k)     = Fsedim(i,k) / dtime
      
-           ENDDO
           ENDDO
+         ENDDO
       
         Fsedim(:,klev+1) = 0.
@@ -988,13 +1063,38 @@
          ENDDO
          ENDDO
+   
          DO iq=1, nqmax
+c  AS: changement
+c  Pourquoi d_tr_vdf(1,1,iq) et tr_seri(1,1,iq)
+c  et pas d_tr_vdf(:,:,iq) tr_seri(:,:,iq)
+c  Je vois pas en quoi cltrac ne prendrait en compte que le traceur à la surface et au point 1 en klon
+c
+
+c  Je garde le source(:,iq) parce que je comprend pas sinon source 
+c   dimension(klon,nqmax) et flux dans cltrac (klon) ???
+
+c             CALL cltrac(dtime,ycoefh,t_seri,
+c     s               tr_seri(1,1,iq),source(:,iq),
+c     e               paprs, pplay,delp,
+c     s               d_tr_vdf(1,1,iq))
+     
              CALL cltrac(dtime,ycoefh,t_seri,
-     s               tr_seri(1,1,iq),source,
+     s               tr_seri(:,:,iq),source(:,iq),
      e               paprs, pplay,delp,
-     s               d_tr_vdf(1,1,iq))
+     s               d_tr_vdf(:,:,iq))
+     
              tr_seri(:,:,iq) = tr_seri(:,:,iq) + d_tr_vdf(:,:,iq)
              d_tr_vdf(:,:,iq)= d_tr_vdf(:,:,iq)/dtime          ! /s
+
+         DO k = 1, klev
+         DO i = 1, klon
+            tr_seri(i,k,iq) = tr_seri(i,k,iq) + d_tr_vdf(i,k,iq)
+            d_tr_vdf(i,k,iq)= d_tr_vdf(i,k,iq)/dtime          ! /s
          ENDDO
-      endif
+         ENDDO
+         
+         ENDDO !nqmax
+
+       endif
 
       IF (if_ebil.ge.2) THEN 
@@ -1084,5 +1184,4 @@
            d_tr_ajs(:,:,:)= d_tr_ajs(:,:,:)/dtime  ! /s
          endif
-
       endif
 
@@ -1104,9 +1203,10 @@
       END IF 
 
+
 c====================================================================
 c RAYONNEMENT
 c====================================================================
 
-c-----------------------------------------------------------------------
+c------------------------------------
 c    . Compute radiative tendencies :
 c------------------------------------
@@ -1118,26 +1218,51 @@
 c     PRINT*,'dtimerad,dtime,radpas',dtimerad,dtime,radpas
             
-c Calcul pour Cp rnew et mmean avec traceurs (Cp independant de T !! )
-       IF(callnlte.or.callthermos) THEN                                 
+
+c------------------------------------
+c    . Compute mean mass, cp and R :
+c------------------------------------
+
+      if(callthermos) then
+         call concentrations2(pplay,t_seri,d_t,tr_seri, nqmax,
+     &                        pdtphys)
+
+      endif
+
+
+cc!!! ADD key callhedin 
+
+      IF(callnlte.or.callthermos) THEN                                 
          call compo_hedin83_mod(pplay,rmu0,   
      &	         co2vmr_gcm,covmr_gcm,ovmr_gcm,n2vmr_gcm,nvmr_gcm)
+
+         IF(ok_chem) then
+ 
+CC  !! GG : Using only mayor species tracers abundances to compute NLTE heating/cooling
+
+CC               Conversion [mmr] ---> [vmr]
+        
+                 co2vmr_gcm(:,:) = tr_seri(1:nlon,1:nlev,i_co2)*
+     &                             mmean(1:nlon,1:nlev)/M_tr(i_co2)
+                 covmr_gcm(:,:)  = tr_seri(1:nlon,1:nlev,i_co)* 
+     &                              mmean(1:nlon,1:nlev)/M_tr(i_co)
+                 ovmr_gcm(:,:)   = tr_seri(1:nlon,1:nlev,i_o)*
+     &                             mmean(1:nlon,1:nlev)/M_tr(i_o)
+                 n2vmr_gcm(:,:)   = tr_seri(1:nlon,1:nlev,i_n2)*
+     &                             mmean(1:nlon,1:nlev)/M_tr(i_n2)
+
+         ENDIF
+
        ENDIF	
 
-      if(callthermos) then
-         call concentrations2(pplay,t_seri,d_t,co2vmr_gcm, n2vmr_gcm,
-     &		covmr_gcm, ovmr_gcm,nvmr_gcm,pdtphys)
-      endif
-
-
-c          NLTE cooling from CO2 emission
-c          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+c
+c   NLTE cooling from CO2 emission
+c   ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
         IF(callnlte) THEN                                 
             if(nltemodel.eq.0.or.nltemodel.eq.1) then
-                CALL nltecool(klon, klev, pplay*9.869e-6, t_seri,
-     $                    co2vmr_gcm,n2vmr_gcm, covmr_gcm, ovmr_gcm, 
-     $                    d_t_nlte)
+                CALL nltecool(klon, klev, nqmax, pplay*9.869e-6, t_seri,
+     $                    tr_seri, d_t_nlte)
             else if(nltemodel.eq.2) then                                
-                CALL nlte_tcool(klon,klev,pplay*9.869e-6,              
+               CALL nlte_tcool(klon,klev,pplay*9.869e-6,              
      $               t_seri,zzlay,co2vmr_gcm, n2vmr_gcm, covmr_gcm, 
      $               ovmr_gcm,d_t_nlte,ierr_nlte,varerr )
@@ -1162,8 +1287,7 @@
      $        CALL nlthermeq(klon, klev, paprs, pplay) 
 
-
-c          LTE radiative transfert / solar / IR matrix
-c          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
+c
+c       LTE radiative transfert / solar / IR matrix
+c       ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
       CALL radlwsw
      e            (dist, rmu0, fract, zzlev,
@@ -1171,12 +1295,11 @@
 
 
-c          CO2 near infrared absorption
-c          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+c       CO2 near infrared absorption
+c      ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
         d_t_nirco2(:,:)=0.
         if (callnirco2) then
-           call nirco2abs (klon, klev, pplay, dist, nqmax, qx,
-     .                 rmu0, fract, d_t_nirco2,
-     .                 co2vmr_gcm, ovmr_gcm)
+           call nirco2abs (klon, klev, pplay, dist, nqmax, tr_seri,
+     .                 rmu0, fract, d_t_nirco2)
         endif
 
@@ -1209,15 +1332,11 @@
         IF (callthermos) THEN
 
-c        call thermosphere(zplev,zplay,dist_sol,
-c     $     mu0,ptimestep,ptime,zday,tsurf,zzlev,zzlay,
-c     &     pt,pq,pu,pv,pdt,pdq,
-c     $     zdteuv,zdtconduc,zdumolvis,zdvmolvis,zdqmoldiff)
-
-           call euvheat(klon, klev, t_seri,paprs,pplay,zzlay,
-     $          rmu0,pdtphys,gmtime,rjourvrai,
-C     $			pq,pdq,zdteuv)
-     $          co2vmr_gcm, n2vmr_gcm, covmr_gcm, 
-     $          ovmr_gcm,nvmr_gcm,d_t_euv )
-
+c           call euvheat(klon, klev,t_seri,paprs,pplay,zzlay,
+c     $          rmu0,pdtphys,gmtime,rjourvrai, co2vmr_gcm, n2vmr_gcm, 
+c     $          covmr_gcm, ovmr_gcm,d_t_euv )
+           call euvheat(klon, klev, nqmax, t_seri,paprs,pplay,zzlay,
+     $         rmu0,pdtphys,gmtime,rjourvrai,
+     $         tr_seri, d_tr, d_t_euv )
+ 		
            DO k=1,klev
               DO ig=1,klon
@@ -1262,5 +1381,5 @@
             call molvis(klon, klev, pdtphys,
      $            pplay,paprs,t_seri,
-     $            v,tsurf,zzlev,zzlay,d_u_molvis)
+     $            v,tsurf,zzlev,zzlay,d_v_molvis)
 
             DO k=1,klev
@@ -1271,6 +1390,31 @@
                ENDDO
             ENDDO
-
-         ENDIF  ! callthermos
+        ENDIF
+
+
+!  --  MOLECULAR DIFFUSION ---
+
+          d_q_moldif(:,:,:)=0
+
+         IF (callthermos .and. ok_chem) THEN
+
+             call moldiff_red(klon, klev, nqmax,
+     &                   pplay,paprs,t_seri, tr_seri, pdtphys,
+     &                   zzlay,d_t_euv,d_t_conduc,d_q_moldif)
+
+
+! --- update tendencies tracers ---
+
+          DO iq = 1, nqmax
+           DO k=1,klev
+              DO ig=1,klon
+                tr_seri(ig,k,iq)= tr_seri(ig,k,iq)+ 
+     &                           d_q_moldif(ig,k,iq)*dtime ! [Kg/kg]?
+              ENDDO
+            ENDDO
+           ENDDO
+           
+
+         ENDIF  ! callthermos & ok_chem
 
 c====================================================================
