Theory QElin_inf
theory QElin_inf
imports LinArith
begin
subsection ‹Quantifier elimination with infinitesimals \label{sec:lin-inf}›
text‹This section formalizes Loos and Weispfenning's quantifier
elimination procedure based on (the simulation of)
infinitesimals~\<^cite>‹"LoosW93"›.›
fun asubst_peps :: "real * real list ⇒ atom ⇒ atom fm" (‹asubst⇩+›) where
"asubst_peps (r,cs) (Less s (d#ds)) =
  (if d=0 then Atom(Less s ds) else
   let u = s - d*r; v = d *⇩s cs + ds; less = Atom(Less u v)
   in if d<0 then less else Or less (Atom(Eq u v)))" |
"asubst_peps rcs (Eq r (d#ds)) = (if d=0 then Atom(Eq r ds) else FalseF)" |
"asubst_peps rcs a = Atom a"
abbreviation subst_peps :: "atom fm ⇒ real * real list ⇒ atom fm" (‹subst⇩+›)
where "subst⇩+ φ rcs ≡ amap⇩f⇩m (asubst⇩+ rcs) φ"
definition "nolb f xs l x = (∀y∈{l<..<x}. y ∉ LB f xs)"
lemma nolb_And[simp]:
  "nolb (And f g) xs l x = (nolb f xs l x ∧ nolb g xs l x)"
apply(clarsimp simp:nolb_def)
apply blast
done
lemma nolb_Or[simp]:
  "nolb (Or f g) xs l x = (nolb f xs l x ∧ nolb g xs l x)"
apply(clarsimp simp:nolb_def)
apply blast
done
context notes [[simp_depth_limit=4]]
begin
lemma innermost_intvl:
  "⟦ nqfree f; nolb f xs l x; l < x; x ∉ EQ f xs; R.I f (x#xs); l < y; y ≤ x⟧
  ⟹ R.I f (y#xs)"
proof(induct f)
  case (Atom a)
  show ?case
  proof (cases a)
    case [simp]: (Less r cs)
    show ?thesis
    proof (cases cs)
      case Nil thus ?thesis using Atom by (simp add:depends⇩R_def)
    next
      case [simp]: (Cons c cs)
      hence "r < c*x + ⟨cs,xs⟩" using Atom by simp
      { assume "c=0" hence ?thesis using Atom by simp }
      moreover
      { assume "c<0"
        hence "x < (r - ⟨cs,xs⟩)/c" (is "_ < ?u") using ‹r < c*x + ⟨cs,xs⟩›
          by (simp add: field_simps)
        have ?thesis
        proof (rule ccontr)
          assume "¬ R.I (Atom a) (y#xs)"
          hence "?u ≤ y" using Atom ‹c<0›
            by (auto simp add: field_simps)
          with ‹x<?u› show False using Atom ‹c<0›
            by(auto simp:depends⇩R_def)
        qed } moreover
      { assume "c>0"
        hence "x > (r - ⟨cs,xs⟩)/c" (is "_ > ?l") using ‹r < c*x + ⟨cs,xs⟩›
          by (simp add: field_simps)
        then have "?l < y" using Atom ‹c>0›
            by (auto simp:depends⇩R_def Ball_def nolb_def)
               (metis linorder_not_le antisym order_less_trans)
        hence ?thesis using ‹c>0› by (simp add: field_simps)
      } ultimately show ?thesis by force
    qed
  next
    case [simp]: (Eq r cs)
    show ?thesis
    proof (cases cs)
      case Nil thus ?thesis using Atom by (simp add:depends⇩R_def)
    next
      case [simp]: (Cons c cs)
      hence "r = c*x + ⟨cs,xs⟩" using Atom by simp
      { assume "c=0" hence ?thesis using Atom by simp }
      moreover
      { assume "c≠0"
        hence ?thesis using ‹r = c*x + ⟨cs,xs⟩› Atom
          by (auto simp: depends⇩R_def split: if_splits) }
      ultimately show ?thesis by force
    qed
  qed
next
  case (And f1 f2) thus ?case by (fastforce)
next
  case (Or f1 f2) thus ?case by (fastforce)
qed simp+
definition "EQ2 = EQ"
lemma EQ2_Or[simp]: "EQ2 (Or f g) xs = (EQ2 f xs ∪ EQ2 g xs)"
by(auto simp:EQ2_def)
lemma EQ2_And[simp]: "EQ2 (And f g) xs = (EQ2 f xs ∪ EQ2 g xs)"
by(auto simp:EQ2_def)
lemma innermost_intvl2:
  "⟦ nqfree f; nolb f xs l x; l < x; x ∉ EQ2 f xs; R.I f (x#xs); l < y; y ≤ x⟧
  ⟹ R.I f (y#xs)"
unfolding EQ2_def by(blast intro:innermost_intvl)
lemma I_subst_peps2:
 "nqfree f ⟹ r+⟨cs,xs⟩ < x ⟹ nolb f xs (r+⟨cs,xs⟩) x
 ⟹ ∀y ∈ {r+⟨cs,xs⟩ <.. x}. R.I f (y#xs) ∧ y ∉ EQ2 f xs
 ⟹ R.I (subst⇩+ f (r,cs)) xs"
proof(induct f)
  case FalseF thus ?case
    by simp (metis antisym_conv1 linorder_neq_iff)
next
  case (Atom a)
  show ?case
  proof(cases "((r,cs),a)" rule:asubst_peps.cases)
    case (1 r cs s d ds)
    { assume "d=0" hence ?thesis using Atom 1 by auto }
    moreover
    { assume "d<0"
      have "s < d*x + ⟨ds,xs⟩" using Atom 1 by simp
      moreover have "d*x < d*(r + ⟨cs,xs⟩)" using ‹d<0› Atom 1
        by (simp add: mult_strict_left_mono_neg)
      ultimately have "s < d * (r + ⟨cs,xs⟩) + ⟨ds,xs⟩" by(simp add:algebra_simps)
      hence ?thesis using 1
        by (auto simp add: iprod_left_add_distrib algebra_simps)
    } moreover
    { let ?L = "(s - ⟨ds,xs⟩) / d" let ?U = "r + ⟨cs,xs⟩"
      assume "d>0"
      hence "?U < x" and "∀y. ?U < y ∧ y < x ⟶ y ≠ ?L"
        and "∀y. ?U < y ∧ y ≤ x ⟶ ?L < y" using Atom 1
        by(simp_all add:nolb_def depends⇩R_def Ball_def field_simps)
      hence "?L < ?U ∨ ?L = ?U"
        by (metis linorder_neqE_linordered_idom order_refl)
      hence ?thesis using Atom 1 ‹d>0›
        by (simp add: iprod_left_add_distrib field_simps)
    } ultimately show ?thesis by force
  next
    case 2 thus ?thesis using Atom
      by (fastforce simp: nolb_def EQ2_def depends⇩R_def field_simps split: if_split_asm)
  qed (insert Atom, auto)
next
  case Or thus ?case by(simp add:Ball_def)(metis order_refl innermost_intvl2)
qed simp_all
end
lemma I_subst_peps:
  "nqfree f ⟹ R.I (subst⇩+ f (r,cs)) xs ⟹
  (∃leps>r+⟨cs,xs⟩. ∀x. r+⟨cs,xs⟩ < x ∧ x ≤ leps ⟶ R.I f (x#xs))"
proof(induct f)
  case TrueF thus ?case by simp (metis less_add_one)
next
  case (Atom a)
  show ?case
  proof (cases "((r,cs),a)" rule: asubst_peps.cases)
    case (1 r cs s d ds)
    { assume "d=0" hence ?thesis using Atom 1 by auto (metis less_add_one) }
    moreover
    { assume "d<0"
      with Atom 1 have "r + ⟨cs,xs⟩ < (s - ⟨ds,xs⟩)/d" (is "?a < ?b")
        by(simp add:field_simps iprod_left_add_distrib)
      then obtain x where "?a < x" "x < ?b" by(metis dense)
      hence " ∀y. ?a < y ∧ y ≤ x ⟶ s < d*y + ⟨ds,xs⟩"
        using ‹d<0› by (simp add:field_simps)
      (metis add_le_cancel_right mult_le_cancel_left order_antisym linear mult.commute xt1(8))
      hence ?thesis using 1 ‹?a<x› by auto
    } moreover
    { let ?a = "s - d * r" let ?b = "⟨d *⇩s cs + ds,xs⟩"
      assume "d>0"
      with Atom 1 have "?a < ?b ∨ ?a = ?b" by auto
      hence ?thesis
      proof
        assume "?a = ?b"
        thus ?thesis using ‹d>0› Atom 1
          by(simp add:field_simps iprod_left_add_distrib)
            (metis add_0_left add_less_cancel_right distrib_left mult.commute mult_strict_left_mono)
      next
        assume "?a < ?b"
        { fix x assume "r+⟨cs,xs⟩ < x ∧ x ≤ r+⟨cs,xs⟩ + 1"
          hence "d*(r + ⟨cs,xs⟩) < d*x"
            using ‹d>0› by(metis mult_strict_left_mono)
          hence "s < d*x + ⟨ds,xs⟩" using ‹d>0› ‹?a < ?b›
            by (simp add:algebra_simps iprod_left_add_distrib)
        }
        thus ?thesis using 1 ‹d>0›
          by(force simp: iprod_left_add_distrib)
      qed
    } ultimately show ?thesis by (metis less_linear)
  qed (insert Atom, auto split:if_split_asm intro: less_add_one)
next
  case And thus ?case
    apply clarsimp
    apply(rule_tac x="min leps lepsa" in exI)
    apply simp
    done
next
  case Or thus ?case by force
qed simp_all
lemma dense_interval:
assumes "finite L" "l ∈ L" "l < x" "P(x::real)"
and dense: "⋀y l. ⟦ ∀y∈{l<..<x}. y ∉ L; l<x; l<y; y≤x ⟧ ⟹ P y"
shows "∃l∈L. l<x ∧ (∀y∈{l<..<x}. y ∉ L) ∧ (∀y. l<y ∧ y≤x ⟶ P y)"
proof -
  let ?L = "{l∈L. l < x}"
  let ?ll = "Max ?L"
  have "?L ≠ {}" using ‹l ∈ L› ‹l<x› by (blast intro:order_less_imp_le)
  hence "∀y. ?ll<y ∧ y<x ⟶ y ∉ L" using ‹finite L› by force
  moreover have "?ll ∈ L"
  proof
    show "?ll ∈ ?L" using ‹finite L› Max_in[OF _ ‹?L ≠ {}›] by simp
    show "?L ⊆ L" by blast
  qed
  moreover have "?ll < x" using ‹finite L› ‹?L ≠ {}› by simp
  ultimately show ?thesis using ‹l < x› ‹?L ≠ {}›
    by(blast intro!:dense greaterThanLessThan_iff[THEN iffD1])
qed
definition
"qe_eps⇩1(f) =
(let as = R.atoms⇩0 f; lbs = lbounds as; ebs = ebounds as
 in list_disj (inf⇩- f # map (subst⇩+ f) lbs @ map (subst f) ebs))"
theorem I_eps1:
assumes "nqfree f" shows "R.I (qe_eps⇩1 f) xs = (∃x. R.I f (x#xs))"
  (is "?QE = ?EX")
proof
  let ?as = "R.atoms⇩0 f" let ?ebs = "ebounds ?as"
  assume ?QE
  { assume "R.I (inf⇩- f) xs"
    hence ?EX using ‹?QE› min_inf[of f xs] ‹nqfree f›
      by(auto simp add:qe_eps⇩1_def amap_fm_list_disj)
  } moreover
  { assume "∀x ∈ EQ f xs. ¬R.I f (x#xs)"
           "¬ R.I (inf⇩- f) xs"
    with ‹?QE› ‹nqfree f› obtain r cs where "R.I (subst⇩+ f (r,cs)) xs"
      by (fastforce simp: qe_eps⇩1_def set_ebounds diff_divide_distrib eval_def I_subst ‹nqfree f›)
    then obtain leps where "R.I f (leps#xs)"
      using I_subst_peps[OF ‹nqfree f›] by fastforce
    hence ?EX .. }
  ultimately show ?EX by blast
next
  let ?as = "R.atoms⇩0 f" let ?ebs = "ebounds ?as"
  assume ?EX
  then obtain x where x: "R.I f (x#xs)" ..
  { assume "R.I (inf⇩- f) xs"
    hence ?QE using ‹nqfree f› by(auto simp:qe_eps⇩1_def)
  } moreover
  { assume "∃rcs ∈ set ?ebs. R.I (subst f rcs) xs"
    hence ?QE by(auto simp:qe_eps⇩1_def) } moreover
  { assume "¬ R.I (inf⇩- f) xs"
    and "∀rcs ∈ set ?ebs. ¬ R.I (subst f rcs) xs"
    hence noE: "∀e ∈ EQ f xs. ¬ R.I f (e#xs)" using ‹nqfree f›
      by (force simp:set_ebounds I_subst diff_divide_distrib eval_def split:if_split_asm)
    hence "x ∉ EQ f xs" using x by fastforce
    obtain l where "l ∈ LB f xs" "l < x"
      using LBex[OF ‹nqfree f› x ‹¬ R.I(inf⇩- f) xs› ‹x ∉ EQ f xs›] ..
    have "∃l∈LB f xs. l<x ∧ nolb f xs l x ∧
            (∀y. l < y ∧ y ≤ x ⟶ R.I f (y#xs))"
      using dense_interval[where P = "λx. R.I f (x#xs)", OF finite_LB ‹l∈LB f xs› ‹l<x› x] x innermost_intvl[OF ‹nqfree f› _ _ ‹x ∉ EQ f xs›]
      by (simp add:nolb_def)
    then obtain r c cs
      where *: "Less r (c#cs) ∈ set(R.atoms⇩0 f) ∧ c>0 ∧
            (r - ⟨cs,xs⟩)/c < x ∧ nolb f xs ((r - ⟨cs,xs⟩)/c) x
            ∧ (∀y. (r - ⟨cs,xs⟩)/c < y ∧ y ≤ x ⟶ R.I f (y#xs))"
      by blast
    then have "R.I (subst⇩+ f (r/c, (-1/c) *⇩s cs)) xs" using noE
      by(auto intro!: I_subst_peps2[OF ‹nqfree f›]
        simp:EQ2_def diff_divide_distrib algebra_simps)
    with * have ?QE
      by(simp add:qe_eps⇩1_def bex_Un set_lbounds) metis
  } ultimately show ?QE by blast
qed
lemma qfree_asubst_peps: "qfree (asubst⇩+ rcs a)"
by(cases "(rcs,a)" rule:asubst_peps.cases) simp_all
lemma qfree_subst_peps: "nqfree φ ⟹ qfree (subst⇩+ φ rcs)"
by(induct φ) (simp_all add:qfree_asubst_peps)
lemma qfree_qe_eps⇩1: "nqfree φ ⟹ qfree(qe_eps⇩1 φ)"
apply(simp add:qe_eps⇩1_def)
apply(rule qfree_list_disj)
apply (auto simp:qfree_min_inf qfree_subst_peps qfree_map_fm)
done
definition "qe_eps = R.lift_nnf_qe qe_eps⇩1"
lemma qfree_qe_eps: "qfree(qe_eps φ)"
by(simp add: qe_eps_def R.qfree_lift_nnf_qe qfree_qe_eps⇩1)
lemma I_qe_eps: "R.I (qe_eps φ) xs = R.I φ xs"
by(simp add:qe_eps_def R.I_lift_nnf_qe qfree_qe_eps⇩1 I_eps1)
end