Index: branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/DOM/phycst.F90
===================================================================
--- branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/DOM/phycst.F90	(revision 5662)
+++ branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/DOM/phycst.F90	(revision 5663)
@@ -68,5 +68,9 @@
    REAL(wp), PUBLIC ::   rcdsn    =    0.31_wp       !: thermal conductivity of snow
    REAL(wp), PUBLIC ::   cpic     = 2067.0_wp        !: specific heat for ice 
+#if defined key_cice
+   REAL(wp), PUBLIC ::   lsub     =    2.835e+6_wp   !: pure ice latent heat of sublimation                   [J/kg]
+#else
    REAL(wp), PUBLIC ::   lsub     =    2.834e+6_wp   !: pure ice latent heat of sublimation                   [J/kg]
+#endif
    REAL(wp), PUBLIC ::   lfus     =    0.334e+6_wp   !: latent heat of fusion of fresh ice                    [J/kg]
    REAL(wp), PUBLIC ::   tmut     =    0.054_wp      !: decrease of seawater meltpoint with salinity
Index: branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_ice.F90
===================================================================
--- branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_ice.F90	(revision 5662)
+++ branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_ice.F90	(revision 5663)
@@ -101,8 +101,12 @@
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   fr_iu              !: ice fraction at NEMO U point
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   fr_iv              !: ice fraction at NEMO V point
-   
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   sstfrz             !: sea surface freezing temperature
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   tsfc_ice           !: sea-ice surface skin temperature (on categories)
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   kn_ice             !: sea-ice surface layer thermal conductivity (on cats)
+
    ! variables used in the coupled interface
    INTEGER , PUBLIC, PARAMETER ::   jpl = ncat
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   u_ice, v_ice          ! jpi, jpj
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  a_p, ht_p ! Meltpond fraction and depth
 #endif
    
@@ -152,14 +156,17 @@
 
 #if defined key_cice
-      ALLOCATE( qla_ice(jpi,jpj,1)    , qlw_ice(jpi,jpj,1)    , qsr_ice(jpi,jpj,1)    , &
+      ALLOCATE( qla_ice(jpi,jpj,ncat) , qlw_ice(jpi,jpj,1)    , qsr_ice(jpi,jpj,1)    , &
                 wndi_ice(jpi,jpj)     , tatm_ice(jpi,jpj)     , qatm_ice(jpi,jpj)     , &
                 wndj_ice(jpi,jpj)     , nfrzmlt(jpi,jpj)      , ss_iou(jpi,jpj)       , &
                 ss_iov(jpi,jpj)       , fr_iu(jpi,jpj)        , fr_iv(jpi,jpj)        , &
                 a_i(jpi,jpj,ncat)     , topmelt(jpi,jpj,ncat) , botmelt(jpi,jpj,ncat) , &
-                STAT= ierr(1) )
-      IF( ln_cpl )   ALLOCATE( u_ice(jpi,jpj)        , fr1_i0(jpi,jpj)       , tn_ice (jpi,jpj,1)    , &
+                sstfrz(jpi,jpj)       , STAT= ierr(1) )
+   ! Alex West: Allocating tn_ice with 5 categories.  When NEMO is used with CICE, this variable
+   ! represents top layer ice temperature, which is multi-category.
+      IF( ln_cpl )   ALLOCATE( u_ice(jpi,jpj)        , fr1_i0(jpi,jpj)       , tn_ice (jpi,jpj,jpl)  , &
          &                     v_ice(jpi,jpj)        , fr2_i0(jpi,jpj)       , alb_ice(jpi,jpj,1)    , &
          &                     emp_ice(jpi,jpj)      , qns_ice(jpi,jpj,1)    , dqns_ice(jpi,jpj,1)   , &
-         &                     STAT= ierr(2) )
+         &                     a_p(jpi,jpj,jpl)      , ht_p(jpi,jpj,jpl)     , tsfc_ice(jpi,jpj,jpl) , &
+         &                     kn_ice(jpi,jpj,jpl) ,    STAT=ierr(2) )
       
 #endif
Index: branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90
===================================================================
--- branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90	(revision 5662)
+++ branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90	(revision 5663)
@@ -46,7 +46,4 @@
    USE p4zflx, ONLY : oce_co2
 #endif
-#if defined key_cice
-   USE ice_domain_size, only: ncat
-#endif
 #if defined key_lim3
    USE limthd_dh       ! for CALL lim_thd_snwblow
@@ -105,5 +102,6 @@
    INTEGER, PARAMETER ::   jpr_e3t1st = 41            ! first T level thickness 
    INTEGER, PARAMETER ::   jpr_fraqsr = 42            ! fraction of solar net radiation absorbed in the first ocean level
-   INTEGER, PARAMETER ::   jprcv      = 42            ! total number of fields received
+   INTEGER, PARAMETER ::   jpr_ts_ice = 43            ! skin temperature of sea-ice (used for melt-ponds)
+   INTEGER, PARAMETER ::   jprcv      = 43            ! total number of fields received
 
    INTEGER, PARAMETER ::   jps_fice   =  1            ! ice fraction sent to the atmosphere
@@ -135,5 +133,10 @@
    INTEGER, PARAMETER ::   jps_e3t1st = 27            ! first level depth (vvl)
    INTEGER, PARAMETER ::   jps_fraqsr = 28            ! fraction of solar net radiation absorbed in the first ocean level
-   INTEGER, PARAMETER ::   jpsnd      = 28            ! total number of fields sended
+   INTEGER, PARAMETER ::   jps_a_p    = 29            ! meltpond fraction  
+   INTEGER, PARAMETER ::   jps_ht_p   = 30            ! meltpond depth (m) 
+   INTEGER, PARAMETER ::   jps_kice   = 31            ! ice surface layer thermal conductivity
+   INTEGER, PARAMETER ::   jps_sstfrz = 32            ! sea-surface freezing temperature
+   INTEGER, PARAMETER ::   jps_fice1  = 33            ! first-order ice concentration (for time-travelling ice coupling)
+   INTEGER, PARAMETER ::   jpsnd      = 33            ! total number of fields sended
 
    !                                                         !!** namelist namsbc_cpl **
@@ -146,8 +149,9 @@
    END TYPE FLD_C
    ! Send to the atmosphere                           !
-   TYPE(FLD_C) ::   sn_snd_temp, sn_snd_alb, sn_snd_thick, sn_snd_crt, sn_snd_co2                        
+   TYPE(FLD_C) ::   sn_snd_temp, sn_snd_alb, sn_snd_thick, sn_snd_crt, sn_snd_co2, sn_snd_cond, sn_snd_mpnd, sn_snd_sstfrz, sn_snd_thick1
+
    ! Received from the atmosphere                     !
    TYPE(FLD_C) ::   sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau, sn_rcv_dqnsdt, sn_rcv_qsr, sn_rcv_qns, sn_rcv_emp, sn_rcv_rnf
-   TYPE(FLD_C) ::   sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2                        
+   TYPE(FLD_C) ::   sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2, sn_rcv_ts_ice 
    ! Other namelist parameters                        !
    INTEGER     ::   nn_cplmodel            ! Maximum number of models to/from which NEMO is potentialy sending/receiving data
@@ -217,7 +221,8 @@
       !!
       NAMELIST/namsbc_cpl/  sn_snd_temp, sn_snd_alb   , sn_snd_thick, sn_snd_crt   , sn_snd_co2,      &
+         &                  sn_snd_cond, sn_snd_mpnd, sn_snd_sstfrz,                                  &
          &                  sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau  , sn_rcv_dqnsdt, sn_rcv_qsr,      &
          &                  sn_rcv_qns , sn_rcv_emp   , sn_rcv_rnf  , sn_rcv_cal   , sn_rcv_iceflx,   &
-         &                  sn_rcv_co2 , nn_cplmodel  , ln_usecplmask
+         &                  sn_rcv_co2 , sn_rcv_ts_ice, nn_cplmodel , ln_usecplmask
       !!---------------------------------------------------------------------
       !
@@ -269,4 +274,8 @@
          WRITE(numout,*)'                      - mesh          = ', sn_snd_crt%clvgrd
          WRITE(numout,*)'      oce co2 flux                    = ', TRIM(sn_snd_co2%cldes   ), ' (', TRIM(sn_snd_co2%clcat   ), ')'
+         WRITE(numout,*)'      ice effective conductivity      = ', TRIM(sn_snd_cond%cldes   ), ' (', TRIM(sn_snd_cond%clcat   ), ')'
+         WRITE(numout,*)'      meltponds fraction & depth      = ', TRIM(sn_snd_mpnd%cldes  ), ' (', TRIM(sn_snd_mpnd%clcat   ), ')'
+         WRITE(numout,*)'      sea surface freezing temp       = ', TRIM(sn_snd_sstfrz%cldes   ), ' (', TRIM(sn_snd_sstfrz%clcat   ), ')'
+
          WRITE(numout,*)'  nn_cplmodel                         = ', nn_cplmodel
          WRITE(numout,*)'  ln_usecplmask                       = ', ln_usecplmask
@@ -383,5 +392,5 @@
       srcv(jpr_snow)%clname = 'OTotSnow'      ! Snow = solid precipitation
       srcv(jpr_tevp)%clname = 'OTotEvap'      ! total evaporation (over oce + ice sublimation)
-      srcv(jpr_ievp)%clname = 'OIceEvap'      ! evaporation over ice = sublimation
+      srcv(jpr_ievp)%clname = 'OIceEvp'      ! evaporation over ice = sublimation
       srcv(jpr_sbpr)%clname = 'OSubMPre'      ! sublimation - liquid precipitation - solid precipitation 
       srcv(jpr_semp)%clname = 'OISubMSn'      ! ice solid water budget = sublimation - solid precipitation
@@ -396,5 +405,7 @@
       CASE default              ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_rcv_emp%cldes' )
       END SELECT
-
+      !Set the number of categories for coupling of sublimation
+      IF ( TRIM( sn_rcv_emp%clcat ) == 'yes' ) srcv(jpr_ievp)%nct = jpl
+      !
       !                                                      ! ------------------------- !
       !                                                      !     Runoffs & Calving     !   
@@ -483,4 +494,17 @@
          srcv(jpr_topm:jpr_botm)%laction = .TRUE.
       ENDIF
+      
+#if defined key_cice && ! defined key_cice4
+      !                                                      ! ----------------------------- !
+      !                                                      !  sea-ice skin temperature     !   
+      !                                                      !  used in meltpond scheme      !
+      !                                                      !  May be calculated in Atm     !
+      !                                                      ! ----------------------------- !
+      srcv(jpr_ts_ice)%clname = 'OTsfIce'
+      IF ( TRIM( sn_rcv_ts_ice%cldes ) == 'ice' ) srcv(jpr_ts_ice)%laction = .TRUE.
+      IF ( TRIM( sn_rcv_ts_ice%clcat ) == 'yes' ) srcv(jpr_ts_ice)%nct = jpl
+      !TODO: Should there be a consistency check here?
+#endif
+
       !                                                      ! ------------------------------- !
       !                                                      !   OPA-SAS coupling - rcv by opa !   
@@ -600,10 +624,10 @@
       !                                                      ! ------------------------- !
       ssnd(jps_toce)%clname = 'O_SSTSST'
-      ssnd(jps_tice)%clname = 'O_TepIce'
+      ssnd(jps_tice)%clname = 'OTepIce'
       ssnd(jps_tmix)%clname = 'O_TepMix'
       SELECT CASE( TRIM( sn_snd_temp%cldes ) )
       CASE( 'none'                                 )       ! nothing to do
       CASE( 'oce only'                             )   ;   ssnd( jps_toce )%laction = .TRUE.
-      CASE( 'oce and ice' , 'weighted oce and ice' )
+      CASE( 'oce and ice' , 'weighted oce and ice' , 'oce and weighted ice')
          ssnd( (/jps_toce, jps_tice/) )%laction = .TRUE.
          IF ( TRIM( sn_snd_temp%clcat ) == 'yes' )  ssnd(jps_tice)%nct = jpl
@@ -634,13 +658,19 @@
 
       !                                                      ! ------------------------- !
-      !                                                      !  Ice fraction & Thickness ! 
+      !                                                      !  Ice fraction & Thickness 
       !                                                      ! ------------------------- !
       ssnd(jps_fice)%clname = 'OIceFrc'
       ssnd(jps_hice)%clname = 'OIceTck'
       ssnd(jps_hsnw)%clname = 'OSnwTck'
+      ssnd(jps_a_p)%clname  = 'OPndFrc'
+      ssnd(jps_ht_p)%clname = 'OPndTck'
+      ssnd(jps_fice1)%clname = 'OIceFrd'
       IF( k_ice /= 0 ) THEN
          ssnd(jps_fice)%laction = .TRUE.                  ! if ice treated in the ocean (even in climato case)
+         ssnd(jps_fice1)%laction = .TRUE.                 ! First-order regridded ice concentration, to be used
+	                                                  ! in producing atmos-to-ice fluxes
 ! Currently no namelist entry to determine sending of multi-category ice fraction so use the thickness entry for now
          IF ( TRIM( sn_snd_thick%clcat ) == 'yes' ) ssnd(jps_fice)%nct = jpl
+         IF ( TRIM( sn_snd_thick1%clcat ) == 'yes' ) ssnd(jps_fice1)%nct = jpl
       ENDIF
       
@@ -657,4 +687,41 @@
       CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_thick%cldes' )
       END SELECT
+
+      !                                                      ! ------------------------- !
+      !                                                      ! Ice Meltponds             !
+      !                                                      ! ------------------------- !
+#if defined key_cice && ! defined key_cice4
+      ! Meltponds only CICE5 
+      ssnd(jps_a_p)%clname = 'OPndFrc'   
+      ssnd(jps_ht_p)%clname = 'OPndTck'   
+      SELECT CASE ( TRIM( sn_snd_mpnd%cldes ) )
+      CASE ( 'none' )
+         ssnd(jps_a_p)%laction = .FALSE.
+         ssnd(jps_ht_p)%laction = .FALSE.
+      CASE ( 'ice only' ) 
+         ssnd(jps_a_p)%laction = .TRUE.
+         ssnd(jps_ht_p)%laction = .TRUE.
+         IF ( TRIM( sn_snd_mpnd%clcat ) == 'yes' ) THEN
+            ssnd(jps_a_p)%nct = jpl
+            ssnd(jps_ht_p)%nct = jpl
+         ELSE
+            IF ( jpl > 1 ) THEN
+               CALL ctl_stop( 'sbc_cpl_init: use weighted ice option for sn_snd_mpnd%cldes if not exchanging category fields' )
+            ENDIF
+         ENDIF
+      CASE ( 'weighted ice' ) 
+         ssnd(jps_a_p)%laction = .TRUE.
+         ssnd(jps_ht_p)%laction = .TRUE.
+         IF ( TRIM( sn_snd_mpnd%clcat ) == 'yes' ) THEN
+            ssnd(jps_a_p)%nct = jpl 
+            ssnd(jps_ht_p)%nct = jpl 
+         ENDIF
+      CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_mpnd%cldes' )
+      END SELECT
+#else
+      IF( TRIM( sn_snd_mpnd%cldes /= 'none' ) THEN
+         CALL ctl_stop('Meltponds can only be used with CICEv5')
+      ENDIF
+#endif
 
       !                                                      ! ------------------------- !
@@ -689,4 +756,36 @@
       !                                                      ! ------------------------- !
       ssnd(jps_co2)%clname = 'O_CO2FLX' ;  IF( TRIM(sn_snd_co2%cldes) == 'coupled' )    ssnd(jps_co2 )%laction = .TRUE.
+      !
+      
+      !                                                      ! ------------------------- !
+      !                                                      ! Sea surface freezing temp !
+      !                                                      ! ------------------------- !
+      ssnd(jps_sstfrz)%clname = 'O_SSTFrz' ; IF( TRIM(sn_snd_sstfrz%cldes) == 'coupled' )  ssnd(jps_sstfrz)%laction = .TRUE.
+      !
+      !                                                      ! ------------------------- !
+      !                                                      !    Ice conductivity       !
+      !                                                      ! ------------------------- !
+      ! Note that ultimately we will move to passing an ocean effective conductivity as well so there
+      ! will be some changes to the parts of the code which currently relate only to ice conductivity
+      ssnd(jps_kice )%clname = 'OIceKn'
+      SELECT CASE ( TRIM( sn_snd_cond%cldes ) )
+      CASE ( 'none' )
+         ssnd(jps_kice)%laction = .FALSE.
+      CASE ( 'ice only' )
+         ssnd(jps_kice)%laction = .TRUE.
+         IF ( TRIM( sn_snd_cond%clcat ) == 'yes' ) THEN
+            ssnd(jps_kice)%nct = jpl
+         ELSE
+            IF ( jpl > 1 ) THEN
+               CALL ctl_stop( 'sbc_cpl_init: use weighted ice option for sn_snd_cond%cldes if not exchanging category fields' )
+            ENDIF
+         ENDIF
+      CASE ( 'weighted ice' )
+         ssnd(jps_kice)%laction = .TRUE.
+         IF ( TRIM( sn_snd_cond%clcat ) == 'yes' ) ssnd(jps_kice)%nct = jpl
+      CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_cond%cldes' )
+      END SELECT
+      !
+      
 
       !                                                      ! ------------------------------- !
@@ -843,5 +942,5 @@
       !!
       LOGICAL  ::   llnewtx, llnewtau      ! update wind stress components and module??
-      INTEGER  ::   ji, jj, jn             ! dummy loop indices
+      INTEGER  ::   ji, jj, jl, jn         ! dummy loop indices
       INTEGER  ::   isec                   ! number of seconds since nit000 (assuming rdttra did not change since nit000)
       REAL(wp) ::   zcumulneg, zcumulpos   ! temporary scalars     
@@ -995,4 +1094,23 @@
       !                                                      ! ================== !
       IF( srcv(jpr_co2)%laction )   atm_co2(:,:) = frcv(jpr_co2)%z3(:,:,1)
+#endif
+
+#if defined key_cice && ! defined key_cice4
+      !  ! Sea ice surface skin temp:
+      IF( srcv(jpr_ts_ice)%laction ) THEN
+        DO jl = 1, jpl
+          DO jj = 1, jpj
+            DO ji = 1, jpi
+              IF (frcv(jpr_ts_ice)%z3(ji,jj,jl) > 0.0) THEN
+                tsfc_ice(ji,jj,jl) = 0.0
+              ELSE IF (frcv(jpr_ts_ice)%z3(ji,jj,jl) < -60.0) THEN
+                tsfc_ice(ji,jj,jl) = -60.0
+              ELSE
+                tsfc_ice(ji,jj,jl) = frcv(jpr_ts_ice)%z3(ji,jj,jl)
+              ENDIF
+            END DO
+          END DO
+        END DO
+      ENDIF
 #endif
 
@@ -1403,6 +1521,33 @@
          zsprecip(:,:) = frcv(jpr_snow)%z3(:,:,1)                  ! May need to ensure positive here
          ztprecip(:,:) = frcv(jpr_rain)%z3(:,:,1) + zsprecip(:,:)  ! May need to ensure positive here
-         zemp_tot(:,:) = frcv(jpr_tevp)%z3(:,:,1) - ztprecip(:,:)
+         zemp_tot(:,:) = frcv(jpr_tevp)%z3(:,:,1) - ztprecip(:,:)         
+#if defined key_cice
+         IF ( TRIM(sn_rcv_emp%clcat) == 'yes' ) THEN
+            ! zemp_ice is the sum of frcv(jpr_ievp)%z3(:,:,1) over all layers - snow
+            zemp_ice(:,:) = - frcv(jpr_snow)%z3(:,:,1)
+            DO jl=1,jpl
+               zemp_ice(:,:   ) = zemp_ice(:,:) + frcv(jpr_ievp)%z3(:,:,jl)
+            ENDDO
+            ! latent heat coupled for each category in CICE
+            qla_ice(:,:,1:jpl) = - frcv(jpr_ievp)%z3(:,:,1:jpl) * lsub
+         ELSE
+            ! If CICE has multicategories it still expects coupling fields for
+            ! each even if we treat as a single field
+            ! The latent heat flux is split between the ice categories according
+            ! to the fraction of the ice in each category
+            zemp_ice(:,:) = frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_snow)%z3(:,:,1)
+            WHERE ( zicefr(:,:) /= 0._wp ) 
+               ztmp(:,:) = 1./zicefr(:,:)
+            ELSEWHERE 
+               ztmp(:,:) = 0.e0
+            END WHERE  
+            DO jl=1,jpl
+               qla_ice(:,:,jl) = - a_i(:,:,jl) * ztmp(:,:) * frcv(jpr_ievp)%z3(:,:,1) * lsub 
+            END DO
+            WHERE ( zicefr(:,:) == 0._wp )  qla_ice(:,:,1) = -frcv(jpr_ievp)%z3(:,:,1) * lsub 
+         ENDIF
+#else         
          zemp_ice(:,:) = frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_snow)%z3(:,:,1)
+#endif                  
             CALL iom_put( 'rain'         , frcv(jpr_rain)%z3(:,:,1)              )   ! liquid precipitation 
          IF( iom_use('hflx_rain_cea') )   &
@@ -1758,4 +1903,15 @@
                CASE default                  ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%clcat' )
                END SELECT
+            CASE( 'oce and weighted ice' )   ;   ztmp1(:,:) =   tsn(:,:,1,jp_tem) + rt0 
+               SELECT CASE( sn_snd_temp%clcat )
+               CASE( 'yes' )   
+        	  ztmp3(:,:,1:jpl) = tn_ice(:,:,1:jpl) * a_i(:,:,1:jpl)
+               CASE( 'no' )
+        	  ztmp3(:,:,:) = 0.0
+        	  DO jl=1,jpl
+                     ztmp3(:,:,1) = ztmp3(:,:,1) + tn_ice(:,:,jl) * a_i(:,:,jl)
+        	  ENDDO
+               CASE default                  ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%clcat' )
+               END SELECT
             CASE( 'mixed oce-ice'        )   
                ztmp1(:,:) = ( ztmp1(:,:) + rt0 ) * zfr_l(:,:) 
@@ -1799,4 +1955,14 @@
          END SELECT
          IF( ssnd(jps_fice)%laction )   CALL cpl_snd( jps_fice, isec, ztmp3, info )
+      ENDIF
+      
+      ! Send ice fraction field (first order interpolation), for weighting UM fluxes to be passed to NEMO
+      IF (ssnd(jps_fice1)%laction) THEN
+         SELECT CASE (sn_snd_thick1%clcat)
+         CASE( 'yes' )   ;   ztmp3(:,:,1:jpl) =  a_i(:,:,1:jpl)
+         CASE( 'no' )    ;   ztmp3(:,:,1) = fr_i(:,:)
+         CASE default    ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_thick1%clcat' )
+	 END SELECT
+         CALL cpl_snd (jps_fice1, isec, ztmp3, info)
       ENDIF
       
@@ -1844,4 +2010,48 @@
          IF( ssnd(jps_hsnw)%laction )   CALL cpl_snd( jps_hsnw, isec, ztmp4, info )
       ENDIF
+      !
+      ! Send meltpond fields 
+      IF( ssnd(jps_a_p)%laction .OR. ssnd(jps_ht_p)%laction ) THEN
+         SELECT CASE( sn_snd_mpnd%cldes) 
+         CASE( 'weighted ice' ) 
+            SELECT CASE( sn_snd_mpnd%clcat ) 
+            CASE( 'yes' ) 
+               ztmp3(:,:,1:jpl) =  a_p(:,:,1:jpl) * a_i(:,:,1:jpl) 
+               ztmp4(:,:,1:jpl) =  ht_p(:,:,1:jpl) * a_i(:,:,1:jpl) 
+            CASE( 'no' ) 
+               ztmp3(:,:,:) = 0.0 
+               ztmp4(:,:,:) = 0.0 
+               DO jl=1,jpl 
+                 ztmp3(:,:,1) = ztmp3(:,:,1) + a_p(:,:,jpl) * a_i(:,:,jpl) 
+                 ztmp4(:,:,1) = ztmp4(:,:,1) + ht_p(:,:,jpl) * a_i(:,:,jpl) 
+               ENDDO 
+            CASE default    ;   CALL ctl_stop( 'sbc_cpl_mpd: wrong definition of sn_snd_mpnd%clcat' ) 
+            END SELECT 
+         CASE( 'ice only' )    
+            ztmp3(:,:,1:jpl) = a_p(:,:,1:jpl) 
+            ztmp4(:,:,1:jpl) = ht_p(:,:,1:jpl) 
+         END SELECT 
+         IF( ssnd(jps_a_p)%laction )   CALL cpl_snd( jps_a_p, isec, ztmp3, info )    
+         IF( ssnd(jps_ht_p)%laction )   CALL cpl_snd( jps_ht_p, isec, ztmp4, info )    
+         !
+         ! Send ice effective conductivity
+         SELECT CASE( sn_snd_cond%cldes)
+         CASE( 'weighted ice' )   
+            SELECT CASE( sn_snd_cond%clcat )
+            CASE( 'yes' )   
+               ztmp3(:,:,1:jpl) =  kn_ice(:,:,1:jpl) * a_i(:,:,1:jpl)
+            CASE( 'no' )
+               ztmp3(:,:,:) = 0.0
+               DO jl=1,jpl
+                 ztmp3(:,:,1) = ztmp3(:,:,1) + kn_ice(:,:,jl) * a_i(:,:,jl)
+               ENDDO
+            CASE default                  ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_cond%clcat' )
+            END SELECT
+         CASE( 'ice only' )   
+           ztmp3(:,:,1:jpl) = kn_ice(:,:,1:jpl)
+         END SELECT
+         IF( ssnd(jps_kice)%laction )   CALL cpl_snd( jps_kice, isec, ztmp3, info )
+      ENDIF
+      !
       !
 #if defined key_cpl_carbon_cycle
@@ -2023,5 +2233,8 @@
       IF( ssnd(jps_rnf   )%laction )  CALL cpl_snd( jps_rnf   , isec, RESHAPE ( rnf , (/jpi,jpj,1/) ), info )
       IF( ssnd(jps_taum  )%laction )  CALL cpl_snd( jps_taum  , isec, RESHAPE ( taum, (/jpi,jpj,1/) ), info )
-
+      
+      ztmp1(:,:) = sstfrz(:,:) + rt0
+      IF( ssnd(jps_sstfrz)%laction )  CALL cpl_snd( jps_sstfrz, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )
+      !
       CALL wrk_dealloc( jpi,jpj, zfr_l, ztmp1, ztmp2, zotx1, zoty1, zotz1, zitx1, zity1, zitz1 )
       CALL wrk_dealloc( jpi,jpj,jpl, ztmp3, ztmp4 )
Index: branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90
===================================================================
--- branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90	(revision 5662)
+++ branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90	(revision 5663)
@@ -15,5 +15,6 @@
    USE dom_oce         ! ocean space and time domain
    USE domvvl
-   USE phycst, only : rcp, rau0, r1_rau0, rhosn, rhoic
+   USE eosbn2, only : eos_fzp ! Function to calculate freezing point of seawater
+   USE phycst, only : rcp, rau0, r1_rau0, rhosn, rhoic, rt0
    USE in_out_manager  ! I/O manager
    USE iom, ONLY : iom_put,iom_use              ! I/O manager library !!Joakim edit
@@ -37,6 +38,6 @@
    USE ice_gather_scatter
    USE ice_calendar, only: dt
+# if defined key_cice4
    USE ice_state, only: aice,aicen,uvel,vvel,vsno,vsnon,vice,vicen
-# if defined key_cice4
    USE ice_flux, only: strax,stray,strocnx,strocny,frain,fsnow,  &
                 strocnxT,strocnyT,                               & 
@@ -45,15 +46,24 @@
                 flatn_f,fsurfn_f,fcondtopn_f,                    &
                 uatm,vatm,wind,fsw,flw,Tair,potT,Qa,rhoa,zlvl,   &
-                swvdr,swvdf,swidr,swidf
+                swvdr,swvdf,swidr,swidf,Tf
    USE ice_therm_vertical, only: calc_Tsfc
 #else
+   USE ice_state, only: aice,aicen,uvel,nt_hpnd,trcrn,vvel,vsno,&
+                vsnon,vice,vicen,nt_Tsfc
    USE ice_flux, only: strax,stray,strocnx,strocny,frain,fsnow,  &
                 strocnxT,strocnyT,                               & 
-                sst,sss,uocn,vocn,ss_tltx,ss_tlty,fsalt_ai,     &
-                fresh_ai,fhocn_ai,fswthru_ai,frzmlt,          &
+                sst,sss,uocn,vocn,ss_tltx,ss_tlty,fsalt_ai,      &
+                fresh_ai,fhocn_ai,fswthru_ai,frzmlt,             &
                 flatn_f,fsurfn_f,fcondtopn_f,                    &
                 uatm,vatm,wind,fsw,flw,Tair,potT,Qa,rhoa,zlvl,   &
-                swvdr,swvdf,swidr,swidf
-   USE ice_therm_shared, only: calc_Tsfc
+                swvdr,swvdf,swidr,swidf,Tf,                      &
+		!! When using NEMO with CICE, this change requires use of 
+		!! one of the following two CICE branches:
+		!! - at CICE5.0,   hadax/r1015_GSI8_with_GSI7
+		!! - at CICE5.1.2, hadax/vn5.1.2_GSI8
+                keffn_top,Tn_top
+
+   USE ice_therm_shared, only: calc_Tsfc, heat_capacity
+   USE ice_shortwave, only: apeffn
 #endif
    USE ice_forcing, only: frcvdr,frcvdf,frcidr,frcidf
@@ -161,4 +171,5 @@
       INTEGER, INTENT( in  ) ::   ksbc                ! surface forcing type
       REAL(wp), DIMENSION(:,:), POINTER :: ztmp1, ztmp2
+      REAL(wp), DIMENSION(:,:,:), POINTER :: ztfrz3d
       REAL(wp) ::   zcoefu, zcoefv, zcoeff            ! local scalar
       INTEGER  ::   ji, jj, jl, jk                    ! dummy loop indices
@@ -173,4 +184,16 @@
       ji_off = INT ( (jpiglo - nx_global) / 2 )
       jj_off = INT ( (jpjglo - ny_global) / 2 )
+
+      ! Set freezing temperatures and ensure consistencey between NEMO and CICE
+      CALL wrk_alloc( jpi,jpj,jpk, ztfrz3d ) 
+      DO jk=1,jpk
+         ztfrz3d(:,:,jk) = eos_fzp( tsn(:,:,jk,jp_sal), fsdept_n(:,:,jk) )
+      ENDDO
+
+      !Ensure that no temperature points are below freezing if not a NEMO restart
+      IF( .NOT. ln_rstart ) THEN
+         tsn(:,:,:,jp_tem) = MAX (tsn(:,:,:,jp_tem),ztfrz3d)
+         tsb(:,:,:,jp_tem) = tsn(:,:,:,jp_tem)
+      ENDIF
 
 #if defined key_nemocice_decomp
@@ -202,9 +225,6 @@
       IF( sbc_ice_cice_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_ice_cice_alloc : unable to allocate cice arrays' )
 
-! Ensure ocean temperatures are nowhere below freezing if not a NEMO restart
-      IF( .NOT. ln_rstart ) THEN
-         tsn(:,:,:,jp_tem) = MAX (tsn(:,:,:,jp_tem),Tocnfrz)
-         tsb(:,:,:,jp_tem) = tsn(:,:,:,jp_tem)
-      ENDIF
+      ! Populate the surface freezing temperature array
+      sstfrz(:,:)=ztfrz3d(:,:,1)
 
       fr_iu(:,:)=0.0
@@ -283,4 +303,5 @@
  
       CALL wrk_dealloc( jpi,jpj, ztmp1, ztmp2 )
+      CALL wrk_dealloc( jpi,jpj,jpk, ztfrz3d ) 
       !
       IF( nn_timing == 1 )  CALL timing_stop('cice_sbc_init')
@@ -343,29 +364,26 @@
          CALL nemo2cice(ztmp,stray,'F', -1. )
 
+
+! Alex West: From configuration GSI8 onwards, when NEMO is used with CICE in
+! HadGEM3 the 'time-travelling ice' coupling approach is used, whereby 
+! atmosphere-ice fluxes are passed as pseudo-local values, formed by dividing
+! gridbox mean fluxes in the UM by future ice concentration obtained through 	
+! OASIS.  This allows for a much more realistic apportionment of energy through
+! the ice - and conserves energy.
+! Therefore the fluxes are now divided by ice concentration in the coupled
+! formulation (jp_purecpl) as well as for jp_flx.  This NEMO branch should only
+! be used at UM10.2 onwards (unless an explicit GSI8 UM branch is included), at
+! which point the GSI8 UM changes were committed.
+
 ! Surface downward latent heat flux (CI_5)
-         IF (ksbc == jp_flx) THEN
+         IF (ksbc == jp_flx .OR. ksbc == jp_purecpl) THEN
             DO jl=1,ncat
                ztmpn(:,:,jl)=qla_ice(:,:,1)*a_i(:,:,jl)
             ENDDO
          ELSE
-! emp_ice is set in sbc_cpl_ice_flx as sublimation-snow
-            qla_ice(:,:,1)= - ( emp_ice(:,:)+sprecip(:,:) ) * Lsub
-! End of temporary code
-            DO jj=1,jpj
-               DO ji=1,jpi
-                  IF (fr_i(ji,jj).eq.0.0) THEN
-                     DO jl=1,ncat
-                        ztmpn(ji,jj,jl)=0.0
-                     ENDDO
-                     ! This will then be conserved in CICE
-                     ztmpn(ji,jj,1)=qla_ice(ji,jj,1)
-                  ELSE
-                     DO jl=1,ncat
-                        ztmpn(ji,jj,jl)=qla_ice(ji,jj,1)*a_i(ji,jj,jl)/fr_i(ji,jj)
-                     ENDDO
-                  ENDIF
-               ENDDO
-            ENDDO
+           !In coupled mode - qla_ice calculated in sbc_cpl for each category
+           ztmpn(:,:,1:ncat)=qla_ice(:,:,1:ncat)
          ENDIF
+
          DO jl=1,ncat
             CALL nemo2cice(ztmpn(:,:,jl),flatn_f(:,:,jl,:),'T', 1. )
@@ -373,5 +391,5 @@
 ! GBM conductive flux through ice (CI_6)
 !  Convert to GBM
-            IF (ksbc == jp_flx) THEN
+            IF (ksbc == jp_flx .OR. ksbc == jp_purecpl) THEN
                ztmp(:,:) = botmelt(:,:,jl)*a_i(:,:,jl)
             ELSE
@@ -382,5 +400,5 @@
 ! GBM surface heat flux (CI_7)
 !  Convert to GBM
-            IF (ksbc == jp_flx) THEN
+            IF (ksbc == jp_flx .OR. ksbc == jp_purecpl) THEN
                ztmp(:,:) = (topmelt(:,:,jl)+botmelt(:,:,jl))*a_i(:,:,jl) 
             ELSE
@@ -442,10 +460,12 @@
       CALL nemo2cice(ztmp,frain,'T', 1. ) 
 
+! Recalculate freezing temperature and send to CICE 
+      sstfrz(:,:)=eos_fzp(sss_m(:,:), fsdept_n(:,:,1)) 
+      CALL nemo2cice(sstfrz,Tf,'T', 1. )
+
 ! Freezing/melting potential
 ! Calculated over NEMO leapfrog timestep (hence 2*dt)
-      nfrzmlt(:,:)=rau0*rcp*fse3t_m(:,:)*(Tocnfrz-sst_m(:,:))/(2.0*dt)
-
-      ztmp(:,:) = nfrzmlt(:,:)
-      CALL nemo2cice(ztmp,frzmlt,'T', 1. )
+      nfrzmlt(:,:)=rau0*rcp*fse3t_m(:,:)*(sstfrz(:,:)-sst_m(:,:))/(2.0*dt) 
+      CALL nemo2cice(nfrzmlt,frzmlt,'T', 1. )
 
 ! SST  and SSS
@@ -453,4 +473,9 @@
       CALL nemo2cice(sst_m,sst,'T', 1. )
       CALL nemo2cice(sss_m,sss,'T', 1. )
+
+! Sea ice surface skin temperature
+      DO jl=1,ncat
+        CALL nemo2cice(tsfc_ice(:,:,jl), trcrn(:,:,nt_tsfc,jl,:),'T',1.)
+      ENDDO 
 
 ! x comp and y comp of surface ocean current
@@ -730,4 +755,32 @@
          CALL cice2nemo(vicen(:,:,jl,:),ht_i(:,:,jl),'T', 1. )
       ENDDO
+
+#if ! defined key_cice4
+! Meltpond fraction and depth
+      DO jl = 1,ncat
+         CALL cice2nemo(apeffn(:,:,jl,:),a_p(:,:,jl),'T', 1. )
+         CALL cice2nemo(trcrn(:,:,nt_hpnd,jl,:),ht_p(:,:,jl),'T', 1. )
+      ENDDO
+#endif
+
+
+! If using multilayers thermodynamics in CICE then get top layer temperature
+! and effective conductivity       
+!! When using NEMO with CICE, this change requires use of 
+!! one of the following two CICE branches:
+!! - at CICE5.0,   hadax/r1015_GSI8_with_GSI7
+!! - at CICE5.1.2, hadax/vn5.1.2_GSI8
+      IF (heat_capacity) THEN
+         DO jl = 1,ncat
+            CALL cice2nemo(Tn_top(:,:,jl,:),tn_ice(:,:,jl),'T', 1. )
+            CALL cice2nemo(keffn_top(:,:,jl,:),kn_ice(:,:,jl),'T', 1. )
+         ENDDO
+! Convert surface temperature to Kelvin
+         tn_ice(:,:,:)=tn_ice(:,:,:)+rt0
+      ELSE
+         tn_ice(:,:,:) = 0.0
+         kn_ice(:,:,:) = 0.0
+      ENDIF       
+
       !
       IF( nn_timing == 1 )  CALL timing_stop('cice_sbc_hadgam')
