Index: trunk/LMDZ.MARS/libf/phymars/conf_phys.F
===================================================================
--- trunk/LMDZ.MARS/libf/phymars/conf_phys.F	(revision 3325)
+++ trunk/LMDZ.MARS/libf/phymars/conf_phys.F	(revision 3333)
@@ -376,6 +376,6 @@
          write(*,*)"albedo_perennialco2 = ",albedo_perennialco2
 
-         write(*,*)"Using lag layer??"
-         lag_layer=.false.
+         write(*,*)"Include water buoyancy effect??"
+         include_waterbuoyancy=.false.
          call getin_p("include_waterbuoyancy",include_waterbuoyancy)
          write(*,*) "include_waterbuoyancy = ", include_waterbuoyancy
Index: trunk/LMDZ.MARS/libf/phymars/paleoclimate_mod.F90
===================================================================
--- trunk/LMDZ.MARS/libf/phymars/paleoclimate_mod.F90	(revision 3325)
+++ trunk/LMDZ.MARS/libf/phymars/paleoclimate_mod.F90	(revision 3333)
@@ -20,6 +20,6 @@
     real,    save                              :: albedo_perennialco2     ! Albedo for perennial co2 ice [1]
     logical, save                              :: lag_layer               ! Does lag layer is present?
-    logical, save                              :: include_waterbuoyancy    ! Include the effect of water buoyancy when computing the sublimation of water ice ?
-!$OMP THREADPRIVATE(h2o_ice_depth,d_coef,lag_co2_ice,albedo_perennialco2,include_waterbuoyancy)
+    logical, save                              :: include_waterbuoyancy   ! Include the effect of water buoyancy when computing the sublimation of water ice ?
+!$OMP THREADPRIVATE(h2o_ice_depth,lag_co2_ice,d_coef,albedo_perennialco2,lag_layer,include_waterbuoyancy)
 
 !=======================================================================
Index: trunk/LMDZ.MARS/libf/phymars/vdif_cd_mod.F90
===================================================================
--- trunk/LMDZ.MARS/libf/phymars/vdif_cd_mod.F90	(revision 3325)
+++ trunk/LMDZ.MARS/libf/phymars/vdif_cd_mod.F90	(revision 3333)
@@ -157,4 +157,6 @@
       pcdh(:,:) = 0.
       z0t(:,:) = 0.
+      fm(:,:) = 0.
+      fh(:,:) = 0.
       f_ri_cd_min = 0.01
 ! this formulation assumes alphah=1., implying betah=betam
@@ -203,9 +205,8 @@
             DO ig=1,ngrid
                ite = 0.
-               residual = abs(pz0tcomp(ig)-pz0t) 
+               residual = 100*tol_iter*pz0(ig) 
                z1z0=pz(ig,1)/pz0(ig)
                cdn(ig)=karman/log(z1z0)
                cdn(ig)=cdn(ig)*cdn(ig)
-            
                DO WHILE((residual .gt. tol_iter*pz0(ig)) .and.  (ite .lt. itemax)) 
 ! Computations of z0T; iterated until z0T converges  
@@ -225,5 +226,4 @@
                       rib_dry(ig,islope) = ric_colaitis
                    ENDIF
-        
 ! Compute the stability functions fm; fh depending on the stability of the surface layer