Index: trunk/LMDZ.MARS/libf/phymars/surflayer_interpol.F
===================================================================
--- trunk/LMDZ.MARS/libf/phymars/surflayer_interpol.F	(revision 338)
+++ trunk/LMDZ.MARS/libf/phymars/surflayer_interpol.F	(revision 339)
@@ -32,4 +32,5 @@
 
 #include "comcstfi.h"
+#include "callkeys.h"
 
 !   Arguments:
@@ -77,5 +78,5 @@
       REAL pz0tcomp(ngrid)     ! computed z0t
       REAL ite
-      REAL residual
+      REAL residual,zcd0,z1
       REAL pcdv(ngrid),pcdh(ngrid)
 ! For output :
@@ -85,7 +86,12 @@
 !   couche de surface:
 !   ------------------
-      zout=z_out 
+
+c Init :
+
+      L_mo(:)=0.
+      ustar(:)=0.
       tstar(:)=0.
-      ustar(:)=0.
+      zout=z_out
+
       reynolds(:)=0.
 
@@ -140,4 +146,10 @@
 !        ENDIF
 
+      if(.not.callrichsl) then
+         rib(ig) = (pg/ph(ig,1))
+     &      *sqrt(pz(ig,1)*pz0(ig))
+     &      *(((log(pz(ig,1)/pz0(ig)))**2)/(log(pz(ig,1)/pz0t)))
+     &      *(ph(ig,1)-pts(ig))/(zu2(ig)+6.)
+      else
           rib(ig) = (pg/ph(ig,1))
 !     &      *pz(ig,1)*pz0(ig)/sqrt(pz(ig,1)*pz0t)
@@ -148,4 +160,6 @@
 !     &  /(MAX(pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1),wmax(ig)**2))
 !     &  /( pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1) + wmax(ig)**2)
+
+      endif
 
          else
@@ -211,9 +225,11 @@
       ENDDO
 
+!      endif !of if callrichsl
+
 ! Large-scale wind at first layer (without gustiness) and
 ! u* theta* computation
       DO ig=1,ngrid
 
-         if (rib(ig) .gt. ric) then
+         if (rib(ig) .ge. ric) then
            ustar(ig)=0.
            tstar(ig)=0.
@@ -230,5 +246,5 @@
      &     /sqrt(pcdv(ig))
 
-           if (tstar(ig) .lt. -50) then
+           if ((tstar(ig) .lt. -50) .and. callrichsl) then
              print*, fh(ig),rib(ig),(ph(ig,1)-pts(ig))
      &               ,log(pz(ig,1)/pz0tcomp(ig)),sqrt(fm(ig))
@@ -329,6 +345,19 @@
       ENDDO
 
+! when using convective adjustment without thermals, a vertical potential temperature
+! profile is assumed up to the thermal roughness length. Hence, theoretically, theta
+! interpolated at any height in the surface layer is theta at the first level.
+
+      IF ((.not.calltherm).and.(calladj)) THEN
+           Teta_out(:)=ph(:,1)
+      ENDIF
+
 ! Usefull diagnostics for the interpolation routine :
 
+!
+!         call WRITEDIAGFI(ngrid,'Ri',
+!     &              'Richardson','m',
+!     &                         2,rib)
+!
 !         call WRITEDIAGFI(ngrid,'z0T',
 !     &              'thermal roughness length','m',