Index: /branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/SBC/sbcblk_core.F90
===================================================================
--- /branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/SBC/sbcblk_core.F90	(revision 2197)
+++ /branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/SBC/sbcblk_core.F90	(revision 2198)
@@ -12,4 +12,5 @@
    !!            3.0  !  2006-06  (G. Madec) sbc rewritting   
    !!            3.2  !  2009-04  (B. Lemaire)  Introduce iom_put
+   !!            3.3  !  2010-10  (S. Masson)  add diurnal cycle
    !!----------------------------------------------------------------------
 
@@ -26,4 +27,5 @@
    USE fldread         ! read input fields
    USE sbc_oce         ! Surface boundary condition: ocean fields
+   USE sbcdcy          ! surface boundary condition: diurnal cycle
    USE iom             ! I/O manager library
    USE in_out_manager  ! I/O manager
@@ -34,5 +36,4 @@
    USE sbc_ice         ! Surface boundary condition: ice fields
 #endif
-
 
    IMPLICIT NONE
@@ -62,8 +63,8 @@
    REAL(wp), PARAMETER ::   Cice =    1.63e-3     ! transfer coefficient over ice
 
-   !                                !!* Namelist namsbc_core : CORE bulk parameters
-   LOGICAL  ::   ln_2m     = .FALSE.     ! logical flag for height of air temp. and hum
-   LOGICAL  ::   ln_taudif = .FALSE.     ! logical flag to use the "mean of stress module - module of mean stress" data
-   REAL(wp) ::   rn_pfac   = 1.          ! multiplication factor for precipitation
+   !                                  !!* Namelist namsbc_core : CORE bulk parameters
+   LOGICAL  ::   ln_2m     = .FALSE.   ! logical flag for height of air temp. and hum
+   LOGICAL  ::   ln_taudif = .FALSE.   ! logical flag to use the "mean of stress module - module of mean stress" data
+   REAL(wp) ::   rn_pfac   = 1.        ! multiplication factor for precipitation
 
    !! * Substitutions
@@ -71,9 +72,8 @@
 #  include "vectopt_loop_substitute.h90"
    !!----------------------------------------------------------------------
-   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009) 
+   !! NEMO/OPA 3.3 , NEMO-consortium (2010) 
    !! $Id$
    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
    !!----------------------------------------------------------------------
-
 CONTAINS
 
@@ -144,8 +144,15 @@
          sn_tdif = FLD_N( 'taudif'  ,    24     ,  'taudif'  ,  .true.    , .false. ,   'yearly'  , ''       , ''         )
          !
-         REWIND( numnam )                    ! ... read in namlist namsbc_core
+         REWIND( numnam )                          ! read in namlist namsbc_core
          READ  ( numnam, namsbc_core )
-         !
-         ! store namelist information in an array
+         !                                         ! check: do we plan to use ln_dm2dc with non-daily forcing?
+         IF( ln_dm2dc .AND. sn_qsr%nfreqh /= 24 )   & 
+            &   CALL ctl_stop( 'sbc_blk_core: ln_dm2dc can be activated only with daily short-wave forcing' ) 
+         IF( ln_dm2dc .AND. sn_qsr%ln_tint ) THEN
+            CALL ctl_warn( 'sbc_blk_core: ln_dm2dc is taking care of the temporal interpolation of daily qsr',   &
+                 &         '              ==> We force time interpolation = .false. for qsr' )
+            sn_qsr%ln_tint = .false.
+         ENDIF
+         !                                         ! store namelist information in an array
          slf_i(jp_wndi) = sn_wndi   ;   slf_i(jp_wndj) = sn_wndj
          slf_i(jp_qsr ) = sn_qsr    ;   slf_i(jp_qlw ) = sn_qlw
@@ -153,11 +160,9 @@
          slf_i(jp_prec) = sn_prec   ;   slf_i(jp_snow) = sn_snow
          slf_i(jp_tdif) = sn_tdif
-         !
-         ! do we use HF tau information?
-         lhftau = ln_taudif
+         !                 
+         lhftau = ln_taudif                        ! do we use HF tau information?
          jfld = jpfld - COUNT( (/.NOT. lhftau/) )
          !
-         ! set sf structure
-         ALLOCATE( sf(jfld), STAT=ierror )
+         ALLOCATE( sf(jfld), STAT=ierror )         ! set sf structure
          IF( ierror > 0 ) THEN
             CALL ctl_stop( 'sbc_blk_core: unable to allocate sf structure' )   ;   RETURN
@@ -167,20 +172,19 @@
             ALLOCATE( sf(ifpr)%fdta(jpi,jpj,1,2) )
          END DO
-         !
-         ! fill sf with slf_i and control print
-         CALL fld_fill( sf, slf_i, cn_dir, 'sbc_blk_core', 'flux formulattion for ocean surface boundary condition', 'namsbc_core' )
+         !                                         ! fill sf with slf_i and control print
+         CALL fld_fill( sf, slf_i, cn_dir, 'sbc_blk_core', 'flux formulation for ocean surface boundary condition', 'namsbc_core' )
          !
       ENDIF
 
-      CALL fld_read( kt, nn_fsbc, sf )                   ! input fields provided at the current time-step
+                       CALL fld_read( kt, nn_fsbc, sf )        ! input fields provided at the current time-step
+
+      IF( ln_dm2dc )   CALL sbc_dcy ( kt , sf(jp_qsr)%fnow )   ! modify now Qsr to include the diurnal cycle
 
 #if defined key_lim3
-      tatm_ice(:,:) = sf(jp_tair)%fnow(:,:,1)
+      tatm_ice(:,:) = sf(jp_tair)%fnow(:,:,1)                  ! LIM3: make Tair available in sea-ice
 #endif
-
-      IF( MOD( kt - 1, nn_fsbc ) == 0 ) THEN
-          CALL blk_oce_core( sf, sst_m, ssu_m, ssv_m )   ! compute the surface ocean fluxes using CLIO bulk formulea
-      ENDIF
-      !                                                  ! using CORE bulk formulea
+      !                                                        ! surface ocean fluxes computed with CLIO bulk formulea
+      IF( MOD( kt - 1, nn_fsbc ) == 0 )   CALL blk_oce_core( sf, sst_m, ssu_m, ssv_m )
+      !
    END SUBROUTINE sbc_blk_core
    
Index: /branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/SBC/sbcdcy.F90
===================================================================
--- /branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/SBC/sbcdcy.F90	(revision 2198)
+++ /branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/SBC/sbcdcy.F90	(revision 2198)
@@ -0,0 +1,214 @@
+MODULE sbcdcy
+   !!======================================================================
+   !!                    ***  MODULE  sbcdcy  ***
+   !! Ocean forcing:  compute the diurnal cycle
+   !!======================================================================
+   !! History : OPA  !  2005-02  (D. Bernie)  Original code
+   !!   NEMO    2.0  !  2006-02  (S. Masson, G. Madec)  adaptation to NEMO
+   !!           3.1  !  2009-07  (J.M. Molines)  adaptation to v3.1
+   !!----------------------------------------------------------------------
+
+   !!----------------------------------------------------------------------
+   !!  sbc_dcy : compute solar flux at kt from daily mean, taking 
+   !!            diurnal cycle into account
+   !!----------------------------------------------------------------------
+   USE oce              ! ocean dynamics and tracers
+   USE phycst           ! ocean physics
+   USE dom_oce          ! ocean space and time domain
+   USE in_out_manager   ! I/O manager
+
+   IMPLICIT NONE
+   PRIVATE
+   INTEGER                      ::   nday_qsr                    ! day when parameters were computed
+   REAL(wp), DIMENSION(jpi,jpj) ::   raa , rbb  , rcc  , rab     ! parameters used to compute the diurnal cycle
+   REAL(wp), DIMENSION(jpi,jpj) ::   rtmd, rdawn, rdusk, rscal   !     -       -         -           -      -
+   REAL(wp), DIMENSION(jpi,jpj) ::   qsr_daily                   ! to hold daily mean QSR
+  
+   PUBLIC   sbc_dcy     ! routine called by sbc
+
+   !!----------------------------------------------------------------------
+   !! NEMO/OPA 3.3 , NEMO-consortium (2010) 
+   !! $Id$ 
+   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
+   !!----------------------------------------------------------------------
+CONTAINS
+
+      SUBROUTINE sbc_dcy( kt, pqsr )
+      !!----------------------------------------------------------------------
+      !!                  ***  ROUTINE sbc_dcy  ***
+      !!
+      !! ** Purpose : introduce a diurnal cycle of qsr from daily values
+      !!
+      !! ** Method  : see Appendix A of Bernie et al. 2007.
+      !!
+      !! ** Action  : redistribute daily QSR on each time step following the diurnal cycle
+      !!
+      !! reference  : Bernie, DJ, E Guilyardi, G Madec, JM Slingo, and SJ Woolnough, 2007
+      !!              Impact of resolving the diurnal cycle in an ocean--atmosphere GCM. 
+      !!              Part 1: a diurnally forced OGCM. Climate Dynamics 29:6, 575-590.
+      !!----------------------------------------------------------------------
+      INTEGER,                        INTENT(in   ) ::   kt     ! ocean time-step index
+      REAL(wp), DIMENSION(jpi,jpj,1), INTENT(inout) ::   pqsr   ! QSR flux with diurnal cycle
+      !!
+      INTEGER  ::   ji, jj                                      ! dummy loop indices
+      REAL(wp) ::   ztwopi, zinvtwopi, zconvrad 
+      REAL(wp) ::   zlo, zup, zlousd, zupusd
+      REAL(wp) ::   zdsws, zdecrad, ztx
+      REAL(wp) ::   ztmp, ztmp1, ztmp2, ztest
+      !---------------------------statement functions------------------------
+      REAL(wp) ::   fintegral, pt1, pt2, paaa, pbbb, pccc       ! dummy statement function arguments
+      fintegral( pt1, pt2, paaa, pbbb, pccc ) =                         &
+         &   paaa * pt2 + zinvtwopi * pbbb * SIN(pccc + ztwopi * pt2)   &
+         & - paaa * pt1 - zinvtwopi * pbbb * SIN(pccc + ztwopi * pt1)
+      !!---------------------------------------------------------------------
+
+      ! Initialization
+      ! --------------
+      ztwopi    = 2. * rpi
+      zinvtwopi = 1. / ztwopi
+      zconvrad  = ztwopi / 360.
+
+      ! When are we during the day (from 0 to 1)
+      zlo = MOD( rdt / rday * REAL( kt-nit000, wp ), 1.)
+      zup = zlo + rdt / rday
+
+      !                                          
+      IF( kt == nit000 ) THEN       ! first time step only               
+         IF(lwp) THEN
+            WRITE(numout,*)
+            WRITE(numout,*) 'sbc_dcy : introduce diurnal cycle from daily mean qsr'
+            WRITE(numout,*) '~~~~~~~'
+            WRITE(numout,*)
+         ENDIF
+         nday_qsr = 0
+         ! Compute rcc needed to compute the time integral of the diurnal cycle
+         rcc(:,:) = zconvrad * glamt(:,:) - rpi
+         ! time of midday
+         rtmd(:,:) = 0.5 - glamt(:,:) / 360.
+         rtmd(:,:) = MOD( (rtmd(:,:) + 1.), 1. )
+      ENDIF
+
+      ! If this is a new day, we have to update the dawn, dusk and scaling function  
+      !----------------------
+    
+      !     2.1 dawn and dusk  
+
+      ! nday is the number of days since the beginning of the current month 
+      IF( nday_qsr /= nday ) THEN 
+         ! save the day of the year and the daily mean of qsr
+         nday_qsr = nday 
+         ! number of days since the previous winter solstice (supposed to be always 21 December)         
+         zdsws = 11 + nday_year
+         ! declination of the earths orbit
+         zdecrad = (-23.5 * zconvrad) * COS( zdsws * ztwopi / REAL(nyear_len(1),wp) )
+         ! save the daily QSR for nest hours of the day
+         qsr_daily(:,:) = pqsr(:,:,1)
+         ! Compute A and B needed to compute the time integral of the diurnal cycle
+        
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               ztmp = zconvrad * gphit(ji,jj)
+               raa(ji,jj) = SIN( ztmp ) * SIN( zdecrad )
+               rbb(ji,jj) = COS( ztmp ) * COS( zdecrad )
+            END DO  
+         END DO  
+
+         ! Compute the time of dawn and dusk
+
+         ! rab to test if the day time is equal to 0, less than 24h of full day        
+         rab(:,:) = -raa(:,:) / rbb(:,:)
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               IF ( ABS(rab(ji,jj)) < 1 ) THEN         ! day duration is less than 24h
+         ! When is it night?
+                  ztx = zinvtwopi * (ACOS(rab(ji,jj)) - rcc(ji,jj))
+                  ztest = -rbb(ji,jj) * SIN( rcc(ji,jj) + ztwopi * ztx )
+         ! is it dawn or dusk?
+                  IF ( ztest > 0 ) THEN
+                     rdawn(ji,jj) = ztx
+                     rdusk(ji,jj) = rtmd(ji,jj) + ( rtmd(ji,jj) - rdawn(ji,jj) )
+                  ELSE
+                     rdusk(ji,jj) = ztx
+                     rdawn(ji,jj) = rtmd(ji,jj) - ( rdusk(ji,jj) - rtmd(ji,jj) )
+                  ENDIF
+               ELSE
+                  rdawn(ji,jj) = rtmd(ji,jj) + 0.5
+                  rdusk(ji,jj) = rdawn(ji,jj)
+               ENDIF
+             END DO  
+         END DO  
+         rdawn(:,:) = MOD((rdawn(:,:) + 1.), 1.)
+         rdusk(:,:) = MOD((rdusk(:,:) + 1.), 1.)
+
+
+         !     2.2 Compute the scalling function:
+         !         S* = the inverse of the time integral of the diurnal cycle from dawm to dusk
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               IF ( ABS(rab(ji,jj)) < 1 ) THEN         ! day duration is less than 24h
+                  IF ( rdawn(ji,jj) < rdusk(ji,jj) ) THEN      ! day time in one part
+                     rscal(ji,jj) = fintegral(rdawn(ji,jj), rdusk(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
+                     rscal(ji,jj) = 1. / rscal(ji,jj)
+                  ELSE                                         ! day time in two parts
+                     rscal(ji,jj) = fintegral(0., rdusk(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj))   &
+                        &         + fintegral(rdawn(ji,jj), 1., raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
+                     rscal(ji,jj) = 1. / rscal(ji,jj)
+                  ENDIF
+               ELSE
+                  IF ( raa(ji,jj) > rbb(ji,jj) ) THEN         ! 24h day
+                     rscal(ji,jj) = fintegral(0., 1., raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
+                     rscal(ji,jj) = 1. / rscal(ji,jj)
+                  ELSE                                          ! No day
+                     rscal(ji,jj) = 0.e0
+                  ENDIF
+               ENDIF
+            END DO  
+         END DO  
+         !
+         ztmp = rday / rdt
+         rscal(:,:) = rscal(:,:) * ztmp
+
+      ENDIF 
+
+         !     3. update qsr with the diurnal cycle
+         !     ------------------------------------
+
+      DO jj = 1, jpj
+         DO ji = 1, jpi
+            IF( ABS(rab(ji,jj)) < 1 ) THEN         ! day duration is less than 24h
+               !
+               IF( rdawn(ji,jj) < rdusk(ji,jj) ) THEN       ! day time in one part
+                  zlousd = MAX(zlo, rdawn(ji,jj))
+                  zlousd = MIN(zlousd, zup)
+                  zupusd = MIN(zup, rdusk(ji,jj))
+                  zupusd = MAX(zupusd, zlo)
+                  ztmp = fintegral(zlousd, zupusd, raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
+                  pqsr(ji,jj,1) = qsr_daily(ji,jj) * ztmp * rscal(ji,jj)
+                  !
+               ELSE                                         ! day time in two parts
+                  zlousd = MIN(zlo, rdusk(ji,jj))
+                  zupusd = MIN(zup, rdusk(ji,jj))
+                  ztmp1 = fintegral(zlousd, zupusd, raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
+                  zlousd = MAX(zlo, rdawn(ji,jj))
+                  zupusd = MAX(zup, rdawn(ji,jj))
+                  ztmp2 = fintegral(zlousd, zupusd, raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
+                  ztmp = ztmp1 + ztmp2
+                  pqsr(ji,jj,1) = qsr_daily(ji,jj) * ztmp * rscal(ji,jj)
+               ENDIF
+            ELSE                                   ! 24h light or 24h night
+               !
+               IF( raa(ji,jj) > rbb(ji,jj) ) THEN         ! 24h day
+                  ztmp = fintegral(zlo, zup, raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
+                  pqsr(ji,jj,1) = qsr_daily(ji,jj) * ztmp * rscal(ji,jj)
+                  !
+               ELSE                                         ! No day
+                  pqsr(ji,jj,1) = 0.e0
+               ENDIF
+            ENDIF
+         END DO  
+      END DO  
+      !
+   END SUBROUTINE sbc_dcy
+
+   !!======================================================================
+END MODULE sbcdcy
Index: /branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/SBC/sbcflx.F90
===================================================================
--- /branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/SBC/sbcflx.F90	(revision 2197)
+++ /branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/SBC/sbcflx.F90	(revision 2198)
@@ -4,28 +4,16 @@
    !! Ocean forcing:  momentum, heat and freshwater flux formulation
    !!=====================================================================
-   !! History :  9.0   !  06-06  (G. Madec)  Original code
+   !! History :  1.0  !  2006-06  (G. Madec)  Original code
+   !!            3.3  !  2010-10  (S. Masson)  add diurnal cycle
    !!----------------------------------------------------------------------
 
    !!----------------------------------------------------------------------
    !!   namflx   : flux formulation namlist
-   !!   sbc_flx  : flux formulation as ocean surface boundary condition
-   !!              (forced mode, fluxes read in NetCDF files)
-   !!----------------------------------------------------------------------
-   !! question diverses
-   !!  *   ajouter un test sur la division entier de freqh et rdttra ???
-   !!  **  ajoute dans namelist: 1 year forcing files
-   !!                         or forcing file starts at the begining of the run
-   !!  *** we assume that the forcing file start and end with the previous
-   !!      year last record and the next year first record (useful for
-   !!      time interpolation, required even if no time interp???)
-   !!  *   ajouter un test sur la division de la frequence en pas de temps
-   !!  ==> daymod ajout de nsec_year = number of second since the begining of the year
-   !!      assumed to be 0 at 0h january the 1st (i.e. 24h december the 31)
-   !!
-   !!  *** regrouper dtatem et dtasal
+   !!   sbc_flx  : flux formulation as ocean surface boundary condition (forced mode, fluxes read in NetCDF files)
    !!----------------------------------------------------------------------
    USE oce             ! ocean dynamics and tracers
    USE dom_oce         ! ocean space and time domain
-   USE sbc_oce         ! Surface boundary condition: ocean fields
+   USE sbc_oce         ! surface boundary condition: ocean fields
+   USE sbcdcy          ! surface boundary condition: diurnal cycle on qsr
    USE phycst          ! physical constants
    USE fldread         ! read input fields
@@ -52,9 +40,8 @@
 #  include "vectopt_loop_substitute.h90"
    !!----------------------------------------------------------------------
-   !!   OPA 9.0 , LOCEAN-IPSL (2006) 
+   !! NEMO/OPA 3.3 , NEMO-consortium (2010) 
    !! $Id$
    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
    !!----------------------------------------------------------------------
-
 CONTAINS
 
@@ -98,28 +85,30 @@
       NAMELIST/namsbc_flx/ cn_dir, sn_utau, sn_vtau, sn_qtot, sn_qsr, sn_emp
       !!---------------------------------------------------------------------
-      !                                         ! ====================== !
-      IF( kt == nit000 ) THEN                   !  First call kt=nit000  !
-         !                                      ! ====================== !
+      !
+      IF( kt == nit000 ) THEN                ! First call kt=nit000  
          ! set file information
          cn_dir = './'        ! directory in which the model is executed
          ! ... default values (NB: frequency positive => hours, negative => months)
-         !              !   file    ! frequency !  variable  ! time intep !  clim   ! 'yearly' or ! weights  ! rotation   !
-         !              !   name    !  (hours)  !   name     !   (T/F)    !  (T/F)  !  'monthly'  ! filename ! pairs      !
-         sn_utau = FLD_N(   'utau'  ,    24     ,  'utau'    ,  .false.   , .false. ,   'yearly'  , ''       , ''         )
-         sn_vtau = FLD_N(   'vtau'  ,    24     ,  'vtau'    ,  .false.   , .false. ,   'yearly'  , ''       , ''         )
-         sn_qtot = FLD_N(   'qtot'  ,    24     ,  'qtot'    ,  .false.   , .false. ,   'yearly'  , ''       , ''         )
-         sn_qsr  = FLD_N(   'qsr'   ,    24     ,  'qsr'     ,  .false.   , .false. ,   'yearly'  , ''       , ''         )
-         sn_emp  = FLD_N(   'emp'   ,    24     ,  'emp'     ,  .false.   , .false. ,   'yearly'  , ''       , ''         )
-
-         REWIND ( numnam )               ! ... read in namlist namflx
+         !              !  file   ! frequency !  variable  ! time intep !  clim   ! 'yearly' or ! weights  ! rotation  !
+         !              !  name   !  (hours)  !   name     !   (T/F)    !  (T/F)  !  'monthly'  ! filename ! pairs     !
+         sn_utau = FLD_N(  'utau' ,    24     ,  'utau'    ,  .false.   , .false. ,   'yearly'  , ''       , ''        )
+         sn_vtau = FLD_N(  'vtau' ,    24     ,  'vtau'    ,  .false.   , .false. ,   'yearly'  , ''       , ''        )
+         sn_qtot = FLD_N(  'qtot' ,    24     ,  'qtot'    ,  .false.   , .false. ,   'yearly'  , ''       , ''        )
+         sn_qsr  = FLD_N(  'qsr'  ,    24     ,  'qsr'     ,  .false.   , .false. ,   'yearly'  , ''       , ''        )
+         sn_emp  = FLD_N(  'emp'  ,    24     ,  'emp'     ,  .false.   , .false. ,   'yearly'  , ''       , ''        )
+         !
+         REWIND ( numnam )                         ! read in namlist namflx
          READ   ( numnam, namsbc_flx ) 
-
-         ! store namelist information in an array
+         !
+         !                                         ! check: do we plan to use ln_dm2dc with non-daily forcing?
+         IF( ln_dm2dc .AND. sn_qsr%nfreqh /= 24 )   &
+            &   CALL ctl_stop( 'sbc_blk_core: ln_dm2dc can be activated only with daily short-wave forcing' ) 
+         !
+         !                                         ! store namelist information in an array
          slf_i(jp_utau) = sn_utau   ;   slf_i(jp_vtau) = sn_vtau
          slf_i(jp_qtot) = sn_qtot   ;   slf_i(jp_qsr ) = sn_qsr 
          slf_i(jp_emp ) = sn_emp
-
-         ! set sf structure
-         ALLOCATE( sf(jpfld), STAT=ierror )
+         !
+         ALLOCATE( sf(jpfld), STAT=ierror )        ! set sf structure
          IF( ierror > 0 ) THEN   
             CALL ctl_stop( 'sbc_flx: unable to allocate sf structure' )   ;   RETURN  
@@ -129,19 +118,16 @@
             ALLOCATE( sf(ji)%fdta(jpi,jpj,1,2) )
          END DO
-
-
-         ! fill sf with slf_i and control print
+         !                                         ! fill sf with slf_i and control print
          CALL fld_fill( sf, slf_i, cn_dir, 'sbc_flx', 'flux formulation for ocean surface boundary condition', 'namsbc_flx' )
          !
       ENDIF
 
-      CALL fld_read( kt, nn_fsbc, sf )           ! Read input fields and provides the
-      !                                          ! input fields at the current time-step
+                       CALL fld_read( kt, nn_fsbc, sf )       ! input fields provided at the current time-step
+      
+      IF( ln_dm2dc )   CALL sbc_dcy( kt , sf(jp_qsr)%fnow )   ! modify now Qsr to include the diurnal cycle
 
-      IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN
-         !
-         ! set the ocean fluxes from read fields
+      IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN                    ! update ocean fluxes at each SBC frequency
 !CDIR COLLAPSE
-         DO jj = 1, jpj
+         DO jj = 1, jpj                                           ! set the ocean fluxes from read fields
             DO ji = 1, jpi
                utau(ji,jj) = sf(jp_utau)%fnow(ji,jj,1)
@@ -152,7 +138,6 @@
             END DO
          END DO
-         
-         ! module of wind stress and wind speed at T-point
-         zcoef = 1. / ( zrhoa * zcdrag ) 
+         !                                                        ! module of wind stress and wind speed at T-point
+         zcoef = 1. / ( zrhoa * zcdrag )
 !CDIR NOVERRCHK
          DO jj = 2, jpjm1
@@ -168,9 +153,7 @@
          CALL lbc_lnk( taum(:,:), 'T', 1. )   ;   CALL lbc_lnk( wndm(:,:), 'T', 1. )
 
-         ! Initialization of emps (when no ice model)
-         emps(:,:) = emp (:,:) 
+         emps(:,:) = emp (:,:)                                    ! Initialization of emps (needed when no ice model)
                   
-         ! control print (if less than 100 time-step asked)
-         IF( nitend-nit000 <= 100 .AND. lwp ) THEN
+         IF( nitend-nit000 <= 100 .AND. lwp ) THEN                ! control print (if less than 100 time-step asked)
             WRITE(numout,*) 
             WRITE(numout,*) '        read daily momentum, heat and freshwater fluxes OK'
Index: /branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/SBC/sbcmod.F90
===================================================================
--- /branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/SBC/sbcmod.F90	(revision 2197)
+++ /branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/SBC/sbcmod.F90	(revision 2198)
@@ -5,6 +5,7 @@
    !!======================================================================
    !! History :  3.0  !  2006-07  (G. Madec)  Original code
-   !!            3.1  !  2008-08  (S. Masson, E. Maisonnave, G. Madec) coupled interface
+   !!            3.1  !  2008-08  (S. Masson, A. Caubel, E. Maisonnave, G. Madec) coupled interface
    !!            3.3  !  2010-04  (M. Leclair, G. Madec)  Forcing averaged over 2 time steps
+   !!            3.3  !  2010-10  (S. Masson)  add diurnal cycle
    !!----------------------------------------------------------------------
 
@@ -13,29 +14,28 @@
    !!   sbc            : surface ocean momentum, heat and freshwater boundary conditions
    !!----------------------------------------------------------------------
-   USE oce             ! ocean dynamics and tracers
-   USE dom_oce         ! ocean space and time domain
-   USE phycst          ! physical constants
-
-   USE sbc_oce         ! Surface boundary condition: ocean fields
-   USE sbc_ice         ! Surface boundary condition: ice fields
-   USE sbcssm          ! surface boundary condition: sea-surface mean variables
-   USE sbcana          ! surface boundary condition: analytical formulation
-   USE sbcflx          ! surface boundary condition: flux formulation
-   USE sbcblk_clio     ! surface boundary condition: bulk formulation : CLIO
-   USE sbcblk_core     ! surface boundary condition: bulk formulation : CORE
-   USE sbcice_if       ! surface boundary condition: ice-if sea-ice model
-   USE sbcice_lim      ! surface boundary condition: LIM 3.0 sea-ice model
-   USE sbcice_lim_2    ! surface boundary condition: LIM 2.0 sea-ice model
-   USE sbccpl          ! surface boundary condition: coupled florulation
+   USE oce              ! ocean dynamics and tracers
+   USE dom_oce          ! ocean space and time domain
+   USE phycst           ! physical constants
+   USE sbc_oce          ! Surface boundary condition: ocean fields
+   USE sbc_ice          ! Surface boundary condition: ice fields
+   USE sbcssm           ! surface boundary condition: sea-surface mean variables
+   USE sbcana           ! surface boundary condition: analytical formulation
+   USE sbcflx           ! surface boundary condition: flux formulation
+   USE sbcblk_clio      ! surface boundary condition: bulk formulation : CLIO
+   USE sbcblk_core      ! surface boundary condition: bulk formulation : CORE
+   USE sbcice_if        ! surface boundary condition: ice-if sea-ice model
+   USE sbcice_lim       ! surface boundary condition: LIM 3.0 sea-ice model
+   USE sbcice_lim_2     ! surface boundary condition: LIM 2.0 sea-ice model
+   USE sbccpl           ! surface boundary condition: coupled florulation
    USE cpl_oasis3, ONLY:lk_cpl      ! are we in coupled mode?
-   USE sbcssr          ! surface boundary condition: sea surface restoring
-   USE sbcrnf          ! surface boundary condition: runoffs
-   USE sbcfwb          ! surface boundary condition: freshwater budget
-   USE closea          ! closed sea
-
-   USE prtctl          ! Print control                    (prt_ctl routine)
-   USE restart         ! ocean restart
-   USE iom
-   USE in_out_manager  ! I/O manager
+   USE sbcssr           ! surface boundary condition: sea surface restoring
+   USE sbcrnf           ! surface boundary condition: runoffs
+   USE sbcfwb           ! surface boundary condition: freshwater budget
+   USE closea           ! closed sea
+
+   USE prtctl           ! Print control                    (prt_ctl routine)
+   USE restart          ! ocean restart
+   USE iom              ! IOM library
+   USE in_out_manager   ! I/O manager
 
    IMPLICIT NONE
@@ -50,9 +50,8 @@
 #  include "domzgr_substitute.h90"
    !!----------------------------------------------------------------------
-   !! NEMO/OPA 3.3 , LOCEAN-IPSL (2010) 
+   !! NEMO/OPA 3.3 , NEMO-consortium (2010) 
    !! $Id$
    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
    !!----------------------------------------------------------------------
-
 CONTAINS
 
@@ -70,6 +69,6 @@
       INTEGER ::   icpt      ! temporary integer
       !!
-      NAMELIST/namsbc/ nn_fsbc, ln_ana, ln_flx, ln_blk_clio, ln_blk_core, ln_cpl,   &
-         &             nn_ice , ln_dm2dc, ln_rnf, ln_ssr, nn_fwb, nn_ico_cpl
+      NAMELIST/namsbc/ nn_fsbc, ln_ana  , ln_flx, ln_blk_clio, ln_blk_core, ln_cpl    ,   &
+         &             nn_ice , ln_dm2dc, ln_rnf, ln_ssr     , nn_fwb     , nn_ico_cpl
       !!----------------------------------------------------------------------
 
@@ -80,23 +79,18 @@
       ENDIF
 
-      REWIND( numnam )                   ! Read Namelist namsbc
+      REWIND( numnam )           ! Read Namelist namsbc
       READ  ( numnam, namsbc )
 
-      ! overwrite namelist parameter using CPP key information
-!!gm here no overwrite, test all option via namelist change: require more incore memory
-!!gm  IF( lk_sbc_cpl       ) THEN   ;   ln_cpl      = .TRUE.   ;   ELSE   ;   ln_cpl      = .FALSE.   ;   ENDIF
-
-      IF( Agrif_Root() ) THEN
-        IF( lk_lim2 )            nn_ice      = 2
-        IF( lk_lim3 )            nn_ice      = 3
-      ENDIF
-      !
-      IF( cp_cfg == 'gyre' ) THEN
+      !                          ! overwrite namelist parameter using CPP key information
+      IF( Agrif_Root() ) THEN                ! AGRIF zoom
+        IF( lk_lim2 )   nn_ice      = 2
+        IF( lk_lim3 )   nn_ice      = 3
+      ENDIF
+      IF( cp_cfg == 'gyre' ) THEN            ! GYRE configuration
           ln_ana      = .TRUE.   
           nn_ice      =   0
       ENDIF
       
-      ! Control print
-      IF(lwp) THEN
+      IF(lwp) THEN               ! Control print
          WRITE(numout,*) '        Namelist namsbc (partly overwritten with CPP key setting)'
          WRITE(numout,*) '           frequency update of sbc (and ice)             nn_fsbc     = ', nn_fsbc
@@ -117,4 +111,5 @@
       ENDIF
 
+      !                          ! Checks:
       IF( .NOT. ln_rnf ) THEN                      ! no specific treatment in vicinity of river mouths 
          ln_rnf_mouth  = .false.                      
@@ -139,5 +134,11 @@
          &   CALL ctl_stop( 'sea-ice model requires a bulk formulation or coupled configuration' )
       
-      ! Choice of the Surface Boudary Condition (set nsbc)
+      IF( ln_dm2dc .AND. .NOT.( ln_flx .OR. ln_blk_core ) )   &
+         &   CALL ctl_stop( 'diurnal cycle into qsr field from daily values requires a flux or core-bulk formulation' )
+      
+      IF( ln_dm2dc .AND. ( ( NINT(rday) / ( nn_fsbc * NINT(rdt) ) )  < 8 ) )   &
+         &   CALL ctl_warn( 'diurnal cycle for qsr: the sampling of the diurnal cycle is too small...' )
+      
+      !                          ! Choice of the Surface Boudary Condition (set nsbc)
       icpt = 0
       IF( ln_ana          ) THEN   ;   nsbc =  1   ; icpt = icpt + 1   ;   ENDIF       ! analytical      formulation
@@ -148,5 +149,5 @@
       IF( cp_cfg == 'gyre') THEN   ;   nsbc =  0                       ;   ENDIF       ! GYRE analytical formulation
       IF( lk_esopa        )            nsbc = -1                                       ! esopa test, ALL formulations
-
+      !
       IF( icpt /= 1 .AND. .NOT.lk_esopa ) THEN
          WRITE(numout,*)
@@ -229,6 +230,4 @@
 
       !                                            !==  Misc. Options  ==!
-
-!!gm  IF( ln_dm2dc       )   CALL sbc_dcy( kt )                 ! Daily mean qsr distributed over the Diurnal Cycle
       
       SELECT CASE( nn_ice )                                     ! Update heat and freshwater fluxes over sea-ice areas
@@ -293,7 +292,7 @@
          CALL iom_put( "emp-rnf" , emp  - rnf )                   ! upward water flux
          CALL iom_put( "emps-rnf", emps - rnf )                   ! c/d water flux
-         CALL iom_put( "qns+qsr" , qns  + qsr )                   ! total heat flux   (caution if ln_dm2dc=true, to be 
-         CALL iom_put( "qns"     , qns        )                   ! solar heat flux    moved after the call to iom_setkt)
-         CALL iom_put( "qsr"     ,       qsr  )                   ! solar heat flux    moved after the call to iom_setkt)
+         CALL iom_put( "qns+qsr" , qns  + qsr )                   ! total heat flux 
+         CALL iom_put( "qns"     , qns        )                   ! solar heat flux
+         CALL iom_put( "qsr"     ,       qsr  )                   ! solar heat flux
          IF( nn_ice > 0 )   CALL iom_put( "ice_cover", fr_i )   ! ice fraction 
       ENDIF
