Index: /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/ice.F90
===================================================================
--- /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/ice.F90	(revision 8060)
+++ /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/ice.F90	(revision 8061)
@@ -281,6 +281,8 @@
    ! MV MP 2016
    !                                     !!** melt pond namelist (namicemp)
-   LOGICAL , PUBLIC ::   ln_limMP         !: activate melt ponds or not
-   INTEGER , PUBLIC ::   nn_limMP         !: type of melt pond implementation (1=full, 2=radiation only, 3=freshwater only)
+   LOGICAL , PUBLIC ::   ln_pnd           !: activate ponds or not
+   INTEGER , PUBLIC ::   nn_pnd_scheme    !: type of melt pond scheme:   =0 prescribed, =1 empirical, =2 topographic
+   INTEGER , PUBLIC ::   nn_pnd_cpl       !: type of melt pond coupling: =0 passive, =1 full, =2 radiation only, =3 freshwater only
+   REAL(wp), PUBLIC ::   rn_apnd          !: prescribed pond fraction (0<rn_apnd<1), only if nn_pnd_scheme = 0
    ! END MV MP 2016
 
Index: /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90
===================================================================
--- /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90	(revision 8060)
+++ /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90	(revision 8061)
@@ -417,13 +417,31 @@
       !--------------------------------------------------------------------
       !
-      IF ( ln_limMP ) THEN
-         a_ip(:,:,:) = 0.1_wp
-         v_ip(:,:,:) = 0.1_wp
-         h_ip(:,:,:) = 0._wp
-         a_ip_frac(:,:,:) = 0._wp
-
-         at_ip(:,:) = 0._wp
-         vt_ip(:,:) = 0._wp
-      ENDIF
+      !IF ( ln_pnd   ) THEN
+      SELECT CASE ( nn_pnd_scheme )
+
+         z1_jpl =  1 / REAL(jpl)
+
+         CASE ( 0 )           !--- Prescribed melt ponds
+
+            DO jl = 1, jpl
+
+               a_ip(:,:,jl) = rn_apnd * z1_jpl * zswitch(:,:)
+               h_ip(:,:,jl) = 0.1 * zswitch(:,:)
+
+            END DO
+
+         CASE ( 1, 2 )        !--- Prognostic melt ponds
+
+            DO jl = 1, jpl
+
+               a_ip(:,:,jl) = 0.1 * zswitch(:,:)
+               h_ip(:,:,jl) = 0.1 * zswitch(:,:)
+
+            END DO
+
+      END SELECT
+
+      v_ip(:,:,:)      = a_ip(:,:,:)  * h_i_p(:,:,:)
+      a_ip_frac(:,:,:) = a_ip(:,:,:) / a_i(:,:,:)
 
       ! END MV MP 2016
@@ -473,5 +491,5 @@
 
       ! MV MP 2016
-      IF ( ln_limMP ) THEN
+      IF ( nn_pnd_scheme > 1 ) THEN
          sxap  (:,:,:) = 0._wp    ; sxvp  (:,:,:) = 0._wp 
          syap  (:,:,:) = 0._wp    ; syvp  (:,:,:) = 0._wp 
Index: /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90
===================================================================
--- /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90	(revision 8060)
+++ /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90	(revision 8061)
@@ -608,5 +608,5 @@
             esrdg (ij) = e_s  (ji,jj,1,jl1) * afrac(ij)
             !MV MP 2016
-            IF ( ln_limMP ) THEN
+            IF ( nn_pnd_scheme > 1 ) THEN
                vprdg(ij) = v_ip(ji,jj, jl1) * afrac(ij)
             ENDIF
@@ -620,5 +620,5 @@
             vsrft (ij) = v_s  (ji,jj,  jl1) * afrft(ij)
             !MV MP 2016
-            IF ( ln_limMP ) THEN
+            IF ( nn_pnd_scheme > 1 ) THEN
                vprft(ij) = v_ip(ji,jj,jl1) * afrft(ij)
             ENDIF
@@ -665,5 +665,5 @@
             !------------------------------------------            
             !  Place part of the melt pond volume into the ocean. 
-            IF ( ln_limMP .AND. ( ( nn_limMP .EQ. 3 ) .OR. ( nn_limMP .EQ. 1 ) ) )  THEN
+            IF ( ( nn_pnd_scheme > 0 ) .AND. ( ( nn_pnd_cpl .EQ. 1 ) .OR. ( nn_pnd_cpl .EQ. 3 ) ) )
                wfx_pnd(ji,jj) = wfx_pnd(ji,jj) + ( rhofw * vprdg(ij) * ( 1._wp - rn_fpondrdg )   & 
                &                                 + rhofw * vprft(ij) * ( 1._wp - rn_fpondrft ) ) * r1_rdtice  ! fresh water source for ocean
@@ -688,5 +688,5 @@
 
             ! MV MP 2016
-            IF ( ln_limMP ) THEN
+            IF ( nn_pnd_scheme > 0 ) THEN
                v_ip (ji,jj,jl1) = v_ip (ji,jj,jl1) - vprdg (ij) - vprft (ij)
             ENDIF
@@ -757,5 +757,5 @@
                   &                                        esrft (ij) * rn_fsnowrft * zswitch(ij) )
                ! MV MP 2016
-               IF ( ln_limMP ) THEN
+               IF ( nn_pnd_scheme > 0 ) THEN
                   v_ip (ji,jj,jl2) = v_ip (ji,jj,jl2)  + ( vprdg (ij) * rn_fpondrdg * fvol(ij)  +  &
                   &                                        vprft (ij) * rn_fpondrft * zswitch(ij) )
Index: /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limitd_th.F90
===================================================================
--- /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limitd_th.F90	(revision 8060)
+++ /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limitd_th.F90	(revision 8061)
@@ -358,5 +358,5 @@
             a_i (ii,ij,1) = a_i(ii,ij,1) * ht_i(ii,ij,1) / rn_himin 
             ! MV MP 2016
-            IF ( ln_limMP ) THEN
+            IF ( nn_pnd_scheme > 0 ) THEN
                a_ip(ii,ij,1) = a_ip(ii,ij,1) * ht_i(ii,ij,1) / rn_himin
             ENDIF
@@ -587,5 +587,5 @@
 
             ! MV MP 2016 
-            IF ( ln_limMP ) THEN
+            IF ( nn_pnd_scheme > 0 ) THEN
             !---------------------
             ! Pond fraction
Index: /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limmp.F90
===================================================================
--- /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limmp.F90	(revision 8060)
+++ /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limmp.F90	(revision 8061)
@@ -83,5 +83,5 @@
       !!-------------------------------------------------------------------
       INTEGER  ::   ios                 ! Local integer output status for namelist read
-      NAMELIST/namicemp/  ln_limMP, nn_limMP
+      NAMELIST/namicemp/  ln_pnd, nn_pnd_scheme, nn_pnd_cpl, rn_apnd
       !!-------------------------------------------------------------------
 
@@ -99,9 +99,13 @@
          WRITE(numout,*) 'lim_mp_init : ice parameters for melt ponds'
          WRITE(numout,*) '~~~~~~~~~~~~'
-         WRITE(numout,*)'    Activate melt ponds                                         ln_limMP      = ', ln_limMP  
-         WRITE(numout,*)'    Type of melt pond implementation (1=all,2=rad, 3=fw, 4=no)  nn_limMP      = ', nn_limMP  
+         WRITE(numout,*)'    Activate melt ponds                                         ln_pnd        = ', ln_pnd
+         WRITE(numout,*)'    Type of melt pond scheme =0 presc, =1 empirical = 2 topo    nn_pnd_scheme = ', nn_pnd_scheme
+         WRITE(numout,*)'    Type of melt pond coupling =0 pass., =1 full, =2 rad, 3=fw  nn_pnd_cpl    = ', nn_pnd_cp
+         WRITE(numout,*)'    Prescribed pond fraction                                    rn_apnd       = ', rn_apnd
       ENDIF
       !
    END SUBROUTINE lim_mp_init
+
+
 
    SUBROUTINE lim_mp( kt )
@@ -128,10 +132,17 @@
       IF( nn_timing == 1 )  CALL timing_start('lim_mp')
 
-      CALL lim_mp_cesm ! test routine with hyper simple melt ponds
-
-      !CALL lim_mp_topo    (at_i, a_i,                                       &
-      !          &          vt_i, v_i, v_s,            t_i, s_i, a_ip_frac,  &
-      !          &          h_ip,     t_su)
-
+      SELECT CASE ( nn_pnd_scheme )
+
+          CASE (1)
+
+             CALL lim_mp_cesm ! empirical melt ponds
+
+          CASE (2)
+
+             CALL lim_mp_topo    (at_i, a_i,                                       &
+                       &          vt_i, v_i, v_s,            t_i, s_i, a_ip_frac,  &
+                       &          h_ip,     t_su)
+
+      END SELECT
 
       ! we should probably not aggregate here since we do it in lim_var_agg
@@ -142,4 +153,6 @@
 
    END SUBROUTINE lim_mp 
+
+
 
    SUBROUTINE lim_mp_cesm
@@ -149,5 +162,7 @@
        !! ** Purpose    : Compute melt pond evolution
        !!
-       !! ** Method     : Accumulation of meltwater and exponential release
+       !! ** Method     : Empirical method. A fraction of meltwater is accumulated 
+       !!                 in pond volume. It is then released exponentially when
+       !!                 surface is freezingAccumulation of meltwater and exponential release
        !!
        !! ** Tunable parameters :
@@ -155,5 +170,5 @@
        !! 
        !! ** Note       : Stolen from CICE for quick test of the melt pond
-       !!                 interface
+       !!                 radiation and freshwater interfaces
        !!
        !! ** References : Holland, M. M. et al (J Clim 2012)
@@ -161,6 +176,6 @@
        !!-------------------------------------------------------------------
 
-       INTEGER, POINTER, DIMENSION(:,:)    :: indxi             ! compressed indices for cells with ice melting
-       INTEGER, POINTER, DIMENSION(:,:)    :: indxj             !
+       INTEGER, POINTER, DIMENSION(:)      :: indxi             ! compressed indices for cells with ice melting
+       INTEGER, POINTER, DIMENSION(:)      :: indxj             !
 
        REAL(wp), POINTER, DIMENSION(:,:)   :: zwfx_mlw          ! available meltwater for melt ponding
@@ -175,17 +190,20 @@
        REAL(wp)                            :: zhs               ! dummy snow depth
        REAL(wp)                            :: zTp               ! reference temperature
-       REAL(wp)                            :: zrexp             ! rate constant to refreeze melt ponds
        REAL(wp)                            :: zdTs              ! dummy temperature difference
        REAL(wp)                            :: z1_rhofw          ! inverse freshwater density
        REAL(wp)                            :: z1_zpnd_aspect    ! inverse pond aspect ratio
 
+       INTEGER                             :: ji, jj, jl, ij    ! loop indices
+       INTEGER                             :: icells            ! size of dummy array
+
        !!-------------------------------------------------------------------
 
-       CALL wrk_alloc( jpi,jpj, zwfx_mlw )
        CALL wrk_alloc( jpi*jpj, indxi, indxj)
+       CALL wrk_alloc( jpi,jpj,     zwfx_mlw )
+       CALL wrk_alloc( jpi,jpj,jpl, zrfrac   )
 
        z1_rhofw       = 1. / rhofw 
        z1_zpnd_aspect = 1. / zpnd_aspect
-       zTp            = -2 + rt0.
+       zTp            = -2. 
 
        !------------------------------------------------------------------
@@ -215,6 +233,6 @@
                IF ( a_i(ji,jj,jl) > epsi10 ) THEN
                   icells = icells + 1
-                  indxi(icells) = i
-                  indxj(icells) = j
+                  indxi(icells) = ji
+                  indxj(icells) = jj
                ENDIF
             END DO                 ! ji
@@ -231,21 +249,33 @@
             IF ( zhi < rn_himin) THEN   !--- Remove ponds on thin ice if ice is too thin
 
-               wfx_pnd(ji,jj)   = wfx_pnd(ji,jj) + v_ip(ji,jj,jl) !--- Give freshwater to the ocean
-
-               a_ip(ji,jj,jl)      = 0._wp                        !--- Dump ponds
+               a_ip(ji,jj,jl)      = 0._wp                               !--- Dump ponds
                v_ip(ji,jj,jl)      = 0._wp
                a_ip_frac(ji,jj,jl) = 0._wp
                h_ip(ji,jj,jl)      = 0._wp
 
+               IF ( ( nn_pnd_cpl .EQ. 1 ) .OR. ( nn_pnd_cpl .EQ. 3 ) ) & !--- Give freshwater to the ocean
+                  wfx_pnd(ji,jj)   = wfx_pnd(ji,jj) + v_ip(ji,jj,jl) 
+
+
             ELSE                        !--- Update pond characteristics
 
                !--- Add retained melt water
-               v_ip(ji,jj,jl) = v_ip(ji,jj,jl) + rfrac(ji,jj,jl) * z1_rhofw * zwfx_mlw(ji,jj) * a_i(ji,jj,jl)
+               v_ip(ji,jj,jl)      = v_ip(ji,jj,jl) + zrfrac(ji,jj,jl) * z1_rhofw * zwfx_mlw(ji,jj) * a_i(ji,jj,jl) * rdt_ice
 
                !--- Shrink pond due to refreezing
-               zdTs = MAX ( zTp - t_su(ji,jj,jl) , 0. )
+               zdTs                = MAX ( zTp - t_su(ji,jj,jl) + rt0 , 0. )
                
-               v_ip(ji,jj,jl) = v_ip(ji,jj,jl) * EXP( rexp * zdTs / zTp )
-           
+               zvpold              = v_ip(ji,jj,jl)
+
+               v_ip(ji,jj,jl)      = v_ip(ji,jj,jl) * EXP( zrexp * zdTs / zTp )
+
+               !--- Dump meltwater due to refreezing ( of course this is wrong
+               !--- but this parameterization is too simple )
+               IF ( ( nn_pnd_cpl .EQ. 1 ) .OR. ( nn_pnd_cpl .EQ. 3 ) ) THEN
+
+                  wfx_pnd(ji,jj)   = wfx_pnd(ji,jj) + rhofw * ( v_ip(ji,jj,jl) - zvpold ) * r1_rdtice
+
+               ENDIF
+
                a_ip_frac(ji,jj,jl) = MIN( 1._wp , SQRT( v_ip(ji,jj,jl) * z1_zpnd_aspect / a_i(ji,jj,jl) ) )
 
@@ -263,7 +293,17 @@
             ENDIF
 
-         END DO
-
-      END DO ! jpl
+          END DO
+
+       END DO ! jpl
+
+       !--- Remove retained meltwater from surface fluxes 
+
+       IF ( ( nn_pnd_cpl .EQ. 1 ) .OR. ( nn_pnd_cpl .EQ. 3 ) ) THEN
+
+           wfx_snw(:,:) = wfx_snw(:,:) *  ( 1. - zrmin - ( zrmax - zrmin ) * at_i(:,:) )
+
+           wfx_sum(:,:) = wfx_sum(:,:) *  ( 1. - zrmin - ( zrmax - zrmin ) * at_i(:,:) )
+
+       ENDIF
 
    END SUBROUTINE lim_mp_cesm
Index: /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limrst.F90
===================================================================
--- /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limrst.F90	(revision 8060)
+++ /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limrst.F90	(revision 8061)
@@ -152,5 +152,5 @@
 
       ! MV MP 2016
-      IF ( ln_limMP ) THEN
+      IF ( nn_pnd_scheme > 0 ) THEN
          DO jl = 1, jpl 
             WRITE(zchar,'(I2.2)') jl
@@ -316,5 +316,5 @@
          END DO
          ! MV MP 2016
-         IF ( ln_limMP ) THEN
+         IF ( nn_pnd_scheme > 0 ) THEN
             DO jl = 1, jpl 
                WRITE(zchar,'(I2.2)') jl
@@ -439,5 +439,5 @@
 
       ! MV MP 2016
-      IF ( ln_limMP ) THEN
+      IF ( nn_pnd_scheme > 0 ) THEN
          DO jl = 1, jpl 
             WRITE(zchar,'(I2.2)') jl
@@ -575,5 +575,5 @@
          END DO
          ! MV MP 2016
-         IF ( ln_limMP ) THEN
+         IF ( nn_pnd_scheme > 0 ) THEN
             DO jl = 1, jpl 
                WRITE(zchar,'(I2.2)') jl
Index: /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limtrp.F90
===================================================================
--- /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limtrp.F90	(revision 8060)
+++ /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limtrp.F90	(revision 8061)
@@ -79,7 +79,7 @@
       ! With melt ponds, we have to diffuse them
       ! We hard code the number of variables to diffuse
-      ! since we can't put an IF ( ln_limMP ) for a declaration
+      ! since we can't put an IF ( nn_pnd_scheme ) for a declaration
       ! ideally, the ihdf_vars should probably be passed as an argument and
-      ! defined somewhere depending on ln_limMP
+      ! defined somewhere depending on nn_pnd_scheme
       ! END MV MP 2016
       INTEGER , PARAMETER                    ::   ihdf_vars  = 8 ! Number of variables in which we apply horizontal diffusion
@@ -217,5 +217,5 @@
                CALL lim_adv_umx( kt, zdt, zudy, zvdx, zcu_box, zcv_box, e_s(:,:,1,jl) )    ! Snow heat content
                ! MV MP 2016
-               IF ( ln_limMP ) THEN
+               IF ( nn_pnd_scheme > 0 ) THEN
                   CALL lim_adv_umx( kt, zdt, zudy, zvdx, zcu_box, zcv_box, a_ip(:,:,jl) )  ! Melt pond fraction
                   CALL lim_adv_umx( kt, zdt, zudy, zvdx, zcu_box, zcv_box, v_ip(:,:,jl) )  ! Melt pond volume
@@ -265,5 +265,5 @@
             END DO
             ! MV MP 2016
-            IF ( ln_limMP ) THEN
+            IF ( nn_pnd_scheme > 0 ) THEN
                z0ap  (:,:,jl)  = a_ip (:,:,jl) * e12t(:,:) ! Melt pond fraction
                z0vp  (:,:,jl)  = v_ip (:,:,jl) * e12t(:,:) ! Melt pond volume
@@ -365,5 +365,5 @@
                   END DO
                   ! MV MP 2016
-                  IF ( ln_limMP ) THEN
+                  IF ( nn_pnd_scheme > 0 ) THEN
                      CALL lim_adv_y( zusnit, v_ice, 1._wp, zarea, z0ap  (:,:,jl), sxap (:,:,jl),   &   !--- melt pond fraction ---
                      &                                         sxxap (:,:,jl), syap (:,:,jl), syyap (:,:,jl), sxyap (:,:,jl)  )
@@ -395,5 +395,5 @@
             END DO
             ! MV MP 2016
-            IF ( ln_limMP ) THEN
+            IF ( nn_pnd_scheme > 0 ) THEN
                a_ip  (:,:,jl)   = z0ap (:,:,jl) * r1_e12t(:,:)
                v_ip  (:,:,jl)   = z0vp (:,:,jl) * r1_e12t(:,:)
@@ -442,5 +442,5 @@
             zhdfptab(:,:,jm)= e_s  (:,:,1,jl); jm = jm + 1
             ! MV MP 2016
-            IF ( ln_limMP ) THEN
+            IF ( nn_pnd_scheme > 0 ) THEN
                zhdfptab(:,:,jm)= a_ip  (:,:,  jl); jm = jm + 1
                zhdfptab(:,:,jm)= v_ip  (:,:,  jl); jm = jm + 1
@@ -481,5 +481,5 @@
             e_s  (:,:,1,jl) = zhdfptab(:,:,jm); jm = jm + 1
             ! MV MP 2016
-            IF ( ln_limMP ) THEN
+            IF ( nn_pnd_scheme > 0 ) THEN
                a_ip (:,:,  jl) = zhdfptab(:,:,jm); jm = jm + 1
                v_ip (:,:,  jl) = zhdfptab(:,:,jm); jm = jm + 1
@@ -539,5 +539,5 @@
 
                         ! MV MP 2016
-                        IF ( ln_limMP ) THEN
+                        IF ( nn_pnd_scheme > 0 ) THEN
                            a_ip (ji,jj,jl)        = rswitch * a_ip (ji,jj,jl)
                            v_ip (ji,jj,jl)        = rswitch * v_ip (ji,jj,jl)
@@ -563,5 +563,5 @@
                a_i (ji,jj,jpl) = v_i(ji,jj,jpl) / MAX( ht_i(ji,jj,jpl) , epsi20 ) * rswitch
                ! MV MP 2016
-               IF ( ln_limMP ) THEN
+               IF ( nn_pnd_scheme > 0 ) THEN
                   ! correct pond fraction to avoid a_ip > a_i
                   a_ip(ji,jj,jpl) = a_ip(ji,jj,jpl) * a_i(ji,jj,jpl) / MAX( za_old , epsi20 ) * rswitch
@@ -598,5 +598,5 @@
 
       ! MV MP 2016 (remove once we get rid of a_i_frac and ht_i)
-      IF ( ln_limMP ) THEN
+      IF ( nn_pnd_scheme > 0 ) THEN
           at_ip(:,:) = SUM( a_ip(:,:,:), dim = 3 )
           vt_ip(:,:) = SUM( v_ip(:,:,:), dim = 3 )
Index: /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limvar.F90
===================================================================
--- /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limvar.F90	(revision 8060)
+++ /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limvar.F90	(revision 8061)
@@ -92,5 +92,5 @@
 
       ! MV MP 2016
-      IF ( ln_limMP ) THEN
+      IF ( nn_pnd_scheme > 0 )
          at_ip(:,:) = SUM( a_ip, dim=3 )
          vt_ip(:,:) = SUM( v_ip, dim=3 )
Index: /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limwri.F90
===================================================================
--- /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limwri.F90	(revision 8060)
+++ /branches/2016/dev_r6859_LIM3_meltponds/NEMOGCM/NEMO/LIM_SRC_3/limwri.F90	(revision 8061)
@@ -208,5 +208,5 @@
 
       ! MV MP 2016
-      IF ( ln_limMP ) THEN
+      IF ( nn_pnd_scheme >= 0 ) THEN ! change once we are sure
          CALL iom_put( "iceamp"  , at_ip  * zswi        )   ! melt pond total fraction
          CALL iom_put( "icevmp"  , vt_ip  * zswi        )   ! melt pond total volume per unit area
@@ -232,5 +232,5 @@
 
       ! MV MP 2016
-      IF ( ln_limMP ) THEN
+      IF ( nn_pnd_scheme >= 0 ) THEN ! change once we are sure that case = 0 works
          IF ( iom_use( "iceamp_cat"  ) )  CALL iom_put( "iceamp_cat"     , a_ip   * zswi2   )       ! melt pond frac for categories
          IF ( iom_use( "icevmp_cat"  ) )  CALL iom_put( "icevmp_cat"     , v_ip   * zswi2   )       ! melt pond frac for categories
@@ -301,5 +301,5 @@
 
       ! MV MP 2016
-      IF ( ln_limMP ) THEN
+      IF ( nn_pnd_scheme >= 0 ) THEN
          CALL histdef( kid, "si_amp", "Melt pond fraction"      , "%"      ,   &
       &         jpi, jpj, kh_i, 1, 1, 1, -99, 32, "inst(x)", rdt, rdt )
@@ -354,5 +354,5 @@
 
       ! MV MP 2016
-      IF ( ln_limMP ) THEN
+      IF ( nn_pnd_scheme >= 0 ) THEN
          CALL histwrite( kid, "si_amp", kt, at_ip         , jpi*jpj, (/1/) )
          CALL histwrite( kid, "si_vmp", kt, vt_ip         , jpi*jpj, (/1/) )
