Index: NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icedyn_adv_pra.F90
===================================================================
--- NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icedyn_adv_pra.F90	(revision 13079)
+++ NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icedyn_adv_pra.F90	(revision 13080)
@@ -44,4 +44,5 @@
    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   sxap , syap , sxxap , syyap , sxyap    ! melt pond fraction
    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   sxvp , syvp , sxxvp , syyvp , sxyvp    ! melt pond volume
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   sxlhp, sylhp, sxxlhp, syylhp, sxylhp   ! melt pond lid thickness
 
    !! * Substitutions
@@ -55,5 +56,5 @@
 
    SUBROUTINE ice_dyn_adv_pra( kt, pu_ice, pv_ice,  &
-      &                        pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, pe_s, pe_i )
+      &                        pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, plh_ip, pe_s, pe_i )
       !!----------------------------------------------------------------------
       !!                **  routine ice_dyn_adv_pra  **
@@ -78,4 +79,5 @@
       REAL(wp), DIMENSION(:,:,:)  , INTENT(inout) ::   pa_ip      ! melt pond fraction
       REAL(wp), DIMENSION(:,:,:)  , INTENT(inout) ::   pv_ip      ! melt pond volume
+      REAL(wp), DIMENSION(:,:,:)  , INTENT(inout) ::   plh_ip     ! melt pond lid thickness
       REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) ::   pe_s       ! snw heat content
       REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) ::   pe_i       ! ice heat content
@@ -88,5 +90,5 @@
       REAL(wp), DIMENSION(jpi,jpj,1)          ::   z0opw
       REAL(wp), DIMENSION(jpi,jpj,jpl)        ::   z0ice, z0snw, z0ai, z0smi, z0oi
-      REAL(wp), DIMENSION(jpi,jpj,jpl)        ::   z0ap , z0vp
+      REAL(wp), DIMENSION(jpi,jpj,jpl)        ::   z0ap , z0vp, z0lhp
       REAL(wp), DIMENSION(jpi,jpj,nlay_s,jpl) ::   z0es
       REAL(wp), DIMENSION(jpi,jpj,nlay_i,jpl) ::   z0ei
@@ -129,4 +131,5 @@
             z0ap(:,:,jl)  = pa_ip(:,:,jl) * e1e2t(:,:)     ! Melt pond fraction
             z0vp(:,:,jl)  = pv_ip(:,:,jl) * e1e2t(:,:)     ! Melt pond volume
+            z0lhp(:,:,jl)  = plh_ip(:,:,jl) * e1e2t(:,:)   ! Melt pond lid thickness
          ENDIF
       END DO
@@ -167,4 +170,6 @@
                CALL adv_x( zdt , pu_ice , 1._wp , zarea , z0vp , sxvp , sxxvp , syvp , syyvp , sxyvp )    !--- melt pond volume
                CALL adv_y( zdt , pv_ice , 0._wp , zarea , z0vp , sxvp , sxxvp , syvp , syyvp , sxyvp ) 
+               CALL adv_x( zdt , pu_ice , 1._wp , zarea , z0lhp, sxlhp, sxxlhp, sylhp, syylhp, sxylhp)    !--- melt pond lid thickness
+               CALL adv_y( zdt , pv_ice , 0._wp , zarea , z0lhp, sxlhp, sxxlhp, sylhp, syylhp, sxylhp)
             ENDIF
          END DO
@@ -202,4 +207,6 @@
                CALL adv_y( zdt , pv_ice , 1._wp , zarea , z0vp , sxvp , sxxvp , syvp , syyvp , sxyvp )    !--- melt pond volume
                CALL adv_x( zdt , pu_ice , 0._wp , zarea , z0vp , sxvp , sxxvp , syvp , syyvp , sxyvp )
+               CALL adv_y( zdt , pv_ice , 1._wp , zarea , z0lhp, sxlhp, sxxlhp, sylhp, syylhp, sxylhp)    !--- melt pond lid thickness
+               CALL adv_x( zdt , pu_ice , 0._wp , zarea , z0lhp, sxlhp, sxxlhp, sylhp, syylhp, sxylhp) 
             ENDIF
          END DO
@@ -225,4 +232,5 @@
             pa_ip(:,:,jl) = z0ap(:,:,jl) * r1_e1e2t(:,:) * tmask(:,:,1)
             pv_ip(:,:,jl) = z0vp(:,:,jl) * r1_e1e2t(:,:) * tmask(:,:,1)
+            plh_ip(:,:,jl) = z0lhp (:,:,jl) * r1_e1e2t(:,:) * tmask(:,:,1)
          ENDIF
       END DO
@@ -231,5 +239,5 @@
       !     Remove negative values (conservation is ensured)
       !     (because advected fields are not perfectly bounded and tiny negative values can occur, e.g. -1.e-20)
-      CALL ice_var_zapneg( zdt, pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, pe_s, pe_i )
+      CALL ice_var_zapneg( zdt, pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, plh_ip, pe_s, pe_i )
       !
       ! --- Ensure snow load is not too big --- !
@@ -653,4 +661,5 @@
          &      sxap(jpi,jpj,jpl)  , syap (jpi,jpj,jpl) , sxxap (jpi,jpj,jpl) , syyap (jpi,jpj,jpl) , sxyap (jpi,jpj,jpl) ,   &
          &      sxvp(jpi,jpj,jpl)  , syvp (jpi,jpj,jpl) , sxxvp (jpi,jpj,jpl) , syyvp (jpi,jpj,jpl) , sxyvp (jpi,jpj,jpl) ,   &
+         &      sxlhp(jpi,jpj,jpl) , sylhp (jpi,jpj,jpl), sxxlhp (jpi,jpj,jpl), syylhp (jpi,jpj,jpl), sxylhp (jpi,jpj,jpl),   &
          !
          &      sxc0 (jpi,jpj,nlay_s,jpl) , syc0 (jpi,jpj,nlay_s,jpl) , sxxc0(jpi,jpj,nlay_s,jpl) , &
@@ -765,4 +774,10 @@
                CALL iom_get( numrir, jpdom_autoglo, 'syyvp', syyvp )
                CALL iom_get( numrir, jpdom_autoglo, 'sxyvp', sxyvp )
+               !                                                     ! melt pond lid thickness
+               CALL iom_get( numrir, jpdom_autoglo, 'sxlhp' , sxlhp  )
+               CALL iom_get( numrir, jpdom_autoglo, 'sylhp' , sylhp  )
+               CALL iom_get( numrir, jpdom_autoglo, 'sxxlhp', sxxlhp )
+               CALL iom_get( numrir, jpdom_autoglo, 'syylhp', syylhp )
+               CALL iom_get( numrir, jpdom_autoglo, 'sxylhp', sxylhp )
             ENDIF
             !
@@ -782,4 +797,5 @@
                sxap  = 0._wp   ;   syap  = 0._wp   ;   sxxap  = 0._wp   ;   syyap  = 0._wp   ;   sxyap  = 0._wp   ! melt pond fraction
                sxvp  = 0._wp   ;   syvp  = 0._wp   ;   sxxvp  = 0._wp   ;   syyvp  = 0._wp   ;   sxyvp  = 0._wp   ! melt pond volume
+               sxlhp  = 0._wp  ;   sylhp  = 0._wp  ;   sxxlhp  = 0._wp  ;   syylhp  = 0._wp  ;   sxylhp  = 0._wp  ! melt pond lid thickness
             ENDIF
          ENDIF
@@ -862,4 +878,10 @@
             CALL iom_rstput( iter, nitrst, numriw, 'syyvp', syyvp )
             CALL iom_rstput( iter, nitrst, numriw, 'sxyvp', sxyvp )
+            !                                                        ! melt pond lid thickness
+            CALL iom_rstput( iter, nitrst, numriw, 'sxlhp' , sxlhp  )
+            CALL iom_rstput( iter, nitrst, numriw, 'sylhp' , sylhp  )
+            CALL iom_rstput( iter, nitrst, numriw, 'sxxlhp', sxxlhp )
+            CALL iom_rstput( iter, nitrst, numriw, 'syylhp', syylhp )
+            CALL iom_rstput( iter, nitrst, numriw, 'sxylhp', sxylhp )
          ENDIF
          !
