examples/verifythis_2017_maximum_sum_submatrix.mlw

   1
   2 (**
   3
   4 {1 VerifyThis @ ETAPS 2017 competition
   5    Challenge 2: Maximum-sum submatrix}
   6
   7 See https://formal.iti.kit.edu/ulbrich/verifythis2017/
   8
   9 Author: Jean-Christophe Filliâtre (CNRS)
  10 *)
  11
  12 (* note: this is a 2D-version of maximum-sum subarray, for which several
  13    verified implementations can be found in maximum_subarray.mlw,
  14    including Kadane's linear algorithm *)
  15
  16 module Kadane2D
  17
  18   use int.Int
  19   use ref.Refint
  20   use int.Sum
  21   use array.Array
  22   use matrix.Matrix
  23
  24   (* maximum-sum subarray problem is assumed *)
  25
  26   function array_sum (a: array int) (l h: int) : int
  27     = sum (fun i -> a[i]) l h
  28
  29   val maximum_subarray (a: array int) : (s: int, lo: int, hi: int)
  30     ensures { 0 <= lo <= hi <= length a /\ s = array_sum a lo hi /\
  31               forall l h. 0 <= l <= h <= length a -> s >= array_sum a l h }
  32
  33   (* sum of a submatrix *)
  34
  35   function col (m: matrix int) (j i: int) : int = m.elts i j
  36
  37   function matrix_sum (m: matrix int) (rl rh cl ch: int) : int
  38     = sum (fun j -> sum (col m j) rl rh) cl ch
  39
  40   let maximum_submatrix (m: matrix int) :
  41       (s: int, rlo: int, rhi: int, clo: int, chi: int)
  42     ensures { (* this is a legal submatrix *)
  43               0 <= rlo <= rhi <= rows m /\
  44               0 <= clo <= chi <= columns m /\
  45               (* s is its sum *)
  46               s = matrix_sum m rlo rhi clo chi /\
  47               (* and it is maximal *)
  48               forall rl rh. 0 <= rl <= rh <= rows m ->
  49               forall cl ch. 0 <= cl <= ch <= columns m ->
  50               s >= matrix_sum m rl rh cl ch }
  51   = let a = Array.make m.columns 0 in
  52     let maxsum = ref 0 in
  53     let rlo = ref 0 in
  54     let rhi = ref 0 in
  55     let clo = ref 0 in
  56     let chi = ref 0 in
  57     for rl = 0 to rows m - 1 do
  58       invariant { 0 <= !rlo <= !rhi <= rows m }
  59       invariant { 0 <= !clo <= !chi <= columns m }
  60       invariant { !maxsum = matrix_sum m !rlo !rhi !clo !chi >= 0 }
  61       invariant { forall rl' rh. 0 <= rl' < rl -> rl' <= rh <= rows m ->
  62                   forall cl ch. 0 <= cl <= ch <= columns m ->
  63                   !maxsum >= matrix_sum m rl' rh cl ch }
  64       fill a 0 (columns m) 0;
  65       for rh = rl + 1 to rows m do
  66         invariant { 0 <= !rlo <= !rhi <= rows m }
  67         invariant { 0 <= !clo <= !chi <= columns m }
  68         invariant { !maxsum = matrix_sum m !rlo !rhi !clo !chi >= 0 }
  69         invariant { forall rl' rh'. 0 <= rl' <= rh' <= rows m ->
  70                     (rl' < rl \/ rl' = rl /\ rh' < rh) ->
  71                     forall cl ch. 0 <= cl <= ch <= columns m ->
  72                     !maxsum >= matrix_sum m rl' rh' cl ch }
  73         invariant { forall j. 0 <= j < columns m ->
  74                     a[j] = sum (col m j) rl (rh - 1) }
  75         (* update array a *)
  76         for c = 0 to columns m -1 do
  77           invariant { forall j. 0 <= j < c ->
  78                       a[j] = sum (col m j) rl rh }
  79           invariant { forall j. c <= j < columns m ->
  80                       a[j] = sum (col m j) rl (rh - 1) }
  81           a[c] <- a[c] + get m (rh - 1) c
  82         done;
  83         (* then use Kadane algorithme on array a *)
  84         let sum, lo, hi = maximum_subarray a in
  85         assert { sum = matrix_sum m rl rh lo hi };
  86         assert { forall cl ch. 0 <= cl <= ch <= columns m ->
  87                  sum >= matrix_sum m rl rh cl ch
  88                  by array_sum a cl ch = matrix_sum m rl rh cl ch };
  89         (* update the maximum if needed *)
  90         if sum > !maxsum then begin
  91           maxsum := sum;
  92           rlo := rl; rhi := rh;
  93           clo := lo; chi := hi
  94         end
  95       done;
  96     done;
  97     !maxsum, !rlo, !rhi, !clo, !chi
  98
  99 end