Module AbLinearize

Require Import Utf8.
Require Import IdealExpr.
Require Import ZArith.
Require Import Coqlib.
Require Import ZIntervals FloatIntervals IdealIntervals.
Require Import Integers.
Require Import AdomLib.
Require Import AbIdealEnv.
Require Import Maps ShareTree Util TreeAl.
Require Import FloatLib.
Require Import LinearQuery.
Require Import FactMsgHelpers.

Transparent Float.add Float.sub Float.mul Float.div Float.neg Float.zero.

Local Existing Instance fitv_gamma_R.
Local Existing Instance fitv_gamma_oR.

Module Tp := Tree_Properties(T).
Module Tsp := ShareTree_Properties(T).

Definition t := unit.
Instance gamma_t : gamma_op t (var -> ideal_num) := fun _ _ => True.

Lemma list_norepet_rev:
  ∀ l, list_norepet (A:=var) l -> list_norepet (rev l).

Definition deletek {A} k (l:list (var * A)) :=
  List.filter (fun kv => negb (eq_dec (fst kv) k)) l.

Lemma deletek_eq:
  ∀ A k l, ~List.In k (List.map fst l) ->
         deletek (A:=A) k l = l.

Lemma deletek_in:
  ∀ A k l k'v',
    List.In k'v' (deletek (A:=A) k l) <->
    List.In k'v' l /\ (fst k'v') <> k.

Lemma deletek_norepet:
  ∀ A k l,
    list_norepet (List.map fst l) ->
    list_norepet (List.map fst (deletek (A:=A) k l)).

Lemma eval_linear_expr_exten:
  ∀ (P1 P2:R->Prop) ρ l1 l2,
    list_norepet (List.map fst l1) ->
    list_norepet (List.map fst l2) ->
    (∀ (k:var) v, List.In (k, v) l1 <-> List.In (k, v) l2) ->
  (P1 ⊆ P2) ->
  ∀ x,
    fold_right
      (fun (kitvs:var * (fitv*fitv)) X =>
         let '(k, (itvα, itvk)) := kitvs in
         fun y =>
           match rρ ρ k with
           | None => False
           | Some kx => (∃ α, α ∈ γ itvα /\ X (y - α * kx)) /\
                        rρ ρ k ∈ γ itvk
           end)
      P1 l1 x ->
    fold_right
      (fun (kitvs:var * (fitv*fitv)) X =>
         let '(k, (itvα, itvk)) := kitvs in
         fun y =>
           match rρ ρ k with
           | None => False
           | Some kx => (∃ α, α ∈ γ itvα /\ X (y - α * kx)) /\
                        rρ ρ k ∈ γ itvk
           end)
      P2 l2 x.

Lemma eval_T_linear_expr_exten:
  ∀ ρ i g1 g2,
    (∀ k, T.get k g1 = T.get k g2) ->
    ∀ x, eval_T_linear_expr ρ (i, g1) x <-> eval_T_linear_expr ρ (i, g2) x.

Lemma eval_T_linear_expr_set_None:
  ∀ ρ i g (k:var) kx,
  ∀ itvs α, α ∈ γ (fst itvs) ->
    rρ ρ k = Some kx ->
    rρ ρ k ∈ γ (snd itvs) ->
    T.get k g = None ->
  ∀ x, eval_T_linear_expr ρ (i, g) x ->
       eval_T_linear_expr ρ (i, T.set k itvs g) (x+α*kx).

Lemma eval_T_linear_expr_rem:
  ∀ ρ i g (k:var) kx, rρ ρ k = Some kx ->
  ∀ itvs, T.get k g = Some itvs ->
  ∀ x, eval_T_linear_expr ρ (i, g) x ->
  (∃ α, α ∈ γ (fst itvs) /\
        eval_T_linear_expr ρ (i, T.remove k g) (x-α*kx)) /\
  rρ ρ k ∈ γ (snd itvs).

Lemma eval_T_linear_expr_set_Some:
  ∀ ρ i g (k:var) kx, rρ ρ k = Some kx ->
  ∀ itvs, T.get k g = Some itvs ->
  ∀ x, eval_T_linear_expr ρ (i, g) x ->
  ∃ α0, α0 ∈ γ (fst itvs) /\
        ∀ itv1 itv1' α1, α1 ∈ γ itv1 -> rρ ρ k ∈ γ itv1' ->
        eval_T_linear_expr ρ (i, T.set k (itv1, itv1') g) (x+(α1-α0)*kx).

Lemma eval_linear_expr_change_cst:
  ∀ ρ i0 g,
  ∀ x, eval_linear_expr ρ (i0, g) x ->
  ∃ α0, α0 ∈ γ i0 /\
        ∀ i1 α1, α1 ∈ γ i1 ->
        eval_linear_expr ρ (i1, g) (x+α1-α0).

Definition add_T_linear_expr (e1 e2:T_linear_expr) : T_linear_expr :=
  (add_itv (fst e1) (fst e2),
   T.combine
     (fun ity1 ity2 =>
        match ity1, ity2 with
        | None, None => None
        | Some i, None => Some i
        | None, Some i => Some i
        | Some (i1, ty), Some (i2, _) => Some (add_itv i1 i2, ty)
        end)
     (snd e1) (snd e2)).

Definition add_T_linear_expr_alt (e1 e2:T_linear_expr) : T_linear_expr :=
  (add_itv (fst e1) (fst e2),
   T.fold
     (fun (g:T.t _) (k:var) ity1 =>
        match T.get k g with
        | None => T.set k ity1 g
        | Some (i2, _) => T.set k (add_itv (fst ity1) i2, snd ity1) g
        end) (snd e1) (snd e2)).

Lemma add_T_linear_expr_alt_eq:
  ∀ e1 e2 k,
    T.get k (snd (add_T_linear_expr_alt e1 e2)) =
    T.get k (snd (add_T_linear_expr e1 e2)).

Lemma add_T_linear_expr_alt_correct:
  ∀ ρ e1 e2 x1 x2,
    eval_T_linear_expr ρ e1 x1 ->
    eval_T_linear_expr ρ e2 x2 ->
    eval_T_linear_expr ρ (add_T_linear_expr_alt e1 e2) (x1+x2).

Lemma add_T_linear_expr_correct:
  ∀ ρ e1 e2 x1 x2,
    eval_T_linear_expr ρ e1 x1 ->
    eval_T_linear_expr ρ e2 x2 ->
    eval_T_linear_expr ρ (add_T_linear_expr e1 e2) (x1+x2).

Definition opp_T_linear_expr (e:T_linear_expr) : T_linear_expr :=
  (opp_itv (fst e), T.map1 (fun ity => (opp_itv (fst ity), snd ity)) (snd e)).

Lemma opp_T_linear_expr_correct:
  ∀ ρ e x,
    eval_T_linear_expr ρ e x ->
    eval_T_linear_expr ρ (opp_T_linear_expr e) (-x).

Definition mult_T_linear_expr (i:fitv) (e:T_linear_expr) : T_linear_expr :=
  (mult_itv i (fst e), T.map1 (fun itvs => (mult_itv i (fst itvs), snd itvs)) (snd e)).

Lemma mult_T_linear_expr_correct:
  ∀ ρ e x i y,
    eval_T_linear_expr ρ e x ->
    y ∈ γ i ->
    eval_T_linear_expr ρ (mult_T_linear_expr i e) (y * x).

Definition inv_itv (i:fitv) : fitv+⊤ :=
  let 'FITV a b := i in
  if Fleb (@B754_finite 53 1024 false 1 (-1074) eq_refl) a ||
     Fleb b (@B754_finite 53 1024 true 1 (-1074) eq_refl) then
    Just (FITV (Bdiv 53 1024 eq_refl eq_refl Float.binop_pl mode_DN one b)
               (Bdiv 53 1024 eq_refl eq_refl Float.binop_pl mode_UP one a))
  else All.

Lemma inv_itv_correct:
  ∀ i x, x ∈ γ i -> (/x) ∈ γ (inv_itv i).

Definition join_itv (i1 i2:fitv) : fitv :=
  match i1, i2 with
  | FITV m1 M1, FITV m2 M2 =>
    FITV (if Fleb m1 m2 || is_nan m2 then m1 else m2) (if Fleb M1 M2 || is_nan M1 then M2 else M1)
  end.

Lemma join_itv_correct:
  ∀ i1 i2 (x:R), x ∈ (γ i1 ∪ γ i2) -> x ∈ γ (join_itv i1 i2).

Definition double_mult_eps_itv : fitv :=
  FITV (@B754_finite 53 1024 false 9007199254740991 (-53) eq_refl)
       (@B754_finite 53 1024 false 4503599627370497 (-52) eq_refl).

Definition double_add_eps : float :=
  @B754_finite 53 1024 false 1 (-1074) eq_refl.

Definition double_add_eps_itv : fitv :=
  FITV (Bopp _ _ Float.neg_pl double_add_eps) double_add_eps.

Lemma double_round_NE_bounds:
  ∀ x,
    let r := round radix2 (Fcore_FLT.FLT_exp (3 - 1024 - 53) 53) (round_mode mode_NE) x in
    (∃ α, α ∈ γ double_mult_eps_itv /\ r = x*α) \/ (r-x) ∈ γ double_add_eps_itv.

Definition double_round_T_linear_expr (e:T_linear_expr) : T_linear_expr :=
  (join_itv (add_itv double_add_eps_itv (fst e))
            (mult_itv double_mult_eps_itv (fst e)),
   T.map1 (fun itvs => (mult_itv double_mult_eps_itv (fst itvs), snd itvs)) (snd e)).

Lemma double_round_T_linear_expr_correct:
  ∀ e ρ x, eval_T_linear_expr ρ e x ->
   let r := round radix2 (Fcore_FLT.FLT_exp (3 - 1024 - 53) 53) (round_mode mode_NE) x in
   eval_T_linear_expr ρ (double_round_T_linear_expr e) r.

Definition intervalize_linear_expr (e:linear_expr) : fitv :=
  List.fold_left
    (fun itv itvs =>
       let '(_, (itvα, itvk)) := itvs in
       add_itv (mult_itv itvα itvk) itv) (snd e) (fst e).

Lemma intervalize_linear_expr_correct :
  ∀ e ρ x,
    eval_linear_expr ρ e x ->
    x ∈ γ (intervalize_linear_expr e).

Definition intervalize_T_linear_expr (e:T_linear_expr) : fitv :=
  T.fold1 (fun itv itvs => add_itv (mult_itv (fst itvs) (snd itvs)) itv) (snd e) (fst e).

Lemma intervalize_T_linear_expr_correct :
  ∀ e ρ x,
    eval_T_linear_expr ρ e x ->
    x ∈ γ (intervalize_T_linear_expr e).

Definition double_round_over_T_linear_expr (e:T_linear_expr) : option T_linear_expr :=
  let le := double_round_T_linear_expr e in
  let 'FITV a b := intervalize_T_linear_expr le in
  if is_finite a && is_finite b then Some le else None.

Lemma double_round_over_T_linear_expr_correct:
  ∀ e ρ x, eval_T_linear_expr ρ e x ->
  let r := round radix2 (Fcore_FLT.FLT_exp (3 - 1024 - 53) 53) (round_mode mode_NE) x in
  ∀ e', double_round_over_T_linear_expr e = Some e' ->
  ∀ (f:float) sgn,
    (if Rlt_bool (Rabs r) (bpow radix2 1024) then (f:R) = r /\ is_finite f = true
     else (f:full_float) = binary_overflow 53 1024 mode_NE sgn) ->
  ∃ y, r_of_inum (INf f) = Some y /\ eval_T_linear_expr ρ e' y.

Definition single_mult_eps_itv : fitv :=
  FITV (@B754_finite 53 1024 false 9007198717870111 (-53) eq_refl)
       (@B754_finite 53 1024 false 4503599895805937 (-52) eq_refl).

Definition single_add_eps : float :=
  @B754_finite 53 1024 false 4503599627370496 (-202) eq_refl.

Definition single_add_eps_itv : fitv :=
  FITV (Bopp _ _ Float.neg_pl single_add_eps) single_add_eps.

Lemma single_round_NE_bounds:
  ∀ x,
    let r := round radix2 (Fcore_FLT.FLT_exp (3 - 128 - 24) 24) (round_mode mode_NE) x in
    (∃ α, α ∈ γ single_mult_eps_itv /\ r = x*α) \/
    (r-x) ∈ γ single_add_eps_itv.

Definition single_round_T_linear_expr (e:T_linear_expr) : T_linear_expr :=
  (join_itv (add_itv single_add_eps_itv (fst e))
            (mult_itv single_mult_eps_itv (fst e)),
   T.map1 (fun itvs => (mult_itv single_mult_eps_itv (fst itvs), snd itvs)) (snd e)).

Lemma single_round_T_linear_expr_correct:
  ∀ e ρ x, eval_T_linear_expr ρ e x ->
   let r := round radix2 (Fcore_FLT.FLT_exp (3 - 128 - 24) 24) (round_mode mode_NE) x in
   eval_T_linear_expr ρ (single_round_T_linear_expr e) r.

Definition single_ub : float :=
  @B754_finite 53 1024 false 9007198986305536 75 eq_refl.

Definition single_lb : float :=
  @B754_finite 53 1024 true 9007198986305536 75 eq_refl.

Definition single_round_over_T_linear_expr (e:T_linear_expr) : option T_linear_expr :=
  let 'FITV a b := intervalize_T_linear_expr e in
  if Fleb a single_lb || Fleb single_ub b then None
  else Some (single_round_T_linear_expr e).

Lemma single_round_over_T_linear_expr_correct:
  ∀ e ρ x, eval_T_linear_expr ρ e x ->
  let r := round radix2 (Fcore_FLT.FLT_exp (3 - 128 - 24) 24) (round_mode mode_NE) x in
  ∀ e', single_round_over_T_linear_expr e = Some e' ->
  ∀ (f:float32) sgn,
    (if Rlt_bool (Rabs r) (bpow radix2 128) then (f:R) = r /\ is_finite f = true
     else (f:full_float) = binary_overflow 24 128 mode_NE sgn) ->
  ∃ y, r_of_inum (INf (Float.of_single f)) = Some y /\ eval_T_linear_expr ρ e' y.

Definition ZofB_T_linear_expr (e:T_linear_expr) : T_linear_expr :=
  (add_itv (fst e) (FITV (Float.neg one) one), snd e).

Lemma Telements_empty:
  T.elements (T.empty (fitv*fitv)) = List.nil.

Lemma Telements_singleton:
  ∀ x y,
    T.elements (T.set x y (T.empty (fitv*fitv))) = (x,y)::List.nil.

Lemma ZofB_T_linear_expr_correct:
  ∀ e ρ x, eval_T_linear_expr ρ e x ->
   eval_T_linear_expr ρ (ZofB_T_linear_expr e) (Ztrunc x).

Fixpoint linearize_expr (c:query_chan) (e:iexpr) : option T_linear_expr :=
  match e with
  | IEvar _ x =>
    match fitv_of_iitv (c.(get_itv) e) with
    | All => None
    | Just (FITV a b as i) =>
      if is_finite a && is_finite b then
        Some (FITV (B754_zero false) (B754_zero false), T.set x (FITV one one, i) (T.empty _))
      else None
    end
  | IEconst (INz z) => Some (fitv_of_z z, T.empty _)
  | IEconst (INf f) =>
    if is_finite f then Some (FITV f f, T.empty _) else None
  | IEZitv itv => Some (fitv_of_zitv itv, T.empty _)
  | IEunop (IOneg | IOnegf) e' =>
    do le' <- linearize_expr c e';
    ret (opp_T_linear_expr le')
  | IEunop IOsingleoffloat e' =>
    do le' <- linearize_expr c e';
    single_round_over_T_linear_expr le'
  | IEunop IOfloatofz e' =>
    do le' <- linearize_expr c e';
    double_round_over_T_linear_expr le'
  | IEunop IOsingleofz e' =>
    do le' <- linearize_expr c e';
    single_round_over_T_linear_expr le'
  | IEunop IOzoffloat e' =>
    do le' <- linearize_expr c e';
    ret (ZofB_T_linear_expr le')
  | IEbinop IOadd e1 e2 =>
    do le1 <- linearize_expr c e1;
    do le2 <- linearize_expr c e2;
    ret (add_T_linear_expr le1 le2)
  | IEbinop IOsub e1 e2 =>
    do le1 <- linearize_expr c e1;
    do le2 <- linearize_expr c e2;
    ret (add_T_linear_expr le1 (opp_T_linear_expr le2))
  | IEbinop IOmul e1 e2 =>
    match e2 with
    | IEconst (INz z) =>
      do le' <- linearize_expr c e1;
      ret (mult_T_linear_expr (fitv_of_z z) le')
    | IEZitv i =>
      do le' <- linearize_expr c e1;
      ret (mult_T_linear_expr (fitv_of_zitv i) le')
    | _ =>
      match c.(get_itv) e1 with
      | NotBot (Just (Az i)) =>
        do le' <- linearize_expr c e2;
        ret (mult_T_linear_expr (fitv_of_zitv i) le')
      | _ => None
      end
    end
  | IEbinop IOdiv e1 e2 =>
    do le' <- linearize_expr c e1;
    match c.(get_itv) e2 with
    | NotBot (Just (Az i)) =>
      match inv_itv (fitv_of_zitv i) with
      | Just i' => Some (ZofB_T_linear_expr (mult_T_linear_expr i' le'))
      | All => None
      end
    | _ => None
    end
  | IEbinop IOaddf e1 e2 =>
    do le1 <- linearize_expr c e1;
    do le2 <- linearize_expr c e2;
    double_round_over_T_linear_expr (add_T_linear_expr le1 le2)
  | IEbinop IOsubf e1 e2 =>
    do le1 <- linearize_expr c e1;
    do le2 <- linearize_expr c e2;
    double_round_over_T_linear_expr (add_T_linear_expr le1 (opp_T_linear_expr le2))
  | IEbinop IOmulf e1 e2 =>
    match e2 with
    | IEconst (INf f) =>
      if is_finite f then
        do le' <- linearize_expr c e1;
        double_round_over_T_linear_expr (mult_T_linear_expr (FITV f f) le')
      else
        None
    | _ =>
      match fitv_of_iitv (c.(get_itv) e1) with
      | Just (FITV a b as i) =>
        do le' <- linearize_expr c e2;
        double_round_over_T_linear_expr (mult_T_linear_expr i le')
      | _ => None
      end
    end
  | IEbinop IOdivf e1 e2 =>
    do le' <- linearize_expr c e1;
    match fitv_of_iitv (c.(get_itv) e2) with
    | Just (FITV a b as i) =>
      match inv_itv i with
      | Just i' => double_round_over_T_linear_expr (mult_T_linear_expr i' le')
      | All => None
      end
    | _ => None
    end
  | _ =>
    match fitv_of_iitv (c.(get_itv) e) with
    | Just i => Some (i, T.empty _)
    | _ => None
    end
  end.

Lemma linearize_expr_correct:
  ∀ c e le, linearize_expr c e = Some le ->
  ∀ ρ x, ρ ∈ γ c -> eval_iexpr ρ e x ->
  ∃ y, r_of_inum x = Some y /\ eval_T_linear_expr ρ le y.

Definition double_min_delta : float :=
  @B754_finite 53 1024 false 1 (-1074) eq_refl.

Lemma double_min_delta_correct :
  ∀ (f1 f2:float),
    is_finite f1 = true -> is_finite f2 = true ->
    (f1 < f2) -> (f1 - f2 <= (Float.neg double_min_delta)).

Definition T_linear_expr_is_top (le:T_linear_expr) : bool :=
  T.fold (fun acc x itvs => acc || fitv_is_top (fst itvs)) (snd le) false.

Definition T_linear_expr_is_nonzero (le:T_linear_expr) : bool :=
  let 'FITV a b := intervalize_T_linear_expr le in
  negb (Fleb a Float.zero && Fleb Float.zero b).

Lemma T_linear_expr_is_nonzero_correct:
  ∀ ρ le, eval_T_linear_expr ρ le 0 -> T_linear_expr_is_nonzero le = false.

Definition process_msg (m:message) (ab:unit * query_chan) :=
  let c := snd ab in
  let new_msg :=
    match m with
    | Fact_msg e b =>
      match simplify_bool_expr e c with
      | (IEbinop ((IOcmp _ | IOcmpf _) as op) e1 e2, b') =>
        do le1 <- linearize_expr c e1;
        do le2 <- linearize_expr c e2;
        let le := add_T_linear_expr le1 (opp_T_linear_expr le2) in
        do i <-
          match op, xorb b b' with
          | IOcmp Ceq, true | IOcmp Cne, false
          | IOcmpf Ceq, true | IOcmpf Cne, false =>
            ret (FITV Float.zero Float.zero)

          | (IOcmp Clt|IOcmpf Clt), false | (IOcmp Cge|IOcmpf Cge), true =>
            ret (FITV (B754_infinity true) Float.zero)
          | IOcmp Cle, false | IOcmp Cgt, true =>
            ret (FITV (B754_infinity true) (Float.neg one))
          | IOcmpf Cle, false | IOcmpf Cgt, true =>
            ret (FITV (B754_infinity true) (Float.neg double_min_delta))

          | (IOcmp Cgt|IOcmpf Cgt), false | (IOcmp Cle|IOcmpf Cle), true =>
            ret (FITV Float.zero (B754_infinity false))
          | IOcmp Cge, false | IOcmp Clt, true =>
            ret (FITV one (B754_infinity false))
          | IOcmpf Cge, false | IOcmpf Clt, true =>
            ret (FITV double_min_delta (B754_infinity false))

          | _, _ => None
          end;
          ret (add_T_linear_expr le (i, T.empty _))
      | _ => None
      end
    | Known_value_msg x v =>
      do itv <-
        match v with
        | INz z => ret (fitv_of_zitv (ZITV z z))
        | INf f => if is_finite f then ret (FITV f f) else None
        end;
      ret (opp_itv itv, T.set x (FITV one one, itv) (T.empty _))
    | Itv_msg x i =>
      match fitv_of_iitv (NotBot (Just i)) with
      | All => None
      | Just itv =>
        ret (opp_itv itv, T.set x (FITV one one, itv) (T.empty _))
      end
    | _ => None
    end
  in
  match new_msg with
  | None => NotBot (tt, m::List.nil)
  | Some m' =>
    if T_linear_expr_is_top m' then NotBot (tt, m::List.nil)
    else
      if T_linear_expr_is_nonzero m' then Bot
      else NotBot (tt, m::(Linear_zero_msg m')::List.nil)
  end.

Transparent Float.cmp.

Lemma process_msg_correct:
  forall m ρ ab,
    ρ ∈ γ ab -> ρ ∈ γ m ->
    ρ ∈ γ (process_msg m ab).

Program Instance abVarExpEq_env : ab_ideal_env t t :=
  {| id_leb := fun _ _ => true;
     id_top := fun _ => NotBot (tt, List.nil);
     id_join := fun _ _ _ => NotBot (tt, List.nil);
     id_init_iter := fun _ => tt;
     id_widen := fun _ _ _ => (tt, NotBot (tt, List.nil));
     assign := fun _ _ _ _ => NotBot (tt, List.nil);
     forget := fun _ _ _ => NotBot (tt, List.nil);
     enrich_query_chan :=
       fun _ c =>
         {| get_var_ty := c.(get_var_ty);
            get_itv := c.(get_itv);
            get_congr := c.(get_congr);
            get_eq_expr := c.(get_eq_expr);
            AbIdealEnv.linearize_expr e :=
              match linearize_expr c e with
              | None => None
              | Some le =>
                if T_linear_expr_is_top le then None else Some le
              end
         |};
     process_msg := process_msg |}.

Opaque t.