Index: /LMDZ5/trunk/libf/phylmd/phys_local_var_mod.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/phys_local_var_mod.F90	(revision 2573)
+++ /LMDZ5/trunk/libf/phylmd/phys_local_var_mod.F90	(revision 2574)
@@ -408,5 +408,6 @@
       allocate(t_seri(klon,klev),q_seri(klon,klev),ql_seri(klon,klev),qs_seri(klon,klev))
       allocate(u_seri(klon,klev),v_seri(klon,klev))
-      allocate(l_mixmin(klon,klev,nbsrf), l_mix(klon,klev,nbtr))
+      allocate(l_mixmin(klon,klev,nbsrf), l_mix(klon,klev,nbsrf))
+      l_mix(:,:,:)=0. ; l_mixmin(:,:,:)=0. ! doit etre initialse car pas toujours remplis
 
       allocate(tr_seri(klon,klev,nbtr))
Index: /LMDZ5/trunk/libf/phylmd/yamada4.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/yamada4.F90	(revision 2573)
+++ /LMDZ5/trunk/libf/phylmd/yamada4.F90	(revision 2574)
@@ -212,7 +212,7 @@
 ! initialize arrays:
 
-  m2(:, :) = 0.0
-  sm(:, :) = 0.0
-  rif(:, :) = 0.0
+  m2(1:ngrid, :) = 0.0
+  sm(1:ngrid, :) = 0.0
+  rif(1:ngrid, :) = 0.0
 
 !------------------------------------------------------------
@@ -264,5 +264,5 @@
  
     DO k = 2, klev
-      q2(:, k) = l(:, k)**2*zz(:, k)
+      q2(1:ngrid, k) = l(1:ngrid, k)**2*zz(1:ngrid, k)
     END DO
 
@@ -349,9 +349,9 @@
 
     DO k = 2, klev - 1
-      km(:, k) = l(:, k)*sqrt(q2(:,k))*sm(:, k)
-      q2(:, k) = q2(:, k) + dt*km(:, k)*m2(:, k)*(1.-rif(:,k))
-      q2(:, k) = min(max(q2(:,k),1.E-10), 1.E4)
-      q2(:, k) = 1./(1./sqrt(q2(:,k))+dt/(2*l(:,k)*b1))
-      q2(:, k) = q2(:, k)*q2(:, k)
+      km(1:ngrid, k) = l(1:ngrid, k)*sqrt(q2(1:ngrid,k))*sm(1:ngrid, k)
+      q2(1:ngrid, k) = q2(1:ngrid, k) + dt*km(1:ngrid, k)*m2(1:ngrid, k)*(1.-rif(1:ngrid,k))
+      q2(1:ngrid, k) = min(max(q2(1:ngrid,k),1.E-10), 1.E4)
+      q2(1:ngrid, k) = 1./(1./sqrt(q2(1:ngrid,k))+dt/(2*l(1:ngrid,k)*b1))
+      q2(1:ngrid, k) = q2(1:ngrid, k)*q2(1:ngrid, k)
     END DO
 
@@ -392,5 +392,5 @@
 
   IF (iflag_pbl>=12) THEN
-    q2(:, 1) = q2(:, 2)
+    q2(1:ngrid, 1) = q2(1:ngrid, 2)
     CALL vdif_q2(dt, g, rconst, ngrid, plev, temp, kq, q2)
   END IF
@@ -803,7 +803,5 @@
  famorti(zzzz, zh_oro, zhlim)=(-1.)*ATAN((zzzz-zh_oro)/(zhlim-zh_oro))*2./3.1416+1.    
 
-
- 
-
+  IF (ngrid==0) RETURN
 
   IF (firstcall) THEN
@@ -811,5 +809,4 @@
     firstcall = .FALSE.
   END IF
-
 
 
@@ -853,5 +850,5 @@
     !......................................................................
 
-    l0(:) = 150.
+    l0(1:ngrid) = 150.
     DO k = 2, klev
       DO ig = 1, ngrid
@@ -868,7 +865,7 @@
 !=================================================================================
 
-   l2(:,:)=0.0
-   l_mixmin(:,:,nsrf)=0.
-   l_mix(:,:,nsrf)=0.
+   l2(1:ngrid,:)=0.0
+   l_mixmin(1:ngrid,:,nsrf)=0.
+   l_mix(1:ngrid,:,nsrf)=0.
 
    IF (nsrf .EQ. 1) THEN
@@ -877,23 +874,19 @@
 !--------------
 
-  extent=2.                                                         ! On ??tend l'impact du relief jusqu'?? extent*h, extent >1.  
-  lmax=150.                                                         ! Longueur de m??lange max dans l'absolu
+     extent=2.                                                         ! On ??tend l'impact du relief jusqu'?? extent*h, extent >1.  
+     lmax=150.                                                         ! Longueur de m??lange max dans l'absolu
 
 ! calculs
 !---------
 
-  DO ig=1,ngrid
+     DO ig=1,ngrid
 
       ! On calcule la hauteur du relief
       !.................................
-
       ! On ne peut pas prendre zstd seulement pour caracteriser le relief sous maille
       ! car sur un terrain pentu mais sans relief, zstd est non nul (comme en Antarctique, C. Genthon)
       ! On corrige donc zstd avec l'ecart type de l'altitude dans la maille sans relief 
       ! (en gros, une maille de taille xcell avec une pente constante zstdslope)
- 
-
       jg=ni(ig) 
-  
 !     IF (zsig(jg) .EQ. 0.) THEN
 !          zstdslope(ig)=0.         
@@ -907,49 +900,30 @@
       h_oro(ig)=zstd(jg)
       hlim(ig)=extent*h_oro(ig)     
-
-     
-
-
-   END DO
-
-  l2limit(1:ngrid)=0.
-
-  DO k=2,klev
-   DO ig=1,ngrid
- 
-      winds(ig,k)=sqrt(u(ig,k)**2+v(ig,k)**2)
-
-      IF (zlev(ig,k) .LE. h_oro(ig)) THEN                                     ! Si on est sous l'orographie
-
-
-      l2strat= kapb*pbl_lmixmin_alpha*winds(ig,k)/sqrt(max(n2(ig,k),1.E-10))              ! si stratifi??, amplitude d'oscillation * kappab (voir Van de Wiel et al 2008)
-                                                
-      l2neutre=kap*zlev(ig,k)*h_oro(ig)/(kap*zlev(ig,k)+h_oro(ig))            ! Dans le cas neutre, formule de blackadar qui tend asymptotiquement vers h
-      l2neutre=MIN(l2neutre,lmax)                                             ! On majore par lmax 
-      l2limit(ig)=MIN(l2neutre,l2strat)                                       ! Calcule de l2 (minimum de la longueur en cas neutre et celle en situation stratifi??e)
-      
-      l2(ig,k)=l2limit(ig)
+     ENDDO
+
+     l2limit(1:ngrid)=0.
+
+     DO k=2,klev
+        DO ig=1,ngrid
+           winds(ig,k)=sqrt(u(ig,k)**2+v(ig,k)**2)
+           IF (zlev(ig,k) .LE. h_oro(ig)) THEN  ! sous l'orographie
+              l2strat= kapb*pbl_lmixmin_alpha*winds(ig,k)/sqrt(max(n2(ig,k),1.E-10))  ! si stratifi??, amplitude d'oscillation * kappab (voir Van de Wiel et al 2008)
+              l2neutre=kap*zlev(ig,k)*h_oro(ig)/(kap*zlev(ig,k)+h_oro(ig))            ! Dans le cas neutre, formule de blackadar. tend asymptotiquement vers h
+              l2neutre=MIN(l2neutre,lmax)                                             ! On majore par lmax 
+              l2limit(ig)=MIN(l2neutre,l2strat)                                       ! Calcule de l2 (minimum de la longueur en cas neutre et celle en situation stratifi??e)
+              l2(ig,k)=l2limit(ig)
                                       
-      ELSE IF (zlev(ig,k) .LE. hlim(ig)) THEN ! Si on est au dessus des montagnes, mais affect?? encore par elles
+           ELSE IF (zlev(ig,k) .LE. hlim(ig)) THEN ! Si on est au dessus des montagnes, mais affect?? encore par elles
 
       ! Au dessus des montagnes, on prend la l2limit au sommet des montagnes 
       ! (la derni??re calcul??e dans la boucle k, vu que k est un indice croissant avec z)
       ! et on multiplie l2limit par une fonction qui d??croit entre h et hlim
-
-         l2(ig,k)=l2limit(ig)*famorti(zlev(ig,k),h_oro(ig), hlim(ig))
-
-      
-      ELSE                                                                    ! Au dessus de extent*h, on prend l2=l0
-
-      l2(ig,k)=0.
-
-
-      END IF
-        
-    END DO
-  END DO
-
-
- END IF                                                                        ! pbl_lmixmin_alpha
+              l2(ig,k)=l2limit(ig)*famorti(zlev(ig,k),h_oro(ig), hlim(ig))
+           ELSE                                                                    ! Au dessus de extent*h, on prend l2=l0 
+              l2(ig,k)=0.
+           END IF
+        ENDDO
+     ENDDO
+   ENDIF                                                                        ! pbl_lmixmin_alpha
 
 !==================================================================================
@@ -961,21 +935,21 @@
  DO k=2,klev
     DO ig=1,ngrid
-
-   lmix(ig,k)=MAX(MAX(l1(ig,k), l2(ig,k)),lmixmin)
-
-   END DO
- END DO
+       lmix(ig,k)=MAX(MAX(l1(ig,k), l2(ig,k)),lmixmin)
+   ENDDO
+ ENDDO
+
+! Diagnostics
 
  DO k=2,klev
-  DO ig=1,ngrid
-   jg=ni(ig)
-   l_mix(jg,k,nsrf)=lmix(ig,k)
-   l_mixmin(jg,k,nsrf)=l2(ig,k)
-  ENDDO
+    DO ig=1,ngrid
+       jg=ni(ig)
+       l_mix(jg,k,nsrf)=lmix(ig,k)
+       l_mixmin(jg,k,nsrf)=l2(ig,k)
+    ENDDO
  ENDDO
-  DO ig=1,ngrid
-   jg=ni(ig)
-      l_mix(jg,1,nsrf)=hlim(ig)
-  ENDDO
+ DO ig=1,ngrid
+    jg=ni(ig)
+    l_mix(jg,1,nsrf)=hlim(ig)
+ ENDDO
 
 
