Index: /trunk/LMDZ.MARS/README
===================================================================
--- /trunk/LMDZ.MARS/README	(revision 1241)
+++ /trunk/LMDZ.MARS/README	(revision 1242)
@@ -2085,2 +2085,6 @@
 - Warnings in conf_phys modified accordingly.
 
+== 10/05/2014 == AS
+- Further reduction (or 'factorisation') of mesoscale includes
+- ustar and tstar are calculated in vdifc. so no need to recalculate
+  those elsewhere. use turb_mod module to store it. 
Index: unk/LMDZ.MARS/libf/phymars/meso_inc/meso_inc_les.F
===================================================================
--- /trunk/LMDZ.MARS/libf/phymars/meso_inc/meso_inc_les.F	(revision 1241)
+++ 	(revision )
@@ -1,77 +1,0 @@
-         if (.not.callrichsl) then
-       
-             DO ig=1,ngrid
-!!! sensible heat flux in W/m2
-!
-!             sensheat(ig) = zflubid(ig)-capcal(ig)*zdtsdif(ig)
-
-!! u star in similarity theory in m/s
-             ustar(ig) = 0.4
-     .               * sqrt( pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1) )
-     .               / log( 1.E+0 + zzlay(ig,1)/z0_default )
-             ENDDO
-
-         else
-
-            DO ig=1,ngrid
-
-!! New SL parametrization, correct formulation for sensheat :
-!
-!            sensheat(ig) = (pplay(ig,1)/(r*pt(ig,1)))*cpp
-!     &        *sqrt(pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1)
-!     &        + (log(1.+0.7*wstar(ig) + 2.3*wstar(ig)**2))**2)
-!     &        *zcdh(ig)*(tsurf(ig)-zh(ig,1))
-
-! New SL parametrization, ustar is more accurately computed in vdif_cd :
-
-            ustar(ig) = sqrt(zcdv(ig)*
-     &   (pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1) 
-     &     + (log(1.+0.7*wstar(ig) + 2.3*wstar(ig)**2))**2)
-     &                     )
-            ENDDO   
-
-         endif  !of if callrichsl
-
-!         write (*,*) 'PHYS HFX cp zdts', sensheat(100), zflubid(100), 
-!     .       capcal(100), 
-!     .       zdtsdif(100)
-!         write (*,*) 'PHYS UST', ustar(100) 
-
-!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-!!! LES LES 
-       IF (turb_resolved) THEN        
-
-         write (*,*) '************************************************' 
-         write (*,*) '** LES mode: the difv part is only used to'
-         write (*,*) '**  provide HFX and UST to the dynamics'
-         write (*,*) '** NB: - dudif, dvdif, dhdif, dqdif are set to 0'
-         write (*,*) '**     - tsurf is updated'     
-         write (*,*) '************************************************'
-
-         IF (lifting .and. doubleq) THEN 
-         !! lifted dust is injected in the first layer. 
-         !! Sedimentation must be called after turbulent mixing, i.e. on next step, after WRF. 
-         !! => lifted dust is not incremented before the sedimentation step.
-            zdqdif(:,1,:)=0.
-            zdqdif(:,1,igcm_dust_number)=-zdqsdif(:,igcm_dust_number)
-            zdqdif(:,1,igcm_dust_mass)=-zdqsdif(:,igcm_dust_mass)
-            zdqdif(:,2:nlayermx,:)=0.
-            DO iq=1, nq
-              IF ((iq .ne. igcm_dust_mass)
-     &        .and. (iq .ne. igcm_dust_number)) THEN
-                zdqsdif(:,iq)=0.
-              ENDIF
-            ENDDO
-         ELSE
-            DO ig=1,ngrid
-             DO l=1,nlayer
-              DO iq=1, nq
-                zdqdif(ig,l,iq) = 0.
-                zdqsdif(ig,iq) = 0. !! sortir de la boucle
-              ENDDO 
-             ENDDO
-            ENDDO
-         ENDIF
-      ENDIF
-!!! LES LES        
-!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
Index: /trunk/LMDZ.MARS/libf/phymars/physiq.F
===================================================================
--- /trunk/LMDZ.MARS/libf/phymars/physiq.F	(revision 1241)
+++ /trunk/LMDZ.MARS/libf/phymars/physiq.F	(revision 1242)
@@ -838,7 +838,4 @@
      &        zdqdif,zdqsdif,wstar,zcdv,zcdh,hfmax_th,sensibFlux)
 
-#ifdef MESOSCALE
-#include "meso_inc/meso_inc_les.F"
-#endif
 
           DO ig=1,ngrid
@@ -871,10 +868,36 @@
            ENDDO
           end if ! of if (tracer)
+         ELSE
+           write (*,*) '******************************************'
+           write (*,*) '** LES mode: the difv part is only used to'
+           write (*,*) '**  - provide HFX and UST to the dynamics'
+           write (*,*) '**  - update TSURF'
+           write (*,*) '******************************************'
+           !! Specific treatment for lifting in turbulent-resolving mode (AC)
+           IF (lifting .and. doubleq) THEN
+             !! lifted dust is injected in the first layer. 
+             !! Sedimentation must be called after turbulent mixing, i.e. on next step, after WRF. 
+             !! => lifted dust is not incremented before the sedimentation step.
+             zdqdif(1:ngrid,1,1:nq)=0.
+             zdqdif(1:ngrid,1,igcm_dust_number) = 
+     .                  -zdqsdif(1:ngrid,igcm_dust_number)
+             zdqdif(1:ngrid,1,igcm_dust_mass) = 
+     .                  -zdqsdif(1:ngrid,igcm_dust_mass)
+             zdqdif(1:ngrid,2:nlayermx,1:nq) = 0.
+             DO iq=1, nq
+              IF ((iq .ne. igcm_dust_mass)
+     &        .and. (iq .ne. igcm_dust_number)) THEN
+                zdqsdif(:,iq)=0.
+              ENDIF
+             ENDDO
+           ELSE
+             zdqdif(1:ngrid,1:nlayer,1:nq) = 0.
+             zdqsdif(1:ngrid,1:nq) = 0.
+           ENDIF
          ENDIF
-
       ELSE    
          DO ig=1,ngrid
             zdtsurf(ig)=zdtsurf(ig)+
-     s      (fluxrad(ig)+fluxgrd(ig))/capcal(ig)
+     s        (fluxrad(ig)+fluxgrd(ig))/capcal(ig)
          ENDDO
          IF (turb_resolved) THEN
@@ -1483,4 +1506,5 @@
      & g,zzlay,zzlev,zu,zv,wstar,hfmax_th,zmax_th,tsurf,zh,z_out,n_out,
      & T_out,u_out,ustar,tstar,L_mo,vhf,vvv)
+                   ! pourquoi ustar recalcule ici? fait dans vdifc.
 
 #ifndef MESOSCALE
@@ -1909,4 +1933,20 @@
       !state  real  RICE      ikj   misc  1  -  h  "RICE"      "ICE RADIUS"                      "m"
       comm_RICE(1:ngrid,1:nlayer) = rice(1:ngrid,1:nlayer)
+
+   
+      !! calculate sensible heat flux in W/m2 for outputs
+      !! -- the one computed in vdifc is not the real one
+      !! -- vdifc must have been called
+      if (.not.callrichsl) then
+        sensibFlux(1:ngrid) = zflubid(1:ngrid)
+     .         - capcal(1:ngrid)*zdtsdif(1:ngrid)
+      else
+        sensibFlux(1:ngrid) = 
+     &   (pplay(1:ngrid,1)/(r*pt(1:ngrid,1)))*cpp
+     &   *sqrt(pu(1:ngrid,1)*pu(1:ngrid,1)+pv(1:ngrid,1)*pv(1:ngrid,1)
+     &         +(log(1.+0.7*wstar(1:ngrid) + 2.3*wstar(1:ngrid)**2))**2)
+     &   *zcdh(1:ngrid)*(tsurf(1:ngrid)-zh(1:ngrid,1))
+      endif
+          
 #else
 #ifndef MESOINI
Index: /trunk/LMDZ.MARS/libf/phymars/turb_mod.F90
===================================================================
--- /trunk/LMDZ.MARS/libf/phymars/turb_mod.F90	(revision 1241)
+++ /trunk/LMDZ.MARS/libf/phymars/turb_mod.F90	(revision 1242)
@@ -6,4 +6,5 @@
   REAL,SAVE,ALLOCATABLE :: ustar(:)
   REAL,SAVE,ALLOCATABLE :: wstar(:)
+  REAL,SAVE,ALLOCATABLE :: tstar(:)
   REAL,SAVE,ALLOCATABLE :: hfmax_th(:)
   REAL,SAVE,ALLOCATABLE :: zmax_th(:)
@@ -25,4 +26,5 @@
     allocate(wstar(ngrid))
     allocate(ustar(ngrid))
+    allocate(tstar(ngrid))
     allocate(hfmax_th(ngrid))
     allocate(zmax_th(ngrid))
Index: /trunk/LMDZ.MARS/libf/phymars/vdifc.F
===================================================================
--- /trunk/LMDZ.MARS/libf/phymars/vdifc.F	(revision 1241)
+++ /trunk/LMDZ.MARS/libf/phymars/vdifc.F	(revision 1242)
@@ -13,5 +13,5 @@
       use surfdat_h, only: watercaptag, frost_albedo_threshold, dryness
       USE comcstfi_h
-      use turb_mod, only: turb_resolved
+      use turb_mod, only: turb_resolved, ustar, tstar
       IMPLICIT NONE
 
@@ -84,5 +84,4 @@
 
       REAL :: pt(ngrid,nlay)
-      REAL ust(ngrid),tst(ngrid)
  
       INTEGER ilev,ig,ilay,nlev
@@ -335,11 +334,11 @@
      &     (log(1.+0.7*wstar(:) + 2.3*wstar(:)**2))**2)
 
-           ust(:)=sqrt(zcdv_true(:))*sqrt(zu2(:)+
+           ustar(:)=sqrt(zcdv_true(:))*sqrt(zu2(:)+
      &     (log(1.+0.7*wstar(:) + 2.3*wstar(:)**2))**2)
 
-           tst(:)=0.
+           tstar(:)=0.
            DO ig=1,ngrid
               IF (zcdh_true(ig) .ne. 0.) THEN      ! When Cd=Ch=0, u*=t*=0
-                 tst(ig)=(ph(ig,1)-ptsrf(ig))*zcdh(ig)/ust(ig)
+                 tstar(ig)=(ph(ig,1)-ptsrf(ig))*zcdh(ig)/ustar(ig)
               ENDIF
            ENDDO
@@ -348,6 +347,6 @@
            zcdv(:)=zcdv_true(:)*sqrt(zu2(:))     ! 1 / bulk aerodynamic momentum conductance 
            zcdh(:)=zcdh_true(:)*sqrt(zu2(:))     ! 1 / bulk aerodynamic heat conductance
-           ust(:)=sqrt(zcdv_true(:))*sqrt(zu2(:))
-           tst(:)=(ph(:,1)-ptsrf(:))*zcdh_true(:)/sqrt(zcdv_true(:))
+           ustar(:)=sqrt(zcdv_true(:))*sqrt(zu2(:))
+           tstar(:)=(ph(:,1)-ptsrf(:))*zcdh_true(:)/sqrt(zcdv_true(:))
         ENDIF
 
@@ -383,5 +382,5 @@
       pt(:,:)=ph(:,:)*ppopsk(:,:)
       CALL yamada4(ngrid,nlay,nq,ptimestep,g,r,pplev,pt
-     s   ,pzlev,pzlay,pu,pv,ph,pq,zcdv_true,pq2,zkv,zkh,zkq,ust
+     s   ,pzlev,pzlay,pu,pv,ph,pq,zcdv_true,pq2,zkv,zkh,zkq,ustar
      s   ,9)
       ENDIF
