Index: /LMDZ5/trunk/libf/phylmd/dyn1d/1D_interp_cases.h
===================================================================
--- /LMDZ5/trunk/libf/phylmd/dyn1d/1D_interp_cases.h	(revision 2306)
+++ /LMDZ5/trunk/libf/phylmd/dyn1d/1D_interp_cases.h	(revision 2307)
@@ -1,2 +1,5 @@
+!
+! $Id$
+!
 !---------------------------------------------------------------------
 ! Interpolation forcing in time and onto model levels
@@ -615,10 +618,10 @@
      &       ,vitw_cas,du_cas,hu_cas,vu_cas           &
      &       ,dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dtrad_cas             &
-     &       ,dq_cas,hq_cas,vq_cas,lat_cas,sens_cas                 &
+     &       ,dq_cas,hq_cas,vq_cas,lat_cas,sens_cas,ustar_cas       &
      &       ,ts_prof_cas,plev_prof_cas,t_prof_cas,q_prof_cas,u_prof_cas,v_prof_cas         &
      &       ,ug_prof_cas,vg_prof_cas,vitw_prof_cas,du_prof_cas,hu_prof_cas,vu_prof_cas     &
      &       ,dv_prof_cas,hv_prof_cas,vv_prof_cas,dt_prof_cas,ht_prof_cas,vt_prof_cas       &
      &       ,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas,lat_prof_cas               &
-     &       ,sens_prof_cas)
+     &       ,sens_prof_cas,ustar_prof_cas)
 
              ts_cur = ts_prof_cas 
@@ -746,4 +749,5 @@
         if (tend_rayo.eq.1) then
            dt_cooling(l) = dtrad_mod_cas(l)
+!          print *,'dt_cooling=',dt_cooling(l)
         else
            dt_cooling(l) = 0.0
Index: /LMDZ5/trunk/libf/phylmd/dyn1d/1D_read_forc_cases.h
===================================================================
--- /LMDZ5/trunk/libf/phylmd/dyn1d/1D_read_forc_cases.h	(revision 2306)
+++ /LMDZ5/trunk/libf/phylmd/dyn1d/1D_read_forc_cases.h	(revision 2307)
@@ -1,2 +1,5 @@
+!
+! $Id$
+!
 !----------------------------------------------------------------------
 ! forcing_les = .T. : Impose a constant cooling 
@@ -738,9 +741,9 @@
      &         ,ug_cas,vg_cas,vitw_cas,du_cas,hu_cas,vu_cas           &
      &         ,dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dtrad_cas             &
-     &         ,dq_cas,hq_cas,vq_cas,lat_cas,sens_cas                 &
+     &         ,dq_cas,hq_cas,vq_cas,lat_cas,sens_cas,ustar_cas       &
      &         ,ts_prof_cas,plev_prof_cas,t_prof_cas,q_prof_cas,u_prof_cas,v_prof_cas         &
      &         ,ug_prof_cas,vg_prof_cas,vitw_prof_cas,du_prof_cas,hu_prof_cas,vu_prof_cas     &
      &         ,dv_prof_cas,hv_prof_cas,vv_prof_cas,dt_prof_cas,ht_prof_cas,vt_prof_cas,dtrad_prof_cas       &
-     &         ,dq_prof_cas,hq_prof_cas,vq_prof_cas,lat_prof_cas,sens_prof_cas)
+     &         ,dq_prof_cas,hq_prof_cas,vq_prof_cas,lat_prof_cas,sens_prof_cas,ustar_prof_cas)
 
 ! vertical interpolation using TOGA interpolation routine:
@@ -778,5 +781,15 @@
        d_u_adv(l) = (hv_mod_cas(l)+vv_mod_cas(l))
       enddo     
-       
+
+! In case fluxes are imposed
+       IF (ok_flux_surf) THEN
+       fsens=sens_prof_cas
+       flat=lat_prof_cas
+       ENDIF
+       IF (ok_prescr_ust) THEN
+       ust=ustar_prof_cas
+       print *,'ust=',ust
+       ENDIF
+
       endif !forcing_case
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
Index: /LMDZ5/trunk/libf/phylmd/mod_1D_cases_read.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/mod_1D_cases_read.F90	(revision 2306)
+++ /LMDZ5/trunk/libf/phylmd/mod_1D_cases_read.F90	(revision 2307)
@@ -1,2 +1,5 @@
+!
+! $Id$
+!
 MODULE mod_1D_cases_read
 
@@ -11,7 +14,10 @@
         real heure_ini_cas
         real day_ju_ini_cas   ! Julian day of case first day
-        parameter (year_ini_cas=2011)
-        parameter (mth_ini_cas=10)
-        parameter (day_ini_cas=1)  ! 10 = 10Juil2006
+!       parameter (year_ini_cas=2011)
+        parameter (year_ini_cas=1969)
+!       parameter (mth_ini_cas=10)
+        parameter (mth_ini_cas=6)
+!       parameter (day_ini_cas=1)  ! 10 = 10Juil2006
+        parameter (day_ini_cas=24)  ! 24 = 24 juin 1969
         parameter (heure_ini_cas=0.) !0h en secondes
         real pdt_cas
@@ -43,5 +49,5 @@
         real, allocatable::  vitw_cas(:,:)
         real, allocatable::  ug_cas(:,:),vg_cas(:,:)
-        real, allocatable::  lat_cas(:),sens_cas(:),ts_cas(:)
+        real, allocatable::  lat_cas(:),sens_cas(:),ts_cas(:),ustar_cas(:)
 
 !champs interpoles
@@ -70,5 +76,5 @@
 
 
-        real lat_prof_cas,sens_prof_cas,ts_prof_cas
+        real lat_prof_cas,sens_prof_cas,ts_prof_cas,ustar_prof_cas
       
 
@@ -143,5 +149,5 @@
         allocate(ug_cas(nlev_cas,nt_cas))
         allocate(vg_cas(nlev_cas,nt_cas))
-        allocate(lat_cas(nt_cas),sens_cas(nt_cas),ts_cas(nt_cas))
+        allocate(lat_cas(nt_cas),sens_cas(nt_cas),ts_cas(nt_cas),ustar_cas(nt_cas))
 
 
@@ -175,5 +181,5 @@
      &     ,ug_cas,vg_cas,vitw_cas,du_cas,hu_cas,vu_cas,dv_cas,hv_cas,vv_cas &
      &     ,dt_cas,dtrad_cas,ht_cas,vt_cas,dq_cas,hq_cas,vq_cas             &
-     &     ,dth_cas,hth_cas,vth_cas,dr_cas,hr_cas,vr_cas,sens_cas,lat_cas,ts_cas)
+     &     ,dth_cas,hth_cas,vth_cas,dr_cas,hr_cas,vr_cas,sens_cas,lat_cas,ts_cas,ustar_cas)
 
 
@@ -203,5 +209,5 @@
         deallocate(ug_cas)
         deallocate(vg_cas)
-        deallocate(lat_cas,sens_cas,ts_cas)
+        deallocate(lat_cas,sens_cas,ts_cas,ustar_cas)
 
 !champs interpoles
@@ -241,5 +247,5 @@
      &     ,zz,pp,temp,qv,rh,theta,rv,u,v,ug,vg,w,                   &
      &     du,hu,vu,dv,hv,vv,dt,dtrad,ht,vt,dq,hq,vq,                     &
-     &     dth,hth,vth,dr,hr,vr,sens,flat,ts)
+     &     dth,hth,vth,dr,hr,vr,sens,flat,ts,ustar)
 
 !program reading forcing of the case study
@@ -265,10 +271,10 @@
       real dth(nlevel,ntime),hth(nlevel,ntime),vth(nlevel,ntime)
       real dr(nlevel,ntime),hr(nlevel,ntime),vr(nlevel,ntime)
-      real flat(ntime),sens(ntime),ts(ntime)
+      real flat(ntime),sens(ntime),ts(ntime),ustar(ntime)
 
 
       integer nid, ierr,rid
       integer nbvar3d
-      parameter(nbvar3d=34)
+      parameter(nbvar3d=39)
       integer var3didin(nbvar3d)
 
@@ -478,4 +484,10 @@
            stop 'ts'
          endif
+
+      ierr=NF_INQ_VARID(nid,"ustar",var3didin(35))
+         if(ierr/=NF_NOERR) then
+           write(*,*) NF_STRERROR(ierr)
+           stop 'ustar'
+         endif
  
 #ifdef NC_DOUBLE
@@ -853,4 +865,16 @@
          endif
 !          write(*,*)'lecture ts ok',ts
+
+#ifdef NC_DOUBLE
+         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(35),ustar)
+#else
+         ierr = NF_GET_VAR_REAL(nid,var3didin(35),ustar)
+#endif
+         if(ierr/=NF_NOERR) then
+            write(*,*) NF_STRERROR(ierr)
+            stop "getvarup"
+         endif
+!         write(*,*)'lecture ustar ok',ustar
+
 
          return 
@@ -862,5 +886,5 @@
      &         ,ug_cas,vg_cas,vitw_cas,du_cas,hu_cas,vu_cas           &
      &         ,dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dtrad_cas   &
-     &         ,dq_cas,hq_cas,vq_cas,lat_cas,sens_cas          &
+     &         ,dq_cas,hq_cas,vq_cas,lat_cas,sens_cas,ustar_cas       &
      &         ,ts_prof_cas,plev_prof_cas,t_prof_cas,q_prof_cas       &
      &         ,u_prof_cas,v_prof_cas,ug_prof_cas,vg_prof_cas         &
@@ -868,5 +892,5 @@
      &         ,dv_prof_cas,hv_prof_cas,vv_prof_cas,dt_prof_cas       &
      &         ,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas    &
-     &         ,hq_prof_cas,vq_prof_cas,lat_prof_cas,sens_prof_cas)
+     &         ,hq_prof_cas,vq_prof_cas,lat_prof_cas,sens_prof_cas,ustar_prof_cas)
           
 
@@ -903,4 +927,5 @@
         real lat_cas(nt_cas)
         real sens_cas(nt_cas)
+        real ustar_cas(nt_cas)
 
 ! outputs:
@@ -915,5 +940,5 @@
         real dtrad_prof_cas(nlev_cas)
         real dq_prof_cas(nlev_cas),hq_prof_cas(nlev_cas),vq_prof_cas(nlev_cas)
-        real lat_prof_cas,sens_prof_cas,ts_prof_cas
+        real lat_prof_cas,sens_prof_cas,ts_prof_cas,ustar_prof_cas
 ! local:
         integer it_cas1, it_cas2,k
@@ -921,20 +946,23 @@
 
 
-        print*,'Check time',day1,day_ini_cas,day_ini_cas+90
+        print*,'Check time',day1,day_ini_cas,day_ini_cas+1
 
         if ((forcing_type.eq.10).and.(1.eq.1)) then
 ! Check that initial day of the simulation consistent with the case:
-       if (annee_ref.ne.2011) then
-        print*,'Pour CINDY, annee_ref doit etre 2011'
+       if (annee_ref.ne.1969) then
+!       print*,'Pour CINDY, annee_ref doit etre 2011'
+        print*,'Pour BOMEX, annee_ref doit etre 1969'
         print*,'Changer annee_ref dans run.def'
         stop
        endif
-       if (annee_ref.eq.2011 .and. day1.lt.day_ini_cas) then
-        print*,'CINDY a débuté le 1 octobre 2011',day1,day_ini_cas
+       if (annee_ref.eq.1969 .and. day1.lt.day_ini_cas) then
+!       print*,'CINDY a debute le 1 octobre 2011',day1,day_ini_cas
+        print*,'BOMEX a debute le 24 juin 1969',day1,day_ini_cas
         print*,'Changer dayref dans run.def'
         stop
        endif
-       if (annee_ref.eq.2011 .and. day1.gt.day_ini_cas+90) then
-        print*,'CINDY a fini le 31 decembre'
+       if (annee_ref.eq.1969 .and. day1.gt.day_ini_cas+1) then
+!       print*,'CINDY a fini le 31 decembre'
+        print*,'BOMEX a fini le 24 juin'
         print*,'Changer dayref ou nday dans run.def',day1,day_ini_cas+90
         stop
@@ -950,5 +978,5 @@
        endif
        if (annee_ref.eq.2006 .and. day1.lt.day_ini_cas) then
-        print*,'AMMA a débuté le 10 juillet 2006',day1,day_ini_cas
+        print*,'AMMA a debute le 10 juillet 2006',day1,day_ini_cas
         print*,'Changer dayref dans run.def'
         stop
@@ -999,6 +1027,8 @@
        sens_prof_cas = sens_cas(it_cas2)                                     &
      &          -frac*(sens_cas(it_cas2)-sens_cas(it_cas1))
-       ts_prof_cas = ts_cas(it_cas2)                                     &
+       ts_prof_cas = ts_cas(it_cas2)                                         &
      &          -frac*(ts_cas(it_cas2)-ts_cas(it_cas1))
+       ustar_prof_cas = ustar_cas(it_cas2)                                   &
+     &          -frac*(ustar_cas(it_cas2)-ustar_cas(it_cas1))
 
        do k=1,nlev_cas
Index: /LMDZ5/trunk/libf/phylmd/pbl_surface_mod.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/pbl_surface_mod.F90	(revision 2306)
+++ /LMDZ5/trunk/libf/phylmd/pbl_surface_mod.F90	(revision 2307)
@@ -1318,5 +1318,5 @@
 !     *sqrt(yu(i,1)*yu(i,1)+yv(i,1)*yv(i,1)))
        ycdragm(i) = ust*ust/(1.+vent)/vent
-       print *,'ycdragm ust yu yv apres=',ycdragm(i),ust,yu(i,1),yv(i,1)
+!      print *,'ycdragm ust yu yv apres=',ycdragm(i),ust,yu(i,1),yv(i,1)
       ENDDO
      ENDIF
@@ -1343,8 +1343,8 @@
         IF (ok_prescr_ust) then
          DO i = 1, knon
-          print *,'ycdragm_x avant=',ycdragm_x(i)
+!         print *,'ycdragm_x avant=',ycdragm_x(i)
           vent= sqrt(yu_x(i,1)*yu_x(i,1)+yv_x(i,1)*yv_x(i,1))
           ycdragm_x(i) = ust*ust/(1.+vent)/vent
-          print *,'ycdragm_x ust yu yv apres=',ycdragm_x(i),ust,yu_x(i,1),yv_x(i,1)
+!         print *,'ycdragm_x ust yu yv apres=',ycdragm_x(i),ust,yu_x(i,1),yv_x(i,1)
          ENDDO
         ENDIF
@@ -1367,11 +1367,11 @@
             ycdragm_w, ycdragh_w, zri1_w, pref_w )
 
-! --- special Dice. JYG+MPL 25112013
+! --- special Dice. JYG+MPL 25112013 puis BOMEX
         IF (ok_prescr_ust) then
          DO i = 1, knon
-          print *,'ycdragm_w avant=',ycdragm_w(i)
+!         print *,'ycdragm_w avant=',ycdragm_w(i)
           vent= sqrt(yu_w(i,1)*yu_w(i,1)+yv_w(i,1)*yv_w(i,1))
           ycdragm_w(i) = ust*ust/(1.+vent)/vent
-          print *,'ycdragm_w ust yu yv apres=',ycdragm_w(i),ust,yu_w(i,1),yv_w(i,1)
+!         print *,'ycdragm_w ust yu yv apres=',ycdragm_w(i),ust,yu_w(i,1),yv_w(i,1)
          ENDDO
         ENDIF
@@ -1747,6 +1747,7 @@
                y_flux_u1, y_flux_v1 )
                
-! Special DICE MPL 05082013
+! Special DICE MPL 05082013 puis BOMEX
        IF (ok_prescr_ust) THEN
+          do j=1,knon
 !         ysnow(:)=0.
 !         yqsol(:)=0.
@@ -1761,6 +1762,7 @@
 !         y_dflux_t(:)=0.
 !         y_dflux_q(:)=0.
-          y_flux_u1(:)=ycdragm(:)*(1.+sqrt(yu(:,1)*yu(:,1)+yv(:,1)*yv(:,1)))*yu(:,1)*ypplay(:,1)/RD/yt(:,1)
-          y_flux_v1(:)=ycdragm(:)*(1.+sqrt(yu(:,1)*yu(:,1)+yv(:,1)*yv(:,1)))*yv(:,1)*ypplay(:,1)/RD/yt(:,1)
+          y_flux_u1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yu(j,1)*ypplay(j,1)/RD/yt(j,1)
+          y_flux_v1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yv(j,1)*ypplay(j,1)/RD/yt(j,1)
+          enddo
       ENDIF
 
@@ -1797,4 +1799,11 @@
           END DO
           ! Martin
+! Special DICE MPL 05082013 puis BOMEX MPL 20150410
+       IF (ok_prescr_ust) THEN
+          DO j=1,knon
+          y_flux_u1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yu(j,1)*ypplay(j,1)/RD/yt(j,1)
+          y_flux_v1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yv(j,1)*ypplay(j,1)/RD/yt(j,1)
+          ENDDO
+      ENDIF
           
        CASE(is_oce)
@@ -1825,4 +1834,11 @@
           print *,'arg de surf_ocean: ytsurf_new ',ytsurf_new
        ENDIF
+! Special DICE MPL 05082013 puis BOMEX MPL 20150410
+       IF (ok_prescr_ust) THEN
+          DO j=1,knon
+          y_flux_u1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yu(j,1)*ypplay(j,1)/RD/yt(j,1)
+          y_flux_v1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yv(j,1)*ypplay(j,1)/RD/yt(j,1)
+          ENDDO
+      ENDIF
           
        CASE(is_sic)
@@ -1844,4 +1860,11 @@
                y_flux_u1, y_flux_v1)
           
+! Special DICE MPL 05082013 puis BOMEX MPL 20150410
+       IF (ok_prescr_ust) THEN
+          DO j=1,knon
+          y_flux_u1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yu(j,1)*ypplay(j,1)/RD/yt(j,1)
+          y_flux_v1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yv(j,1)*ypplay(j,1)/RD/yt(j,1)
+          ENDDO
+      ENDIF
 
        CASE DEFAULT
@@ -1886,6 +1909,8 @@
 !
           IF (iflag_split .eq.0) THEN
-             Kech_h(:) = ycdragh(:) * (1.0+SQRT(yu(:,1)**2+yv(:,1)**2)) * &
-                  ypplay(:,1)/(RD*yt(:,1))
+             do j=1,knon
+             Kech_h(j) = ycdragh(j) * (1.0+SQRT(yu(j,1)**2+yv(j,1)**2)) * &
+                  ypplay(j,1)/(RD*yt(j,1))
+             enddo
           ENDIF ! (iflag_split .eq.0)
 
@@ -1895,9 +1920,13 @@
           ENDDO
 
-          y_d_ts(:) = ytsurf_new(:) - yts(:) 
+          do j=1,knon
+          y_d_ts(j) = ytsurf_new(j) - yts(j) 
+          enddo
 
         ELSE ! (ok_flux_surf)
-          y_flux_t1(:) =  yfluxsens(:)
-          y_flux_q1(:) = -yevap(:)
+          do j=1,knon
+          y_flux_t1(j) =  yfluxsens(j)
+          y_flux_q1(j) = -yevap(j)
+          enddo
         ENDIF
 
