Index: /NEMO/branches/UKMO/NEMO_4.0_add_pond_lids_prints/src/OCE/SBC/sbc_ice.F90
===================================================================
--- /NEMO/branches/UKMO/NEMO_4.0_add_pond_lids_prints/src/OCE/SBC/sbc_ice.F90	(revision 12369)
+++ /NEMO/branches/UKMO/NEMO_4.0_add_pond_lids_prints/src/OCE/SBC/sbc_ice.F90	(revision 12370)
@@ -93,5 +93,6 @@
    
    ! already defined in ice.F90 for SI3
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  a_i
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  a_i               ! Sea ice fraction on categories
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  a_i_last_couple   ! Sea ice fraction on categories at the last coupling point
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  h_i, h_s
 
Index: /NEMO/branches/UKMO/NEMO_4.0_add_pond_lids_prints/src/OCE/SBC/sbccpl.F90
===================================================================
--- /NEMO/branches/UKMO/NEMO_4.0_add_pond_lids_prints/src/OCE/SBC/sbccpl.F90	(revision 12369)
+++ /NEMO/branches/UKMO/NEMO_4.0_add_pond_lids_prints/src/OCE/SBC/sbccpl.F90	(revision 12370)
@@ -48,4 +48,8 @@
    USE lib_mpp        ! distribued memory computing library
    USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
+
+#if defined key_oasis3 
+   USE mod_oasis, ONLY : OASIS_Sent, OASIS_ToRest, OASIS_SentOut, OASIS_ToRestOut 
+#endif 
 
    IMPLICIT NONE
@@ -152,5 +156,5 @@
    INTEGER, PARAMETER ::   jps_wlev   = 32   ! water level 
    INTEGER, PARAMETER ::   jps_fice1  = 33   ! first-order ice concentration (for semi-implicit coupling of atmos-ice fluxes)
-   INTEGER, PARAMETER ::   jps_a_p    = 34   ! meltpond area
+   INTEGER, PARAMETER ::   jps_a_p    = 34   ! meltpond area fraction
    INTEGER, PARAMETER ::   jps_ht_p   = 35   ! meltpond thickness
    INTEGER, PARAMETER ::   jps_kice   = 36   ! sea ice effective conductivity
@@ -159,4 +163,12 @@
 
    INTEGER, PARAMETER ::   jpsnd      = 38   ! total number of fields sent 
+
+#if ! defined key_oasis3 
+   ! Dummy variables to enable compilation when oasis3 is not being used 
+   INTEGER                    ::   OASIS_Sent        = -1 
+   INTEGER                    ::   OASIS_SentOut     = -1 
+   INTEGER                    ::   OASIS_ToRest      = -1 
+   INTEGER                    ::   OASIS_ToRestOut   = -1 
+#endif 
 
    !                                  !!** namelist namsbc_cpl **
@@ -184,4 +196,5 @@
    LOGICAL     ::   ln_usecplmask         !  use a coupling mask file to merge data received from several models
                                           !   -> file cplmask.nc with the float variable called cplmask (jpi,jpj,nn_cplmodel)
+   LOGICAL     ::   ln_scale_ice_fluxes   ! Scale sea ice fluxes by the sea ice fractions at the previous coupling point
    TYPE ::   DYNARR     
       REAL(wp), POINTER, DIMENSION(:,:,:) ::   z3   
@@ -256,5 +269,5 @@
          &                  sn_rcv_iceflx, sn_rcv_co2   , nn_cplmodel , ln_usecplmask, sn_rcv_mslp ,   &
          &                  sn_rcv_icb   , sn_rcv_isf   , sn_rcv_wfreq , sn_rcv_tauw, nn_cats_cpl  ,   &
-         &                  sn_rcv_ts_ice
+         &                  sn_rcv_ts_ice, ln_scale_ice_fluxes
 
       !!---------------------------------------------------------------------
@@ -330,4 +343,5 @@
          WRITE(numout,*)'  ln_usecplmask                       = ', ln_usecplmask
          WRITE(numout,*)'  nn_cats_cpl                         = ', nn_cats_cpl
+         WRITE(numout,*)'  ln_scale_ice_fluxes                 = ', ln_scale_ice_fluxes
       ENDIF
 
@@ -815,4 +829,8 @@
       END SELECT
 
+      ! Initialise ice fractions from last coupling time to zero
+       a_i_last_couple(:,:,:) = 0._wp
+
+
       !                                                      ! ------------------------- ! 
       !                                                      !      Ice Meltponds        ! 
@@ -1643,4 +1661,5 @@
       REAL(wp), DIMENSION(jpi,jpj)     ::   zemp_tot, zemp_ice, zemp_oce, ztprecip, zsprecip  , zevap_oce, zdevap_ice
       REAL(wp), DIMENSION(jpi,jpj)     ::   zqns_tot, zqns_oce, zqsr_tot, zqsr_oce, zqprec_ice, zqemp_oce, zqemp_ice
+      REAL(wp), DIMENSION(jpi,jpj)     ::   zevap_ice_total, ztmp1
       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zqns_ice, zqsr_ice, zdqns_ice, zqevap_ice, zevap_ice    !!gm , zfrqsr_tr_i
       !!----------------------------------------------------------------------
@@ -1663,10 +1682,9 @@
          ztprecip(:,:) =   frcv(jpr_rain)%z3(:,:,1) + zsprecip(:,:)  ! May need to ensure positive here
          zemp_tot(:,:) =   frcv(jpr_tevp)%z3(:,:,1) - ztprecip(:,:)
-         zemp_ice(:,:) = ( frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_snow)%z3(:,:,1) ) * picefr(:,:)
       CASE( 'oce and ice'   )   ! received fields: jpr_sbpr, jpr_semp, jpr_oemp, jpr_ievp
          zemp_tot(:,:) = ziceld(:,:) * frcv(jpr_oemp)%z3(:,:,1) + picefr(:,:) * frcv(jpr_sbpr)%z3(:,:,1)
          zemp_ice(:,:) = frcv(jpr_semp)%z3(:,:,1) * picefr(:,:)
          zsprecip(:,:) = frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_semp)%z3(:,:,1)
-         ztprecip(:,:) = frcv(jpr_semp)%z3(:,:,1) - frcv(jpr_sbpr)%z3(:,:,1) + zsprecip(:,:)
+         ztprecip(:,:) = frcv(jpr_semp)%z3(:,:,1) - frcv(jpr_sbpr)%z3(:,:,1) + zsprecip(:,:)         
       CASE( 'none'      )       ! Not available as for now: needs additional coding below when computing zevap_oce 
       !                         ! since fields received are not defined with none option
@@ -1675,4 +1693,26 @@
 
 #if defined key_si3
+
+      ! --- evaporation over ice (kg/m2/s) --- !
+      zevap_ice_total(:,:) = 0._wp   
+      IF (sn_rcv_emp%clcat == 'yes') THEN
+         DO jl=1,jpl
+            IF (ln_scale_ice_fluxes) THEN
+               zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,jl) * a_i_last_couple(:,:,jl)
+            ELSE
+               zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,jl)
+            ENDIF
+            zevap_ice_total(:,:) = zevap_ice_total(:,:) + zevap_ice(:,:,jl)
+         ENDDO
+      ELSE
+         zevap_ice(:,:,1) = frcv(jpr_ievp)%z3(:,:,1 )
+         zevap_ice_total(:,:) = zevap_ice(:,:,1)
+      ENDIF
+
+      IF ( TRIM( sn_rcv_emp%cldes ) == 'conservative' ) THEN
+         ! For conservative case zemp_ice has not been defined yet. Do it now.
+         zemp_ice(:,:) = zevap_ice_total(:,:) - frcv(jpr_snow)%z3(:,:,1) * picefr(:,:)
+      END IF
+
       ! zsnw = snow fraction over ice after wind blowing (=picefr if no blowing)
       zsnw(:,:) = 0._wp   ;   CALL ice_thd_snwblow( ziceld, zsnw )
@@ -1683,11 +1723,8 @@
 
       ! --- evaporation over ocean (used later for qemp) --- !
-      zevap_oce(:,:) = frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:)
-
-      ! --- evaporation over ice (kg/m2/s) --- !
-      DO jl=1,jpl
-         IF (sn_rcv_emp%clcat == 'yes') THEN   ;   zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,jl)
-         ELSE                                  ;   zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,1 )   ;   ENDIF
-      ENDDO
+      zevap_oce(:,:) = frcv(jpr_tevp)%z3(:,:,1) - zevap_ice_total(:,:)
+
+      
+
 
       ! since the sensitivity of evap to temperature (devap/dT) is not prescribed by the atmosphere, we set it to 0
@@ -1727,5 +1764,11 @@
          sprecip (:,:)   = zsprecip (:,:)
          tprecip (:,:)   = ztprecip (:,:)
-         evap_ice(:,:,:) = zevap_ice(:,:,:)
+         IF (ln_scale_ice_fluxes) THEN
+            ! Convert from grid box means to sea ice means
+            WHERE( a_i(:,:,:) > 0.0_wp ) evap_ice(:,:,:) = zevap_ice(:,:,:) / a_i(:,:,:)
+            WHERE( a_i(:,:,:) <= 0.0_wp ) evap_ice(:,:,:) = 0.0
+         ELSE
+            evap_ice(:,:,:) = zevap_ice(:,:,:)
+         ENDIF
          DO jl = 1, jpl
             devap_ice(:,:,jl) = zdevap_ice(:,:)
@@ -1774,7 +1817,7 @@
       IF( iom_use('snow_ao_cea') )  CALL iom_put( 'snow_ao_cea' , sprecip(:,:) * ( 1._wp - zsnw(:,:) )                  )  ! Snow over ice-free ocean  (cell average)
       IF( iom_use('snow_ai_cea') )  CALL iom_put( 'snow_ai_cea' , sprecip(:,:) *           zsnw(:,:)                    )  ! Snow over sea-ice         (cell average)
-      IF( iom_use('subl_ai_cea') )  CALL iom_put( 'subl_ai_cea' , frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) * tmask(:,:,1) )  ! Sublimation over sea-ice (cell average)
+      IF( iom_use('subl_ai_cea') )  CALL iom_put( 'subl_ai_cea' , zevap_ice_total(:,:) * tmask(:,:,1) )  ! Sublimation over sea-ice (cell average)
       IF( iom_use('evap_ao_cea') )  CALL iom_put( 'evap_ao_cea' , ( frcv(jpr_tevp)%z3(:,:,1)  &
-         &                                                        - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) ) * tmask(:,:,1) )  ! ice-free oce evap (cell average)
+         &                                                        - zevap_ice_total(:,:) ) * tmask(:,:,1) )  ! ice-free oce evap (cell average)
       ! note: runoff output is done in sbcrnf (which includes icebergs too) and iceshelf output is done in sbcisf
       !
@@ -1784,4 +1827,7 @@
       CASE( 'oce only' )         ! the required field is directly provided
          zqns_tot(:,:) = frcv(jpr_qnsoce)%z3(:,:,1)
+         ! For Met Office sea ice non-solar fluxes are already delt with by JULES so setting to zero
+         ! here so the only flux is the ocean only one.
+         zqns_ice(:,:,:) = 0._wp 
       CASE( 'conservative' )     ! the required fields are directly provided
          zqns_tot(:,:) = frcv(jpr_qnsmix)%z3(:,:,1)
@@ -1847,6 +1893,4 @@
          &             +   zsprecip(:,:)                   * ( 1._wp - zsnw ) * ( zcptsnw (:,:) - rLfus )  ! solid precip over ocean + snow melting
       zqemp_ice(:,:) =     zsprecip(:,:)                   * zsnw             * ( zcptsnw (:,:) - rLfus )  ! solid precip over ice (qevap_ice=0 since atm. does not take it into account)
-!!    zqemp_ice(:,:) = -   frcv(jpr_ievp)%z3(:,:,1)        * picefr(:,:)      *   zcptsnw (:,:)   &        ! ice evap
-!!       &             +   zsprecip(:,:)                   * zsnw             * zqprec_ice(:,:) * r1_rhos  ! solid precip over ice
       
       ! --- total non solar flux (including evap/precip) --- !
@@ -1900,6 +1944,6 @@
       IF ( srcv(jpr_icb)%laction       ) CALL iom_put('hflx_icb_cea'    , - frcv(jpr_icb)%z3(:,:,1) * rLfus )                      ! latent heat from icebergs melting
       IF ( iom_use('hflx_rain_cea')    ) CALL iom_put('hflx_rain_cea'   , ( tprecip(:,:) - sprecip(:,:) ) * zcptrain(:,:) )        ! heat flux from rain (cell average)
-      IF ( iom_use('hflx_evap_cea')    ) CALL iom_put('hflx_evap_cea'   , ( frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1) &
-           &                                                              * picefr(:,:) ) * zcptn(:,:) * tmask(:,:,1) )            ! heat flux from evap (cell average)
+      IF ( iom_use('hflx_evap_cea')    ) CALL iom_put('hflx_evap_cea'   , ( frcv(jpr_tevp)%z3(:,:,1) - zevap_ice_total(:,:)  )  &
+                                                                          * zcptn(:,:) * tmask(:,:,1) )            ! heat flux from evap (cell average)
       IF ( iom_use('hflx_snow_cea')    ) CALL iom_put('hflx_snow_cea'   , sprecip(:,:) * ( zcptsnw(:,:) - rLfus )  )               ! heat flux from snow (cell average)
       IF ( iom_use('hflx_snow_ao_cea') ) CALL iom_put('hflx_snow_ao_cea', sprecip(:,:) * ( zcptsnw(:,:) - rLfus ) &
@@ -1914,4 +1958,7 @@
       CASE( 'oce only' )
          zqsr_tot(:,:  ) = MAX( 0._wp , frcv(jpr_qsroce)%z3(:,:,1) )
+         ! For Met Office sea ice solar fluxes are already delt with by JULES so setting to zero
+         ! here so the only flux is the ocean only one.
+         zqsr_ice(:,:,:) = 0._wp
       CASE( 'conservative' )
          zqsr_tot(:,:  ) = frcv(jpr_qsrmix)%z3(:,:,1)
@@ -1992,4 +2039,6 @@
             ENDDO
          ENDIF
+      CASE( 'none'      ) 
+         zdqns_ice(:,:,:) = 0._wp
       END SELECT
       
@@ -2007,6 +2056,13 @@
       !                                                      ! ========================= !
       CASE ('coupled')
-         qml_ice(:,:,:) = frcv(jpr_topm)%z3(:,:,:)
-         qcn_ice(:,:,:) = frcv(jpr_botm)%z3(:,:,:)
+         IF (ln_scale_ice_fluxes) THEN
+            WHERE( a_i(:,:,:) > 0.0_wp ) qml_ice(:,:,:) = frcv(jpr_topm)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:)
+            WHERE( a_i(:,:,:) <= 0.0_wp ) qml_ice(:,:,:) = 0.0_wp
+            WHERE( a_i(:,:,:) > 0.0_wp ) qcn_ice(:,:,:) = frcv(jpr_botm)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:)
+            WHERE( a_i(:,:,:) <= 0.0_wp ) qcn_ice(:,:,:) = 0.0_wp
+         ELSE
+            qml_ice(:,:,:) = frcv(jpr_topm)%z3(:,:,:)
+            qcn_ice(:,:,:) = frcv(jpr_botm)%z3(:,:,:)
+         ENDIF
       END SELECT
       !
@@ -2098,5 +2154,5 @@
                SELECT CASE( sn_snd_temp%clcat ) 
                CASE( 'yes' )    
-                  ztmp3(:,:,1:jpl) = tn_ice(:,:,1:jpl) * a_i(:,:,1:jpl) 
+                  ztmp3(:,:,1:jpl) = t_i(:,:,1,1:jpl) * a_i(:,:,1:jpl) 
                CASE( 'no' ) 
                   ztmp3(:,:,:) = 0.0 
@@ -2192,4 +2248,14 @@
          IF( ssnd(jps_fice)%laction )   CALL cpl_snd( jps_fice, isec, ztmp3, info )
       ENDIF
+
+      ! If this coupling was successful then save ice fraction for use between coupling points. 
+      ! This is needed for some calculations where the ice fraction at the last coupling point 
+      ! is needed. 
+      IF( info == OASIS_Sent     .OR. info == OASIS_ToRest .OR.   & 
+                     & info == OASIS_SentOut  .OR. info == OASIS_ToRestOut ) THEN 
+         IF ( sn_snd_thick%clcat == 'yes' ) THEN 
+           a_i_last_couple(:,:,1:jpl) = a_i(:,:,1:jpl)
+         ENDIF
+      ENDIF    
 
       IF( ssnd(jps_fice1)%laction ) THEN
@@ -2250,5 +2316,5 @@
       !                                                      !      Ice melt ponds       ! 
       !                                                      ! ------------------------- !
-      ! needed by Met Office
+      ! needed by Met Office - 1) fraction of ponded ice; 2) local/actual pond depth
       IF( ssnd(jps_a_p)%laction .OR. ssnd(jps_ht_p)%laction ) THEN 
          SELECT CASE( sn_snd_mpnd%cldes)  
@@ -2256,12 +2322,12 @@
             SELECT CASE( sn_snd_mpnd%clcat )  
             CASE( 'yes' )  
-               ztmp3(:,:,1:jpl) =  a_ip(:,:,1:jpl)
-               ztmp4(:,:,1:jpl) =  v_ip(:,:,1:jpl)  
+               ztmp3(:,:,1:jpl) =  a_ip_frac(:,:,1:jpl)
+               ztmp4(:,:,1:jpl) =  h_ip(:,:,1:jpl)
             CASE( 'no' )  
                ztmp3(:,:,:) = 0.0  
                ztmp4(:,:,:) = 0.0  
                DO jl=1,jpl  
-                 ztmp3(:,:,1) = ztmp3(:,:,1) + a_ip(:,:,jpl)  
-                 ztmp4(:,:,1) = ztmp4(:,:,1) + v_ip(:,:,jpl) 
+                 ztmp3(:,:,1) = ztmp3(:,:,1) + a_ip_frac(:,:,jpl)
+                 ztmp4(:,:,1) = ztmp4(:,:,1) + h_ip(:,:,jpl)
                ENDDO  
             CASE default   ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_mpnd%clcat' )  
