Index: trunk/mars/README
===================================================================
--- trunk/mars/README	(revision 91)
+++ trunk/mars/README	(revision 117)
@@ -622,2 +622,12 @@
 >> new version version for aeroptproperties.F in phymars to limit uncertainties and be able to play with ngau
 >> this was coded by JBM in his personal reference version but not transmitted to the team reference version
+
+== 10/05/2011 == AS + JF
+>> in newsedim.F used for mesoscale computations, spurious values close to the surface
+   --> this was related to 1 - exp(-x) calculated as zero in w(ig,l) if x is very small
+   --> fix: when this happens, replace exp(-x) by 1 - x since x ~ 0
+>> in newsedim.F, "if (dztop.gt.epaisseur(ig,l)) then" was closed too soon by an "endif"
+   --> hence basically the simple method was never used 
+       and useless calculations with the complex method were carried out
+   --> fix by moving the closing "endif" [in addition to corrections mentioned in the previous point]
+
Index: trunk/mars/libf/phymars/newsedim.F
===================================================================
--- trunk/mars/libf/phymars/newsedim.F	(revision 91)
+++ trunk/mars/libf/phymars/newsedim.F	(revision 117)
@@ -67,5 +67,4 @@
 
 
-
 c    ** un petit test de coherence
 c       --------------------------
@@ -138,16 +137,32 @@
              k=0
 
+            w(ig,l) = 0. !! JF+AS ajout initialisation
 c **************************************************************
 c            Simple Method
-             w(ig,l) =
-     &       (1- exp(-dztop*g/(r*pt(ig,l))))*pplev(ig,l) / g
-cc           write(*,*) 'OK simple method l,w =', l, w(ig,l)
-cc           write(*,*) 'OK simple method dztop =', dztop
+
+cc             w(ig,l) =
+cc     &       (1.- exp(-dztop*g/(r*pt(ig,l))))*pplev(ig,l) / g
+cccc           write(*,*) 'OK simple method l,w =', l, w(ig,l)
+cccc           write(*,*) 'OK simple method dztop =', dztop
+
+             w(ig,l) = 1. - exp(-dztop*g/(r*pt(ig,l)))
+             !!! Diagnostic: JF. Fix: AS. Date: 05/11 
+             !!! Probleme arrondi avec la quantite ci-dessus
+             !!! ---> vaut 0 pour -dztop*g/(r*pt(ig,l)) trop petit
+             !!! ---> dans ce cas on utilise le developpement limite !
+             !!! ---> exp(-x) = 1 - x lorsque x --> 0 avec une erreur de x^2 / 2            
+             IF ( w(ig,l) .eq. 0. ) THEN
+                w(ig,l) = ( dztop*g/(r*pt(ig,l)) ) * pplev(ig,l) / g
+             ELSE
+                w(ig,l) = w(ig,l) * pplev(ig,l) / g  
+             ENDIF
+
+
 c **************************************************************
 cccc         Complex method :
-            if (dztop.gt.epaisseur(ig,l)) then
+            if (dztop.gt.epaisseur(ig,l)) then                !!!if on traverse plus d'une couche
 cccc            Cas ou on "epuise" la couche l : On calcule le flux
 cccc            Venant de dessus en tenant compte de la variation de Vstokes
-
+c **************************************************************
                Ep= epaisseur(ig,l)
                Stra= traversee(ig,l)
@@ -159,8 +174,19 @@
                enddo 
                Ep = Ep - epaisseur(ig,l+k)
-             end if
-             ptop=pplev(ig,l+k)*exp(-(dztop-Ep)*g/(r*pt(ig,l+k)))
-             w(ig,l) = (pplev(ig,l) -Ptop)/g
-c
+             !ptop=pplev(ig,l+k)*exp(-(dztop-Ep)*g/(r*pt(ig,l+k)))
+
+             !!! JF+AS 05/11 Probleme arrondi potentiel, meme solution que ci-dessus
+             ptop=exp(-(dztop-Ep)*g/(r*pt(ig,l+k)))
+             IF ( ptop .eq. 1. ) THEN
+                PRINT*, 'exposant trop petit ', ig, l
+                ptop=pplev(ig,l+k) * ( 1. - (dztop-Ep)*g/(r*pt(ig,l+k)))
+             ELSE
+                ptop=pplev(ig,l+k) * ptop 
+             ENDIF
+             w(ig,l) = (pplev(ig,l) - Ptop)/g
+
+            endif                !!!!!if complex method
+
+
 cc           write(*,*) 'OK new    method l,w =', l, w(ig,l)
 cc           write(*,*) 'OK new    method dztop =', dztop
@@ -171,4 +197,6 @@
 cc            if(l.eq.8)write(*,*)'l=8,k, w',k, w(1,l)
 c **************************************************************
+
+
         end do
       end do
