Index: WRFV3/phys/module_lmdz_phys.F
===================================================================
--- WRFV3/phys/module_lmdz_phys.F	(revision 1)
+++ WRFV3/phys/module_lmdz_phys.F	(revision 6)
@@ -91,5 +91,5 @@
          & wi_start,wi_end,wj_start,wj_end
     INTEGER, INTENT(IN)                                  :: wparfirstscal, wnum3dm
-    INTEGER, INTENT(IN)                                  :: wjulday
+    REAL, INTENT(IN)                                     :: wjulday
     REAL, INTENT(IN)                                     :: wtime_step
     INTEGER, INTENT(IN)                                  :: wrf_qvid, wrf_qcid,      &
@@ -235,5 +235,5 @@
 ! wnum_tiles: number of tiles
 ! wjulday: WRF julian day
-! wgmt: WRF gmt hour
+! wgmt: WRF gmt hour (0, 1)
 ! wtime_step: WRF time-step (in seconds)
 ! wrf_dbg: level of debug from WRF
@@ -375,9 +375,4 @@
       ' wrf_xtime: ',wrf_xtime,' wrf_rst: ', wrf_isrestart,' wrf_lwbdy: ',           &
       wrf_islowbdyin
-    PRINT *,'     '//TRIM(Spt)//': k p ph ___________'
-    DO iz=1,dimz
-      PRINT  *,iz,wrf_perP(ip,iz,jp),wrf_perGeopot(ip,iz,jp)
-    END DO
-    PRINT *,'dimx: ',dimx,' dimy: ',dimy,' dimz: ',dimz,' ddimx: ',ddimx,' ddimy: ',ddimy
 
     phykeys=0.
@@ -390,5 +385,4 @@
     lmdz_lafin = .FALSE.
 
-!    PRINT *,'Lluis: i j ij lon lat lonvals latvals_______'
     DO iy=1,ddimy
       DO ix=1,ddimx
@@ -400,6 +394,4 @@
         LatValues(ddimx*(iy-1)+ix) = wrf_Lat(ix,iy)*Pi/180.
         PsfctendValues(ddimx*(iy-1)+ix) = wrf_psfctend(ix,iy)
-!        PRINT *,ix,iy,(ddimx)*(iy-1)+ix,wrf_Lon(ix,iy),LonValues(ddimx*(iy-1)+ix),   &
-!          wrf_Lat(ix,iy),LatValues(ddimx*(iy-1)+ix)
       END DO
     END DO
@@ -533,8 +525,4 @@
           wrf_grid%loce,wrf_grid%lsic)
 
-        PRINT *,'    Lluis ltksoil. '// Spt //' ter: ',wrf_grid%lter(ip,jp),' lic: ',&
-          wrf_grid%llic(ip,jp),' loce: ', wrf_grid%loce(ip,jp),' sic: ',             &
-          wrf_grid%lsic(ip,jp)
-
 ! LANDMASK in LMDZ as fractions
         wrf_grid%lmsk = wrf_grid%lter + wrf_grid%llic
@@ -565,6 +553,4 @@
 ! Value is given at the first time-step, not any more, since it will be continuously 
 !   used by LMDZ
-         PRINT *,'  Lluis before '//TRIM(Spt)//' SMOIS: ',wrf_grid%smois(ip,:,jp),   &
-           ' dzs: ',wrf_grid%dzs,' lwsol: ',wrf_grid%lwsol(ip,jp)
 ! l. Fita, LMD. January 2014
 ! This should be in this way, but, because we can not use any scheme DZS=0., so 
@@ -574,6 +560,4 @@
 !            wrf_grid%soil_layers(:,iz,:)/100.
 !        END DO
-        PRINT *,'  Lluis after: '//TRIM(Spt)//' SMOIS: ',wrf_grid%smois(ip,:,jp),    &
-           ' lwsol: ',wrf_grid%lwsol(ip,jp)      
         DO iz=1,4
           wrf_grid%ltksoil(:,iz,:)=wrf_grid%st(:,1,:)
@@ -626,9 +610,4 @@
         ORCHIDEElevels = (/ 0.098, 0.1955, 0.5865, 1.3685, 2.9326, 6.0606, 12.3167,  &
           24.8289, 49.8534, 99.9022, 200. /)
-        PRINT *,'  Lluis interpolating:',SIZE(wrf_grid%soil_layers(ip,:,jp)),' : ',  &
-          ' soil: ',                                                                 &
-          wrf_grid%soil_layers(ip,SIZE(wrf_grid%soil_layers(ip,:,jp)):1:-1,jp),      &
-          ' st: ',wrf_grid%st(ip,:,jp), 11, ' | ', ORCHIDEElevels,' soiltop: ',      &
-          wrf_grid%soilctop(ip,:,jp),' soilbot: ',wrf_grid%soilcbot(ip,:,jp)
 
         datasetname = soilcat
@@ -671,9 +650,4 @@
 !L. Fita, LMD. January 2014 Soil layers are in cm!
 
-            IF ( (ix == ip) .AND. (iy == jp)) PRINT *,' Lluis '//TRIM(Spt)//         &
-              ' soiltop: ', wrf_grid%soilctop(ix,:,iy),' soilbot: ',                 &
-              wrf_grid%soilcbot(ix,:,iy),' layers: ',                                &
-              wrf_grid%soil_layers(ix,Nsoillayers:1:-1,iy)/100.,' sm: ',             &
-              wrf_grid%sm(ix,:,iy)
 
 ! NOTE: wrf_dx, wrf_dy should be matrix !!!
@@ -684,9 +658,4 @@
               soilTmass, wrf_grid%lwsol(ix,iy))
 
-            IF ( (ix == ip) .AND. (iy == jp)) PRINT *,'TOTAL soil mass at '//        &
-              TRIM(Spt)//': ', soilTmass,' [kg] water content: ',                    &
-              wrf_grid%lwsol(ix,iy), ' [kg] equiv. water depth: ',                   &
-              wrf_grid%lwsol(ix,iy)/(wrf_dx*wrf_dy*1000.)
-
             wrf_grid%lwsol(ix,iy) =  wrf_grid%lwsol(ix,iy)/(wrf_dx*wrf_dy*1000.)
 
@@ -694,7 +663,4 @@
               wrf_grid%soil_layers(ix,SIZE(wrf_grid%soil_layers(ix,:,iy)):1:-1,iy),  &
               wrf_grid%st(ix,:,iy),11, ORCHIDEElevels,wrf_grid%lotter(ix,:,iy))
-            IF ( (ix == ip) .AND. (iy == jp)) PRINT *,' Lluis orchidee temperature'//&
-              ' at '// TRIM(Spt)//': ',wrf_grid%lotter(ix,:,iy)
-!            IF ( (ix == ip) .AND. (iy == jp)) STOP
           END DO
         END DO
@@ -852,10 +818,4 @@
         wrf_grid%lzrugsrel = wrf_grid%lrugKsoil(:,1,:)
 
-        PRINT *,'Lluis Topography values at',ip, jp, ' zmea: ',                      &
-          wrf_grid%lzmea(ip,jp),' zstd: ',wrf_grid%lzstd(ip,jp),' zsig: ',           &
-          wrf_grid%lzsig(ip,jp),' zgam: ',wrf_grid%lzgam(ip,jp),' zthe: ',           &
-          wrf_grid%lzthe(ip,jp),' zpic: ',wrf_grid%lzpic(ip,jp),' zval: ',           &
-          wrf_grid%lzval(ip,jp),' rugoro: ',wrf_grid%lzrugsrel(ip,jp)
-
 ! R = 8.3144621 ([J mol/K] http://en.wikipedia.org/wiki/Gas_constant)
 ! mu = 28.97 (Molecular Mass of Air [kg/kmol] from http://www.engineeringtoolbox.com/molecular-mass-air-d_679.html)
@@ -867,5 +827,4 @@
 
 !  Grids do not start/end at the same index
-!!      PRINT *,'  Lluis: Creating area/cu/cv...'
       wrf_area = 0.
       wrf_cu = 0.
@@ -899,10 +858,4 @@
       CALL lmdz_vars_init(ddimxy) 
 
-!      PRINT *,'  Lluis before loading 5 11 ltksoil: ',wrf_grid%ltksoil(5,:,11),      &
-!        ' ftsol: ',ftsol(315,:)
-!      PRINT *,'  Lluis before loading 5 11 qsurf: ',qsurf(315,:),' qsurf_rst: ',     &
-!        qsurf_rst(315,:),' smois :',wrf_grid%smois(5,:,11),' lqksoil: ',             &
-!        wrf_grid%lqksoil(5,:,11)
-
       CALL load_lmdz(wrf_grid,ddimx,ddimy,dimz,wbdyw,ddimxy,4,11,rlon,rlat,zmasq,    &
         pctsrf,ftsol,ftsoil_rst,qsurf_rst,qsol_rst,falb1,falb2,evap_rst,snow_rst,    &
@@ -916,26 +869,4 @@
         rugos_rst, agesno_rst, ftsoil_rst)
       CALL fonte_neige_init(restart_runoff)
-
-      PRINT *,'  Lluis on module restart_runoff: ',UBOUND(restart_runoff)
-      PRINT *,'    ffonte_global: ',ALLOCATED(ffonte_global)
-      PRINT *,'    fqfonte_global: ',ALLOCATED(fqfonte_global)
-      PRINT *,'    fqcalving_global: ',ALLOCATED(fqcalving_global)
-
-      PRINT *,'  Lluis after loading '//TRIM(Spt)//':  qsurf_rst: ', qsurf_rst(lp,:),&
-        ' smois :',wrf_grid%smois(ip,:,jp),' lqksoil: ', wrf_grid%lqksoil(ip,:,jp)
-
-      PRINT *,'  Lluis after loading '//TRIM(Spt)//' ltksoil: ',                     &
-        wrf_grid%ltksoil(ip,:,jp),' ftsol: ',ftsol(lp,:)
-
-      PRINT *,'  Lluis after loading '//TRIM(Spt)//' lwsol: ',wrf_grid%lwsol(ip,jp), &
-        ' qsol_rst: ',qsol_rst(lp)
-
-      PRINT *,'  Lluis lmdz pt= ',lp,' k qsurf_rst qsol_rst evap_rst snow_rst ' //   &
-        'fder_rst agesno_rst rugos_rst'
-      DO iz=1,4
-        PRINT *,'    * ',iz,qsurf_rst(lp,iz),qsol_rst(lp), evap_rst(lp,iz),          &
-          snow_rst(lp,iz),fder_rst(lp), agesno_rst(lp,iz),rugos_rst(lp,iz)
-        PRINT *,'      ftsoil: ',ftsoil_rst(lp,:,iz)
-      END DO
 
       DO ix=1,ddimx
@@ -959,19 +890,11 @@
       CALL NOread_limit_init(ddimxy, klon_glo)
 
-      PRINT *,'   Lluis klon_glo,', klon_glo,' klon: ',klon,'lp= ',lp,' sst: ',sst(lp)
       sst_glo = sst(1:ddimxy)
 ! In ocean_forced_mod is used as tsurf_limit = sst(knindex(1:knon))
       tsurf_limit = ftsol(:,3)
-      PRINT *,'   Lluis pt tsurf_limit_check_____'
-      DO ix=-10,10
-        PRINT *,lp+ix,tsurf_limit(lp+ix),sst(lp+ix)
-      END DO
 
 ! In surf_land_bucket_mod is used as z0_new=rugos(knindex(1:knon),1) [z0_new == z0_limit]
       z0_limit = rugos_rst(:,1)
       alb_limit = albedo
-      PRINT *,'  Lluis pre-limit_values ftsol: ',ftsol(lp,:),' pctsrf: ',pctsrf(lp,:)
-      PRINT *,'  Lluis limit_values lp= ',lp,' tsurf_limit: ',tsurf_limit(lp),       &
-        ' z0_limit: ', z0_limit(lp),' alb_limit: ',alb_limit(lp)
 
       rug_glo = rugos_rst(:,is_ter)
@@ -981,5 +904,4 @@
 !      temps=wjulday
 
-!!      PRINT *,'  Lluis: calling iniphysiq'
       CALL iniphysiq(ddimx*ddimy, dimz, 86400., wjulday, wtime_step,                    &
         & LatValues, LonValues, wrf_area, CUValues, CVValues, 1./reradius, g, r_d, cp,  &
@@ -999,19 +921,4 @@
       PRINT *,'  Lluis: END of firststep !!!!!'
     END IF firststep
-
-    ip = 12
-    jp = 7
-    PRINT *,'  ',wrf_grid%itimestep,' Lluis at 12 7: ',wrf_grid%lter(ip,jp), wrf_grid%llic(ip,jp), &
-      wrf_grid%loce(ip,jp), wrf_grid%lsic(ip,jp), ' rugksoil: ',                     &
-      wrf_grid%lrugksoil(ip,1,jp), wrf_grid%lrugksoil(ip,2,jp),                      &
-      wrf_grid%lrugksoil(ip,3,jp), wrf_grid%lrugksoil(ip,4,jp),' rug: ',             &
-      wrf_grid%lrug(ip,jp)
-    ip = 6
-    jp = 9
-    PRINT *,'  ',wrf_grid%itimestep,' Lluis at 6 9: ',wrf_grid%lter(ip,jp), wrf_grid%llic(ip,jp), &
-      wrf_grid%loce(ip,jp), wrf_grid%lsic(ip,jp), ' rugksoil: ',                     &
-      wrf_grid%lrugksoil(ip,1,jp), wrf_grid%lrugksoil(ip,2,jp),                      &
-      wrf_grid%lrugksoil(ip,3,jp), wrf_grid%lrugksoil(ip,4,jp),' rug: ',             &
-      wrf_grid%lrug(ip,jp)
 
     ip=icheck_p
@@ -1067,7 +974,4 @@
       PRINT *,'  ' // subname // ': Problem allocating WRF moisture array!', ierr
     END IF
-
-!!    PRINT *,'  Lluis wrf MOIST:',LBOUND(wrf_MOIST),' x ',UBOUND(wrf_MOIST)
-!!    PRINT *,'Mix ratio: ',LBOUND(MixingRatioValues),' x ', UBOUND(MixingRatioValues)
 
     IF (ALLOCATED(MixingRatiotendValues)) DEALLOCATE(MixingRatiotendValues)
@@ -1230,23 +1134,4 @@
     END IF
 
-    PRINT *,'  Lluis before physiq '//TRIM(Spt)//' ltksoil: ',                       &
-      wrf_grid%ltksoil(ip,:,jp), ' ftsol: ',ftsol(lp,:)
-
-    PRINT *,'  Lluis beofre physiq '//TRIM(Spt)//' qsurf: ',zxqsurf(lp),             &
-        ' qsurf_rst: ',qsurf_rst(lp,:),' smois :',wrf_grid%smois(ip,:,jp),           &
-        ' lqksoil: ', wrf_grid%lqksoil(ip,:,jp)
-
-    PRINT *,'  Lluis before physiq: '//TRIM(Spt)//' FullPressValues HalfPressValues '&
-       // 'GeopotValues UValues VValues TValues: '
-    DO iz=1,dimz
-      PRINT *,FullPressValues(lp,iz), HalfPressValues(lp,iz), GeopotValues(lp,iz),   &
-        UValues(lp,iz), VValues(lp,iz), TValues(lp,iz)
-    END DO
-
-    ip = 23
-    jp = 18
-    PRINT *,'  ',wrf_grid%itimestep,' Lluis before physiq at ',ip,', ', jp,' sst: ', &
-      wrf_grid%sst(ip,jp), sst((jp-1)*ddimx+ip),' ftsol 3: ',ftsol((jp-1)*ddimx+ip,is_oce)
-
     ip=icheck_p
     jp=jcheck_p
@@ -1281,21 +1166,4 @@
        &  lmdz_dudyn, lmdz_PVtheta)
 
-    PRINT *,'  Lluis after physiq '//TRIM(Spt)//' ltksoil: ',                        &
-      wrf_grid%ltksoil(ip,:,jp), ' ftsol: ',ftsol(315,:)
-    PRINT *,'  Lluis after physiq '//TRIM(Spt)//' qsurf: ',zxqsurf(lp),              &
-        ' qsurf_rst: ', qsurf_rst(lp,:),' smois :',wrf_grid%smois(ip,:,jp),          &
-        ' lqksoil: ', wrf_grid%lqksoil(ip,:,jp)
-    PRINT *,'  Lluis after physiq: '//TRIM(Spt)//' FullPressValues HalfPressValues ' &
-       // 'GeopotValuesUValues VValues TValues: '
-    DO iz=1,dimz
-      PRINT *,FullPressValues(lp,iz), HalfPressValues(lp,iz), GeopotValues(lp,iz),   &
-        UValues(lp,iz), VValues(lp,iz), TValues(lp,iz)
-    END DO
-
-    ip = 23
-    jp = 18
-    PRINT *,'  ',wrf_grid%itimestep,' Lluis after physiq at ',ip,', ', jp,' sst: ', &
-      wrf_grid%sst(ip,jp), sst((jp-1)*ddimx+ip),' ftsol 3: ',ftsol((jp-1)*ddimx+ip,is_oce)
-
     ip=icheck_p
     jp=jcheck_p
@@ -1370,9 +1238,4 @@
     wrf_Vtend(wime,:,wjme) = wrfVtenddestagg(wime-1,:,wjme-1)
 
-    PRINT *,'  Lluis ' // TRIM(Spt) // ' wrf_Utend: ',wrf_Utend(ip,1,jp),            &
-      ' wrfUtenddestagg a: ', wrfUtenddestagg(ip-1,1,jp),' wrfUtenddestagg b: ',     &
-      wrfUtenddestagg(ip,1,jp), ' wrfUtenddestagg c: ',wrfUtenddestagg(ip+1,1,jp),   &
-      ' UtendValues: ', UtendValues(lp,1),' lp: ',lp
-    PRINT *,'dim UtendValues: ',1-wbdyw,':',dimx-1-wbdyw+1,', ',dimz+1,', ',1-wbdyw,':',dimy-1-wbdyw+1
 !!    wrf_Ttend = wrf_Temptend
 
@@ -1467,20 +1330,4 @@
     LMDZvarmethod = 'prod'
 
-    PRINT *,'  AFTER physic and before recomputing!!!!!'
-    ip = 12
-    jp = 7
-    PRINT *,'  ',wrf_grid%itimestep,' Lluis at 12 7: ',wrf_grid%lter(ip,jp), wrf_grid%llic(ip,jp), &
-      wrf_grid%loce(ip,jp), wrf_grid%lsic(ip,jp), ' rugksoil: ',                   &
-      wrf_grid%lrugksoil(ip,1,jp), wrf_grid%lrugksoil(ip,2,jp),                    &
-      wrf_grid%lrugksoil(ip,3,jp), wrf_grid%lrugksoil(ip,4,jp),' rug: ',           &
-      wrf_grid%lrug(ip,jp)
-    ip = 6
-    jp = 9
-    PRINT *,'  ',wrf_grid%itimestep,' Lluis at 6 9: ',wrf_grid%lter(ip,jp), wrf_grid%llic(ip,jp), &
-      wrf_grid%loce(ip,jp), wrf_grid%lsic(ip,jp), ' rugksoil: ',                   &
-      wrf_grid%lrugksoil(ip,1,jp), wrf_grid%lrugksoil(ip,2,jp),                    &
-      wrf_grid%lrugksoil(ip,3,jp), wrf_grid%lrugksoil(ip,4,jp),' rug: ',           &
-      wrf_grid%lrug(ip,jp)
-
     CALL wrf_varKsoil(wims,wime,wjms,wjme,dimx,dimy,wbdyw,LMDZvarmethod,             &
       wrf_grid%lter,wrf_grid%llic,wrf_grid%loce,wrf_grid%lsic,                       &
@@ -1488,25 +1335,4 @@
       wrf_grid%lrugksoil(:,4,:),wrf_grid%lrug)
 
-    PRINT *,'  AFTER physic and recomputing!!!!!'
-    ip = 12
-    jp = 7
-    PRINT *,'  ',wrf_grid%itimestep,' Lluis at 12 7: ',wrf_grid%lter(ip,jp), wrf_grid%llic(ip,jp), &
-      wrf_grid%loce(ip,jp), wrf_grid%lsic(ip,jp), ' rugksoil: ',                   &
-      wrf_grid%lrugksoil(ip,1,jp), wrf_grid%lrugksoil(ip,2,jp),                    &
-      wrf_grid%lrugksoil(ip,3,jp), wrf_grid%lrugksoil(ip,4,jp),' rug: ',           &
-      wrf_grid%lrug(ip,jp)
-    ip = 6
-    jp = 9
-    PRINT *,'  ',wrf_grid%itimestep,' Lluis at 6 9: ',wrf_grid%lter(ip,jp), wrf_grid%llic(ip,jp), &
-      wrf_grid%loce(ip,jp), wrf_grid%lsic(ip,jp), ' rugksoil: ',                   &
-      wrf_grid%lrugksoil(ip,1,jp), wrf_grid%lrugksoil(ip,2,jp),                    &
-      wrf_grid%lrugksoil(ip,3,jp), wrf_grid%lrugksoil(ip,4,jp),' rug: ',           &
-      wrf_grid%lrug(ip,jp)
-
-    ip = 23
-    jp = 18
-    PRINT *,'  ',wrf_grid%itimestep,' Lluis at ',ip,', ', jp,' sst : ',              &
-      wrf_grid%sst(ip,jp), sst((jp-1)*ddimx+ip),' ftsol 3: ',ftsol((jp-1)*ddimx+ip,is_oce)
-  
     CALL wrf_varKsoil(wims,wime,wjms,wjme,dimx,dimy,wbdyw,LMDZvarmethod,             &
       wrf_grid%lter,wrf_grid%llic,wrf_grid%loce,wrf_grid%lsic,                       &
@@ -1553,9 +1379,4 @@
 !       & wrf_grid%seed2)
 
-    ip = 23
-    jp = 18
-    PRINT *,'  ',wrf_grid%itimestep,' Lluis after put at ',ip,', ', jp,' sst: ',     &
-      wrf_grid%sst(ip,jp), sst((jp-1)*ddimx+ip),' ftsol 3: ',ftsol((jp-1)*ddimx+ip,is_oce)
-
     ip=icheck_p
     jp=jcheck_p
@@ -1654,5 +1475,4 @@
 !       & wrf_grid%lweakinv, wrf_grid%lwind10m, wrf_grid%lwind10max, wrf_grid%lzmax_th)
 
-    PRINT *, '   Lluis: wrf_grid%la_th', UBOUND(wrf_grid%la_th),' fraca: ',UBOUND(fraca)
 !    DO iz=wkms, wkme-2
 !      CALL get_lmdz_out3D_z_i(wims, wime, wjms, wjme, iz, ddimx, ddimy, dimz, wbdyw, &
@@ -1820,9 +1640,4 @@
 !    END DO
 
-    PRINT *,'  Lluis Loading done! iz ldtlwr _______________'
-    DO iz=1,dimz
-      PRINT *,iz,wrf_grid%ldtlwr(ip,iz,jp)
-    END DO
-
 !    CALL get_lmdz_out3D_ii(wims, wime, wjms, wjme, wkms, wkme, ddimx, ddimy, dimz,   &
 !      & wbdyw, ddimxy, nbsrf, nsoilmx, lmdzmixingratios, 1, 2, wrf_l_thermals,       &
