(**************************************************************************)
(*                                                                        *)
(* Proof of the Bresenham line drawing algorithm.                         *)
(*                                                                        *)
(* Jean-Christophe Filliâtre (LRI, Université Paris Sud)                  *)
(* May 2001                                                               *)
(*                                                                        *)
(**************************************************************************)

(*  Parameters. 
    Without loss of generality, we can take [x1=0] and [y1=0]. 
    Thus the line to draw joins [(0,0)] to [(x2,y2)]
    and we have [deltax = x2] and [deltay = y2]. 
    Moreover we assume being in the first octant, i.e.
    [0 <= y2 <= x2] (see the Coq file). The seven other cases can be easily
    deduced by symmetry. *)

logic x2,y2 : int

axiom first_octant : 0 <= y2 <= x2

(*s Global variables of the program. *)

parameter x,y,e : int ref

(*s The code.
    [(best x y)] expresses that the point [(x,y)] is the best
    possible point i.e. the closest to the real line (see the Coq file).
    The invariant relates [x], [y] and [e] and
    gives lower and upper bound for [e] (see the Coq file). *)

logic abs : int -> int

axiom abs_def : 
  forall x:int. (x >= 0 and abs(x) = x) or (x <= 0 and abs(x) = -x)

predicate best(x:int, y:int) =
  forall y':int. abs (x2 * y - x * y2) <= abs (x2 * y' - x * y2)

predicate Invariant(x:int, y:int, e:int) =
  e = 2 * (x + 1) * y2 - (2 * y + 1) * x2 and
  2 * (y2 - x2) <= e <= 2 * y2

axiom invariant_is_ok : forall x,y,e:int. Invariant(x,y,e) -> best(x,y)

let bresenham () =
  begin
    x := 0;
    y := 0;
    e := 2 * y2 - x2;
    while !x <= x2 do
      { invariant 0 <= x <= x2 + 1 and Invariant(x, y, e)
        variant   x2 + 1 - x }
      (* here we would do (plot x y) *) 
      assert { best(x, y) };
      if !e < 0 then
	e := !e + 2 * y2
      else begin
	y := !y + 1;
	e := !e + 2 * (y2 - x2)
      end;
      x := !x + 1
    done
  end

    
(*
Local Variables: 
compile-command: "gwhy-bin -split-user-conj bresenham.mlw"
End: 
*)