Index: /NEMO/branches/2020/SI3-05_MP/src/ICE/icethd_pnd.F90
===================================================================
--- /NEMO/branches/2020/SI3-05_MP/src/ICE/icethd_pnd.F90	(revision 13810)
+++ /NEMO/branches/2020/SI3-05_MP/src/ICE/icethd_pnd.F90	(revision 13811)
@@ -53,4 +53,15 @@
       !!               
       !! ** Purpose :   change melt pond fraction and thickness
+      !!
+      !! Note: Melt ponds affect currently affect radiative transfer 
+      !!       They carry no heat, and the melt water they carry is not
+      !!       exchanged with the ocean
+      !!
+      !! 						This means freshwater is directly released after surface & bottom melt in ice_thd_dh
+      !!
+      !!       A wfx_pnd has been coded for diagnostic purposes
+      !!       Each time wfx_pnd is updated, wfx_sum / wfx_snw_sum must be updated
+      !!
+      !!!      The current diagnostic lacks a contribution from drainage
       !!                
       !!-------------------------------------------------------------------
@@ -226,6 +237,5 @@
             !--- mass flux ---!
             IF( zdv_mlt > 0._wp ) THEN
-               ! MV add comment on what that mass flux means
-               ! water removed from fw flux due to melt pond growth
+             ! MV add comment on what that mass flux means
                zfac = zdv_mlt * rhow * r1_rdtice                        ! melt pond mass flux < 0 [kg.m-2.s-1]
                wfx_pnd_1d(ji) = wfx_pnd_1d(ji) - zfac
@@ -277,5 +287,4 @@
             
             ! Calculate the permeability of the ice (Assur 1958, see Flocco 2010)
-            ! zsbr = - ( t_i_1d(ji,jk) - rt0 ) / rTmlt
             DO jk = 1, nlay_i
                zsbr = - 1.2_wp                                  &
@@ -293,4 +302,6 @@
             zdv_flush   = MAX( zdv_flush, -v_ip_1d(ji) )
             v_ip_1d(ji) = v_ip_1d(ji) + zdv_flush
+            
+            ! MV --- why pond drainage does not give back water into freshwater flux ?
             
             !--- Set new pond area and depth ---! assuming linear relation between h_ip and a_ip_frac
@@ -407,4 +418,7 @@
 
        REAL (wp), PARAMETER :: &
+          ! MV ouate de phoque!!! constants hard coded ????
+          ! 917 = rhoice
+          ! 0.334 =Lfus!!!!
           zr1_rlfus = 1._wp / 0.334e+6 / 917._wp , & ! (J/m^3)
           zTp       = -0.15_wp, & ! pond freezing temperature (C)
@@ -474,8 +488,8 @@
 
       !js 23/04/19: rewired back to a fraction with a_i
-      zrfrac(:,:,:) = rn_pnd_fracmin + ( rn_pnd_fracmax - rn_pnd_fracmin ) * aicen(:,:,:)
+      !!! zrfrac(:,:,:) = rn_pnd_fracmin + ( rn_pnd_fracmax - rn_pnd_fracmin ) * aicen(:,:,:)
+      ! MV2020
+      zrfrac(:,:,:) = rn_apnd_min + ( rn_apnd_max - rn_apnd_min ) * at_i_1d(:) !  = ( 1 - r ) = fraction of melt water that is not flushed
       zwfx_tmp(:,:) = 0._wp
-
-! MV ---> use expression from level-ice ponds, clarify what is needed
 
       !--- Add retained melt water to melt ponds
@@ -496,8 +510,15 @@
                                     * ( dh_i_pnd(:,:,jl)*rhoic + dh_s_pnd(:,:,jl)*rhosn )          & 
                                     * z1_rhow * a_i(:,:,jl)
+! MV ---> use expression from level ice melt ponds (dv_mlt)
+                                            
       END DO
 
-! MV ---> use expression from level ice melt ponds (dv_mlt)
-
+
+!!! MV ---> rewrite this accounting for level ice melt ponds 
+!!! 1) wfx_pnd_in is obsolete, ln_pnd_fw as well
+!!! 2) 
+!!! wfx_sum
+      !!!
+      !!!
       !js 03/05/19: we truncate negative values after calculating zvolp, in a
       ! similar manner to the subroutine ice_thd_pnd_cesm. Variation dh_i_pnd and
@@ -531,4 +552,6 @@
          wfx_snw_sum(:,:) = ( 1._wp - zwfx_tmp(:,:) ) * wfx_snw_sum(:,:)
       ENDIF
+      
+      
 
        !-----------------------------------------------------------------
@@ -568,4 +591,7 @@
       END DO                     ! jj
 
+
+
+!!! MV sounds like this should be replaced by lid
 
        DO ij = 1, icells
