Theory BenOr_Kozen_Reif.BKR_Decision

theory BKR_Decision
  imports BKR_Algorithm
    "Berlekamp_Zassenhaus.Factorize_Rat_Poly"
    "Algebraic_Numbers.Real_Roots"
     BKR_Proofs
    "HOL.Deriv"
begin

section "Algorithm"

subsection "Parsing"
  (* Formula type *)
datatype 'a fml =
  And "'a fml" "'a fml"
  | Or "'a fml" "'a fml"
  | Gt 'a   (* 'a > 0 *)
  | Geq 'a  (* 'a ≥ 0 *)
  | Lt 'a   (* 'a < 0 *)
  | Leq 'a  (* 'a ≤ 0 *)
  | Eq 'a   (* 'a = 0 *)
  | Neq 'a  (* 'a ≠ 0 *)

(* Evaluating a formula over a lookup semantics where 'a is nat *)
primrec lookup_sem :: "nat fml ⇒ ('a::linordered_field list) ⇒ bool"
  where
    "lookup_sem (And l r) ls = (lookup_sem l ls ∧ lookup_sem r ls)"
  | "lookup_sem (Or l r) ls = (lookup_sem l ls ∨ lookup_sem r ls)"
  | "lookup_sem (Gt p) ls = (ls ! p > 0)"  
  | "lookup_sem (Geq p) ls = (ls ! p ≥ 0)" 
  | "lookup_sem (Lt p) ls = (ls ! p < 0)" 
  | "lookup_sem (Leq p) ls = (ls ! p ≤ 0)" 
  | "lookup_sem (Eq p) ls = (ls ! p = 0)" 
  | "lookup_sem (Neq p) ls = (ls ! p ≠ 0)" 

(* (compute) all polynomials mentioned in a formula *)
primrec poly_list :: "'a fml ⇒ 'a list"
  where
    "poly_list (And l r) = poly_list l @ poly_list r"  
  | "poly_list (Or l r) = poly_list l @ poly_list r"  
  | "poly_list (Gt p) = [p]"
  | "poly_list (Geq p) = [p]"
  | "poly_list (Lt p) = [p]"
  | "poly_list (Leq p) = [p]"
  | "poly_list (Eq p) = [p]"
  | "poly_list (Neq p) = [p]"

primrec index_of_aux :: "'a list ⇒ 'a ⇒ nat ⇒ nat" where
  "index_of_aux [] y n = n"
| "index_of_aux (x#xs) y n =
  (if x = y then n else index_of_aux xs y (n+1))"

definition index_of :: "'a list ⇒ 'a ⇒ nat" where
  "index_of xs y = index_of_aux xs y 0"

definition convert :: "'a fml ⇒ (nat fml × 'a list)"
  where
    "convert fml = (
    let ps = remdups (poly_list fml)
    in
    (map_fml (index_of ps) fml, ps)
  )"


subsection "Factoring"

(* Makes sure the result of factorize_rat_poly is monic *)
definition factorize_rat_poly_monic :: "rat poly ⇒ (rat × (rat poly × nat) list)"
  where
    "factorize_rat_poly_monic p = (
    let (c,fs) =  factorize_rat_poly p ;
        lcs = prod_list (map (λ(f,i). (lead_coeff f) ^ Suc i) fs) ;
        fs = map (λ(f,i). (normalize f, i)) fs
    in
    (c * lcs,fs)
  )"

(* Factoring an input list of polynomials *)
definition factorize_polys :: "rat poly list ⇒ (rat poly list × (rat × (nat × nat) list) list)"
  where
    "factorize_polys ps = (
    let fact_ps = map factorize_rat_poly_monic ps;
        factors = remdups (map fst (concat (map snd fact_ps))) ;
        data = map (λ(c,fs). (c, map (λ(f,pow). (index_of factors f, pow) ) fs)) fact_ps
    in
      (factors,data)
  )"

(* After turning a polynomial into factors,
  this turns a sign condition on the factors
  into a sign condition for the polynomial *)
definition undo_factorize :: "rat × (nat × nat) list ⇒ rat list ⇒ rat"
  where
    "undo_factorize cfs signs =
   squash
    (case cfs of (c,fs) ⇒
    (c * prod_list (map (λ(f,pow). (signs ! f) ^ Suc pow) fs)))
  "

definition undo_factorize_polys :: "(rat × (nat × nat) list) list ⇒ rat list ⇒ rat list"
  where
    "undo_factorize_polys ls signs = map (λl. undo_factorize l signs) ls"

subsection "Auxiliary Polynomial"
definition crb:: "real poly ⇒ int" where
  "crb p = ceiling (2 + max_list_non_empty (map (λ i. norm (coeff p i)) [0 ..< degree p]) 
    / norm (lead_coeff p))"

(* Because we are using prod_list instead of lcm, it's important that this is called 
    when ps is pairwise coprime. *)
definition coprime_r :: "real poly list ⇒ real poly"
  where
    "coprime_r ps = pderiv (prod_list ps) * ([:-(crb (prod_list ps)),1:]) * ([:(crb (prod_list ps)),1:])"


subsection "Setting Up the Procedure"
  (* 0 indexed *)
definition insertAt :: "nat ⇒ 'a  ⇒ 'a list ⇒ 'a list" where
  "insertAt n x ls = take n ls @ x # (drop n ls)"

(* 0 indexed *)
definition removeAt :: "nat ⇒ 'a list ⇒ 'a list" where
  "removeAt n ls = take n ls @ (drop (n+1) ls)"

definition find_sgas_aux:: "real poly list ⇒ rat list list"
  where "find_sgas_aux in_list =
  concat (map (λi.
    map (λv. insertAt i 0 v) (find_consistent_signs_at_roots (in_list ! i) (removeAt i in_list))
  ) [0..<length in_list])"

(* For an input list of real polynomials, apply BKR to all positions *)
definition find_sgas :: "real poly list ⇒ rat list list"
  where
    "find_sgas ps = ( 
    let r = coprime_r ps in
    find_consistent_signs_at_roots r ps @ find_sgas_aux ps
  )"

(* Putting the sign condition preprocessing together with BKR *)
definition find_consistent_signs :: "rat poly list ⇒ rat list list"
  where
    "find_consistent_signs ps = (
     let (fs,data) = factorize_polys ps;
         sgas = find_sgas (map (map_poly of_rat) fs);
         rsgas = map (undo_factorize_polys data) sgas
    in
    (if fs = [] then [(map (λx. if poly x 0 < 0 then -1 else if poly x 0 = 0 then 0 else 1) ps)] else rsgas)
  )"


subsection "Deciding Univariate Problems"
definition decide_universal :: "rat poly fml ⇒ bool"
  where [code]:
    "decide_universal fml = (
    let (fml_struct,polys) = convert fml;
        conds = find_consistent_signs polys
    in
     list_all (lookup_sem fml_struct) conds
  )"

definition decide_existential :: "rat poly fml ⇒ bool"
  where [code]:
    "decide_existential fml = (
    let (fml_struct,polys) = convert fml;
        conds = find_consistent_signs polys
    in
      find (lookup_sem fml_struct) conds ≠ None
  )"

section "Proofs"
subsection "Parsing and Semantics"
  (* Evaluating a formula where 'a is a real poly *)
primrec real_sem :: "real poly fml ⇒ real ⇒ bool"
  where
    "real_sem (And l r) x = (real_sem l x ∧ real_sem r x)"
  | "real_sem (Or l r) x = (real_sem l x ∨ real_sem r x)"
  | "real_sem (Gt p) x = (poly p x > 0)" 
  | "real_sem (Geq p) x = (poly p x ≥ 0)" 
  | "real_sem (Lt p) x = (poly p x < 0)" 
  | "real_sem (Leq p) x = (poly p x ≤ 0)" 
  | "real_sem (Eq p) x = (poly p x = 0)" 
  | "real_sem (Neq p) x = (poly p x ≠ 0)" 

(* Evaluating a formula where 'a is a rat poly *)
primrec fml_sem :: "rat poly fml ⇒ real ⇒ bool"
  where
    "fml_sem (And l r) x = (fml_sem l x ∧ fml_sem r x)"
  | "fml_sem (Or l r) x = (fml_sem l x ∨ fml_sem r x)"
  | "fml_sem (Gt p) x = (rpoly p x > 0)" 
  | "fml_sem (Geq p) x = (rpoly p x ≥ 0)" 
  | "fml_sem (Lt p) x = (rpoly p x < 0)" 
  | "fml_sem (Leq p) x = (rpoly p x ≤ 0)" 
  | "fml_sem (Eq p) x = (rpoly p x = 0)" 
  | "fml_sem (Neq p) x = (rpoly p x ≠ 0)" 

lemma poly_list_set_fml:
  shows "set (poly_list fml) = set_fml fml"
  apply (induction) by auto

lemma convert_semantics_lem:
  assumes "⋀p. p ∈ set (poly_list fml) ⟹
    ls ! (index_of ps p) = rpoly p x"
  shows "fml_sem fml x = lookup_sem (map_fml (index_of ps) fml) ls"
  using assms apply (induct fml)
  by auto

lemma index_of_aux_more:
  shows "index_of_aux ls p n ≥ n"
  apply (induct ls arbitrary: n)
  apply auto
  using Suc_leD by blast

lemma index_of_aux_lookup:
  assumes "p ∈ set ls"
  shows "(index_of_aux ls p n) - n < length ls"
    "ls ! ((index_of_aux ls p n) - n) = p"
  using assms apply (induct ls arbitrary: n)
  apply auto
  apply (metis Suc_diff_Suc index_of_aux_more lessI less_Suc_eq_0_disj less_le_trans)
  by (metis Suc_diff_Suc index_of_aux_more lessI less_le_trans nth_Cons_Suc)

lemma index_of_lookup:
  assumes "p ∈ set ls"
  shows "index_of ls p < length ls"
    "ls ! (index_of ls p) = p"
  apply (metis assms index_of_aux_lookup(1) index_of_def minus_nat.diff_0)
  by (metis assms index_of_aux_lookup(2) index_of_def minus_nat.diff_0)

lemma convert_semantics:
  shows "fml_sem fml x = lookup_sem (fst (convert fml)) (map (λp. rpoly p x) (snd (convert fml)))"
  unfolding convert_def Let_def apply simp
  apply (intro convert_semantics_lem)
  by (simp add: index_of_lookup(1) index_of_lookup(2))

lemma convert_closed:
  shows "⋀i. i ∈ set_fml (fst (convert fml)) ⟹ i < length (snd (convert fml))"
  unfolding convert_def Let_def
  apply (auto simp add: fml.set_map)
  by (simp add: index_of_lookup(1) poly_list_set_fml)

(* Rational sign vector of polynomials qs with rational coefficients at x *)
definition sign_vec::"rat poly list ⇒ real ⇒ rat list"
  where "sign_vec qs x ≡
    map (squash ∘ (λp. rpoly p x)) qs"

(* The set of all rational sign vectors for qs wrt the set S
  When S = UNIV, then this quantifies over all reals *)
definition consistent_sign_vectors::"rat poly list ⇒ real set ⇒ rat list set"
  where "consistent_sign_vectors qs S = (sign_vec qs) ` S"

lemma sign_vec_semantics:
  assumes "⋀i. i ∈ set_fml fml ⟹ i < length ls"
  shows "lookup_sem fml (map (λp. rpoly p x) ls) = lookup_sem fml (sign_vec ls x)"
  using assms apply (induction)
  by (auto simp add: sign_vec_def squash_def)

(* The universal and existential decision procedure is easy if we know the consistent sign vectors *)
lemma universal_lookup_sem:
  assumes "⋀i. i ∈ set_fml fml ⟹ i < length qs"
  assumes "set signs = consistent_sign_vectors qs UNIV"
  shows "(∀x::real. lookup_sem fml (map (λp. rpoly p x) qs)) ⟷
    list_all (lookup_sem fml) signs"
  using assms(2) unfolding consistent_sign_vectors_def list_all_iff
  by (simp add: assms(1) sign_vec_semantics)

lemma existential_lookup_sem:
  assumes "⋀i. i ∈ set_fml fml ⟹ i < length qs"
  assumes "set signs = consistent_sign_vectors qs UNIV"
  shows "(∃x::real. lookup_sem fml (map (λp. rpoly p x) qs)) ⟷
    find (lookup_sem fml) signs ≠ None"
  using assms(2) unfolding consistent_sign_vectors_def find_None_iff
  by (simp add: assms(1) sign_vec_semantics)

subsection "Factoring Lemmas"
  (*definition real_factorize_list:: "rat poly list ⇒ real poly list"
where "real_factorize_list qs = map (map_poly of_rat) (fst(factorize_polys qs))"
*)
interpretation of_rat_poly_hom: map_poly_comm_semiring_hom of_rat..
interpretation of_rat_poly_hom: map_poly_comm_ring_hom of_rat..
interpretation of_rat_poly_hom: map_poly_idom_hom of_rat..

lemma finite_prod_map_of_rat_poly_hom:
  shows "poly (real_of_rat_poly (∏(a,b)∈s. f a b)) y = (∏(a,b)∈s. poly (real_of_rat_poly (f a b)) y)"
  apply (simp add: of_rat_poly_hom.hom_prod poly_prod)
  by (simp add: case_prod_app prod.case_distrib)

lemma sign_vec_index_of:
  assumes "f ∈ set ftrs"
  shows "sign_vec ftrs x ! (index_of ftrs f) = squash (rpoly f x)"
  by (simp add: assms index_of_lookup(1) index_of_lookup(2) sign_vec_def)

lemma squash_idem:
  shows "squash (squash x) = squash x"
  unfolding squash_def by auto

lemma squash_mult:
  shows "squash ((a::real) * b) = squash a * squash b"
  unfolding squash_def apply auto
  using less_not_sym mult_neg_neg apply blast
  using mult_less_0_iff by blast

lemma squash_prod_list:
  shows "squash (prod_list (ls::real list)) = prod_list (map squash ls)"
  apply (induction ls)
  unfolding squash_def apply auto
  apply (simp add: mult_less_0_iff)
  by (simp add: zero_less_mult_iff)

lemma squash_pow:
  shows "squash ((x::real) ^ (y::nat)) = (squash x) ^ y"
  unfolding squash_def apply auto
  by (auto simp add: zero_less_power_eq)

lemma squash_real_of_rat[simp]:
  shows "squash (real_of_rat x) = squash x"
  unfolding squash_def by auto

lemma factorize_rat_poly_monic_irreducible_monic:
  assumes "factorize_rat_poly_monic f = (c,fs)"
  assumes "(fi,i) ∈ set fs"
  shows "irreducible fi ∧ monic fi"
proof -
  obtain c' fs' where cfs: "factorize_rat_poly f = (c',fs')"
    by (meson surj_pair)
  then have fs: "fs = map (λ(f,i). (normalize f, i)) fs'"
    using factorize_rat_poly_monic_def assms by auto
  obtain "fi'" where "(fi',i) ∈ set fs'" "fi = normalize fi'"
    using assms(2) unfolding fs by auto
  thus ?thesis using factorize_rat_poly irreducible_normalize_iff
    by (metis cfs monic_normalize not_irreducible_zero)
qed

lemma square_free_normalize:
  assumes "square_free p"
  shows "square_free (normalize p)"
  by (metis assms square_free_multD(3) unit_factor_mult_normalize)

lemma coprime_normalize:
  assumes "coprime a b"
  shows "coprime (normalize a) b"
  using assms by auto

lemma undo_normalize:
  shows "a = Polynomial.smult (unit_factor (lead_coeff a)) (normalize a)"
  by (metis add.right_neutral mult_pCons_right mult_zero_right normalize_mult_unit_factor pCons_0_hom.hom_zero unit_factor_poly_def)

lemma finite_smult_distr:
  assumes "distinct fs"
  shows "(∏(x,y)∈set fs. Polynomial.smult ((f x y)::rat) (g x y)) = 
    Polynomial.smult (∏(x,y)∈set fs. f x y) (∏(x,y)∈set fs. g x y)"
  using assms
proof (induction fs)
  case Nil
  then show ?case by auto
next
  case (Cons a fs)
  then show ?case apply auto
    using mult.commute mult_smult_right prod.case_distrib smult_smult split_cong split_conv
    by (simp add: Groups.mult_ac(2) split_beta)
qed

lemma normalize_coprime_degree:
  assumes "normalize (f::rat poly) = normalize g"
  assumes "coprime f g"
  shows "degree f = 0"
proof -
  have "f dvd g" by (simp add: assms(1) associatedD2) 
  then have "f dvd 1"
    using assms(2) associatedD1 by auto
  thus ?thesis
    using Missing_Polynomial_Factorial.is_unit_field_poly by blast
qed

lemma factorize_rat_poly_monic_square_free_factorization:
  assumes res: "factorize_rat_poly_monic f = (c,fs)"
  shows "square_free_factorization f (c,fs)"
proof (unfold square_free_factorization_def split, intro conjI impI allI)
  obtain c' fs' where cfs: "factorize_rat_poly f = (c',fs')"
    by (meson surj_pair)
  then have fs: "fs = map (λ(f,i). (normalize f, i)) fs'"
    using factorize_rat_poly_monic_def assms by auto
  have sq: "square_free_factorization f (c',fs')"
    using cfs factorize_rat_poly(1) by blast
  obtain lcs where lcs: "lcs = prod_list (map (λ(f,i). lead_coeff f ^ Suc i) fs')" by force
  have c: "c = c' * lcs"  using assms unfolding factorize_rat_poly_monic_def cfs Let_def lcs by auto
  show "f = 0 ⟹ c = 0" using c cfs by auto                       
  show "f = 0 ⟹ fs = []" using fs cfs by auto
  have dist: "distinct fs'" using sq square_free_factorizationD(5) by blast
  show dist2: "distinct fs" unfolding fs
    unfolding distinct_conv_nth apply auto
  proof -
    fix i j
    assume ij: "i < length fs'" "j < length fs'" "i ≠ j"
    assume eq: "(case fs' ! i of
            (f, x) ⇒ (normalize f, x)) =
           (case fs' ! j of
            (f, x) ⇒ (normalize f, x))"
    obtain f a where fa: "fs' ! i = (f,a)" by force
    obtain g where g: "fs' ! j = (g,a)" "normalize f = normalize g"
      using eq fa  apply auto
      by (metis case_prod_conv prod.collapse prod.inject)
    have "f ≠ g" using dist ij fa g
      using nth_eq_iff_index_eq by fastforce
    then have "coprime f g"
      using square_free_factorizationD(3)[OF sq, of f a g a] fa g ij
      apply auto
      using nth_mem by force
    then have "degree f = 0"
      by (simp add: g(2) normalize_coprime_degree)
    thus False
      using fa ij(1) nth_mem sq square_free_factorizationD'(3) by fastforce 
  qed
  have ceq: "c = c' * (∏(a, i)∈set fs'. (lead_coeff a) ^ Suc i)" using c lcs
    by (simp add: dist prod.distinct_set_conv_list) 
  have fseq: " (∏(a, i)∈set fs. a ^ Suc i) =  (∏(a, i)∈set fs'. (normalize a) ^ Suc i)"
    apply (subst prod.distinct_set_conv_list[OF dist])
    apply (subst prod.distinct_set_conv_list[OF dist2])
    unfolding fs apply (auto simp add: o_def )
    by (metis (no_types, lifting) case_prod_conv old.prod.exhaust)

  have "f = Polynomial.smult c' (∏(a, i)∈set fs'. a ^ Suc i)"  using sq square_free_factorizationD(1) by blast
  moreover have "... =  Polynomial.smult c' (∏(a, i)∈set fs'. (Polynomial.smult ((unit_factor (lead_coeff a))) (normalize a)) ^ Suc i)"
    apply (subst undo_normalize[symmetric]) by auto
  moreover have "... =  Polynomial.smult c'
    (∏(a, i)∈set fs'. (Polynomial.smult ((lead_coeff a) ^ Suc i) ((normalize a) ^ Suc i)))"
    apply (subst smult_power) by auto
  moreover have "... =  Polynomial.smult c'
    (Polynomial.smult (∏(a, i)∈set fs'. ((lead_coeff a) ^ Suc i))
        (∏(a, i)∈set fs'. (normalize a) ^ Suc i))"
    apply (subst finite_smult_distr) by (auto simp add: dist)
  moreover have "... =  Polynomial.smult (c' * (∏(a, i)∈set fs'. (lead_coeff a) ^ Suc i))
        (∏(a, i)∈set fs'. (normalize a) ^ Suc i)"
    using smult_smult by blast
  moreover have "... = Polynomial.smult c  (∏(a, i)∈set fs. a ^ Suc i)"
    unfolding ceq fseq by auto
  ultimately show "f =  Polynomial.smult c  (∏(a, i)∈set fs. a ^ Suc i)" by auto
  fix a i
  assume ai: "(a,i) ∈ set fs"
  obtain a' where a': "(a',i) ∈ set fs'" "a = normalize a'" using ai unfolding fs by auto
  show "square_free a" using square_free_normalize a'
    using sq square_free_factorizationD(2) by blast
  show "0 < degree a" using degree_normalize a'
    using sq square_free_factorizationD'(3) by fastforce
  fix  b j
  assume bj: "(b,j) ∈ set fs" "(a,i) ≠ (b,j)"
  obtain b' where b': "(b',j) ∈ set fs'" "b = normalize b'" using bj unfolding fs by auto
  show "algebraic_semidom_class.coprime a b" using a' b' apply auto
    using bj(2) sq square_free_factorizationD(3) by fastforce
qed  

lemma undo_factorize_correct:
  assumes "factorize_rat_poly_monic p = (c,fs)"
  assumes "⋀f p. (f,p) ∈ set fs ⟹ f ∈ set ftrs"
  shows "undo_factorize (c,map (λ(f,pow). (index_of ftrs f, pow)) fs) (sign_vec ftrs x) = squash (rpoly p x)"
proof -
  have p: "p = smult c (∏(a, i)∈ set fs. a ^ Suc i)"
    using assms(1) factorize_rat_poly_monic_square_free_factorization square_free_factorizationD(1) by blast
  have fs: "distinct fs"
    using assms(1) factorize_rat_poly_monic_square_free_factorization square_free_factorizationD(5) by blast
  have "rpoly p x = ((real_of_rat c) * rpoly (∏(a, i)∈ set fs. a ^ Suc i) x)"
    using p by (simp add: of_rat_hom.map_poly_hom_smult)
  moreover have "... = ((real_of_rat c) * rpoly (∏ai∈ set fs. case ai of (a,i) ⇒ a ^ Suc i) x)"
    by blast
  moreover have "... = ((real_of_rat c) * (∏ai∈ set fs. case ai of (a,i) ⇒ rpoly (a ^ Suc i) x))"
    by (simp add: finite_prod_map_of_rat_poly_hom)
  moreover have "... = ((real_of_rat c) * (∏ai∈ set fs. case ai of (a,i) ⇒ (rpoly a x) ^ Suc i))"
    by (metis (mono_tags, lifting) of_rat_poly_hom.hom_power poly_hom.hom_power split_cong)
  moreover have "...  = ((real_of_rat c) * (prod_list (map (λai. case ai of (a,i) ⇒ (rpoly a x) ^ Suc i) fs)))"
    by (simp add: fs prod.distinct_set_conv_list)
  ultimately have "rpoly p x = ((real_of_rat c) * (prod_list (map (λai. case ai of (a,i) ⇒ (rpoly a x) ^ Suc i) fs)))" by auto

  then have "squash (rpoly p x) = squash c * prod_list (map squash (map (λai. case ai of (a,i) ⇒ (rpoly a x) ^ Suc i) fs))"
    by (auto simp add: squash_mult squash_prod_list o_def)
  moreover have "... = squash c * prod_list (map (λai. case ai of (a,i) ⇒ squash ((rpoly a x) ^ Suc i)) fs)"
    apply (simp add: o_def)
    by (simp add: prod.case_distrib)
  ultimately have rp:"squash(rpoly p x) = squash c * prod_list (map (λai. case ai of (a,i) ⇒ squash (rpoly a x) ^ Suc i) fs)"
    using squash_pow
    by presburger
  have "undo_factorize
     (c, map (λ(f, pow).(index_of ftrs f, pow)) fs) (sign_vec ftrs x) =
    squash
     (c * (∏xa←fs. case xa of (f, y) ⇒ sign_vec ftrs x ! index_of ftrs f ^ Suc y))"
    unfolding undo_factorize_def apply (auto simp add: o_def)
    by (metis (mono_tags, lifting) case_prod_conv old.prod.exhaust)
  moreover have "... =  squash
     (c * (∏xa←fs. case xa of (f, y) ⇒ (squash (rpoly f x)) ^ Suc y))"
    using assms(2) sign_vec_index_of
      map_eq_conv split_cong
    apply (auto) 
    by (smt map_eq_conv split_cong)
  ultimately show ?thesis using rp
    by (metis (mono_tags, lifting) of_rat_hom.hom_mult squash_idem squash_mult squash_real_of_rat)
qed

lemma length_sign_vec[simp]:
  shows "length (sign_vec ps x) = length ps" unfolding sign_vec_def by auto

lemma factorize_polys_has_factors:
  assumes "factorize_polys ps = (ftrs,data)"
  assumes "p ∈ set ps"
  assumes "factorize_rat_poly_monic p = (c,fs)"
  shows "set (map fst fs) ⊆ set ftrs"
  using assms unfolding factorize_polys_def Let_def apply auto
  by (metis UN_iff fst_conv image_eqI snd_conv)

lemma factorize_polys_undo_factorize_polys:
  assumes "factorize_polys ps = (ftrs,data)"
  shows "undo_factorize_polys data (sign_vec ftrs x) = sign_vec ps x"
  unfolding list_eq_iff_nth_eq undo_factorize_polys_def apply auto
proof -
  show leq:"length data = length ps"
    using assms unfolding factorize_polys_def by (auto simp add: Let_def)
  fix i
  assume il:"i < length data"
  obtain c fs where cfs: "factorize_rat_poly_monic (ps ! i) = (c,fs)"
    by (meson surj_pair)
  then have fsts:"set (map fst fs) ⊆ set ftrs"
    using assms factorize_polys_has_factors il leq nth_mem by fastforce
  have *:"data ! i = (c,map (λ(f,pow). (index_of ftrs f, pow)) fs)"
    using assms unfolding factorize_polys_def
    using cfs il by (auto simp add: Let_def cfs)
  have "undo_factorize (data ! i) (sign_vec ftrs x) = squash (rpoly (ps ! i) x)" unfolding *
    apply (subst  undo_factorize_correct[of "ps ! i"])
    apply (auto simp add: cfs)
    using fsts by auto
  thus "undo_factorize (data ! i) (sign_vec ftrs x) =  sign_vec ps x ! i" 
    using leq il sign_vec_def by auto
qed

lemma factorize_polys_irreducible_monic:
  assumes "factorize_polys ps = (fs,data)"
  shows "distinct fs" "⋀f. f ∈ set fs ⟹ irreducible f ∧ monic f"
  using assms unfolding factorize_polys_def Let_def apply auto
  using factorize_rat_poly_monic_irreducible_monic
  apply (metis prod.collapse) 
  using factorize_rat_poly_monic_irreducible_monic
  by (metis prod.collapse)

lemma factorize_polys_square_free:
  assumes "factorize_polys ps = (fs,data)"
  shows "⋀f. f ∈ set fs ⟹ square_free f"
  using assms factorize_polys_irreducible_monic(2) irreducible_imp_square_free by blast

lemma irreducible_monic_coprime:
  assumes f: "monic f" "irreducible (f::rat poly)" 
  assumes g: "monic g" "irreducible (g::rat poly)"
  assumes "f ≠ g"
  shows "coprime f g"
  by (metis (no_types, lifting) assms(5) coprime_0(2) coprime_def' f(1) f(2) g(1) g(2) irreducible_normalized_divisors normalize_dvd_iff normalize_idem normalize_monic)

lemma factorize_polys_coprime:
  assumes "factorize_polys ps = (fs,data)"
  shows "⋀f g. f ∈ set fs ⟹ g ∈ set fs ⟹ f ≠ g ⟹ coprime f g"
  using assms factorize_polys_irreducible_monic(2) irreducible_monic_coprime by auto

lemma coprime_rat_poly_real_poly:
  assumes "coprime p (q::rat poly)"
  shows "coprime (real_of_rat_poly p) ((real_of_rat_poly q)::real poly)"
  by (metis assms gcd_dvd_1 of_rat_hom.map_poly_gcd of_rat_poly_hom.hom_dvd_1)

lemma coprime_rat_poly_iff_coprimereal_poly:
  shows "coprime p (q::rat poly) ⟷ coprime (real_of_rat_poly p) ((real_of_rat_poly q)::real poly)"
proof -
  have forward: "coprime p (q::rat poly) ⟶ coprime (real_of_rat_poly p) ((real_of_rat_poly q)::real poly)"
    using coprime_rat_poly_real_poly by auto
  have backward: "coprime (real_of_rat_poly p) ((real_of_rat_poly q)::real poly) ⟹ coprime p (q::rat poly)"
  proof -
    assume copr_real: "comm_monoid_mult_class.coprime (real_of_rat_poly p) (real_of_rat_poly q)"
    have "degree (gcd p (q::rat poly)) > 0 ⟹ False" 
    proof -
      assume deg: "degree (gcd p (q::rat poly)) > 0"
      then have "∃y. y dvd p ∧ y dvd q ∧ degree y > 0"
        by blast
      then obtain y where yprop: "y dvd p ∧ y dvd q ∧ degree y > 0"
        by auto
      then have "(real_of_rat_poly y) dvd (real_of_rat_poly p) ∧
        (real_of_rat_poly y ) dvd (real_of_rat_poly q) ∧ degree y > 0"
        by simp
      then show "False" 
        using copr_real apply (auto)
        by fastforce 
    qed
    then show "comm_monoid_mult_class.coprime p (q::rat poly)" 
      using comm_monoid_gcd_class.gcd_dvd_1
      by (metis Missing_Polynomial_Factorial.is_unit_field_poly copr_real gcd_zero_iff' neq0_conv of_rat_poly_hom.hom_zero)
  qed
  show ?thesis
    using forward backward by auto
qed

lemma factorize_polys_map_distinct:
  assumes "factorize_polys ps = (fs,data)"
  assumes "fss =  map real_of_rat_poly fs"
  shows "distinct fss"
  using factorize_polys_irreducible_monic[OF assms(1)]
  unfolding assms(2) 
  apply (simp add: distinct_conv_nth)
  by (metis of_rat_eq_iff of_rat_hom.coeff_map_poly_hom poly_eqI)

lemma factorize_polys_map_square_free:
  assumes "factorize_polys ps = (fs,data)"
  assumes "fss =  map real_of_rat_poly fs"
  shows "⋀f. f ∈ set fss ⟹ square_free f"
  using factorize_polys_square_free[OF assms(1)]
  using assms(2) field_hom_0'.square_free_map_poly of_rat_hom.field_hom_0'_axioms by auto 

lemma factorize_polys_map_coprime:
  assumes "factorize_polys ps = (fs,data)"
  assumes "fss =  map real_of_rat_poly fs"
  shows "⋀f g. f ∈ set fss ⟹ g ∈ set fss ⟹ f ≠ g ⟹ coprime f g"
  using factorize_polys_coprime[OF assms(1)] coprime_rat_poly_real_poly unfolding assms(2)
  by auto

lemma coprime_prod_list:
  assumes "⋀p. p ∈ set ps ⟹ p ≠ 0"
  assumes "coprime (prod_list ps) (q::real poly)"  
  shows "⋀p. p ∈ set ps ⟹ coprime p q"
proof -
  fix p
  assume "p ∈ set ps"
  then obtain r where r: "prod_list ps = r * p"
    using remove1_retains_prod by blast
  show "coprime p q"
    apply (rule coprime_prod[of r 1])
    using assms r apply auto
    by blast
qed

(* basically copied from square_free_factorizationD' *)
lemma factorize_polys_square_free_prod_list:
  assumes "factorize_polys ps = (fs,data)"
  shows "square_free (prod_list fs)"
proof (rule square_freeI)
  from factorize_polys_coprime[OF assms]
  have coprime: "⋀p q. p ∈ set fs ⟹ q ∈ set fs ⟹ p ≠ q ⟹ coprime p q" .
  from factorize_polys_square_free[OF assms]
  have sq: "⋀p. p ∈ set fs ⟹ square_free p" .
  thus "prod_list fs ≠ 0" unfolding prod_list_zero_iff
    using square_free_def by blast
  fix q
  assume "degree q > 0" "q * q dvd prod_list fs"
  from irreducible⇩d_factor[OF this(1)] this(2) obtain q where 
    irr: "irreducible q" and dvd: "q * q dvd prod_list fs" unfolding dvd_def by auto
  hence dvd': "q dvd prod_list fs" unfolding dvd_def by auto
  from irreducible_dvd_prod_list[OF irr dvd'] obtain b where 
    mem: "b ∈ set fs" and dvd1: "q dvd b" by auto
  from dvd1 obtain k where b: "b = q * k" unfolding dvd_def by auto
  from split_list[OF mem] b obtain bs1 bs2 where bs: "fs = bs1 @ b # bs2" by auto
  from irr have q0: "q ≠ 0" and dq: "degree q > 0" unfolding irreducible⇩d_def by auto
  have "square_free (q * k)" using sq b mem by auto
  from this[unfolded square_free_def, THEN conjunct2, rule_format, OF dq] 
  have qk: "¬ q dvd k" by simp
  from dvd[unfolded bs b] have "q * q dvd q * (k * prod_list (bs1 @ bs2))"
    by (auto simp: ac_simps)
  with q0 have "q dvd k * prod_list (bs1 @ bs2)" by auto
  with irr qk have "q dvd prod_list (bs1 @ bs2)" by auto
  from irreducible_dvd_prod_list[OF irr this] obtain b' where 
    mem': "b' ∈ set (bs1 @ bs2)" and dvd2: "q dvd b'" by fastforce
  from dvd1 dvd2 have "q dvd gcd b b'" by auto
  with dq is_unit_iff_degree[OF q0] have cop: "¬ coprime b b'" by force
  from mem' have "b' ∈ set fs" unfolding bs by auto
  have b': "b' = b" using coprime
    using ‹b' ∈ set fs› cop mem by blast
  with mem' bs factorize_polys_irreducible_monic(1)[OF assms] show False by auto
qed

lemma factorize_polys_map_square_free_prod_list:
  assumes "factorize_polys ps = (fs,data)"
  assumes "fss =  map real_of_rat_poly fs"
  shows "square_free (prod_list fss)"
  using  factorize_polys_square_free_prod_list[OF assms(1)] unfolding assms(2)
  by (simp add: of_rat_hom.square_free_map_poly)

lemma factorize_polys_map_coprime_pderiv:
  assumes "factorize_polys ps = (fs,data)"
  assumes "fss = map real_of_rat_poly fs"
  shows "⋀f. f ∈ set fss ⟹ coprime f (pderiv (prod_list fss))"
proof -
  fix f
  assume f: "f ∈ set fss"
  from factorize_polys_map_square_free[OF assms]
  have sq: "⋀p. p ∈ set fss ⟹ square_free p" .
  have z: "⋀p. p ∈ set fss ⟹ p ≠ 0" using sq square_free_def by blast
  have c: "coprime (prod_list fss) (pderiv (prod_list fss))"
    apply (simp add: separable_def[symmetric] square_free_iff_separable[symmetric])
    using factorize_polys_map_square_free_prod_list[OF assms] .
  from coprime_prod_list[OF z c f]
  show "coprime f (pderiv (prod_list fss))" by auto
qed

definition pairwise_coprime_list:: "rat poly list ⇒ bool"
  where "pairwise_coprime_list qs = 
    (∀m < length qs. ∀ n < length qs.
     m ≠ n ⟶ coprime (qs ! n) (qs ! m))"

(* Restating factorize_polys_map_coprime to match later definitions *)
lemma coprime_factorize:
  fixes qs:: "rat poly list"
  shows "pairwise_coprime_list (fst(factorize_polys qs))"
proof -
  let ?fs = "fst(factorize_polys qs)"
  have "(∀m < length ?fs. ∀ n < length ?fs.
     m ≠ n ⟶ coprime (?fs ! n) (?fs ! m))"
  proof clarsimp 
    fix m n
    assume "m < length (fst (factorize_polys qs))"
    assume "n < length (fst (factorize_polys qs))"
    assume "m ≠ n"
    show " algebraic_semidom_class.coprime (fst (factorize_polys qs) ! n)
            (fst (factorize_polys qs) ! m)"
      by (metis ‹m < length (fst (factorize_polys qs))› ‹m ≠ n› ‹n < length (fst (factorize_polys qs))› coprime_iff_coprime distinct_conv_nth factorize_polys_coprime factorize_polys_def factorize_polys_irreducible_monic(1) fstI nth_mem)
  qed
  then show ?thesis unfolding pairwise_coprime_list_def by auto
qed

lemma squarefree_factorization_degree:
  assumes "square_free_factorization p (c,fs)"
  shows "degree p = sum_list (map (λ(f,c). (c+1) * degree f) fs)"
proof -
  have "p =
    Polynomial.smult c
     (∏(a, i)∈set fs. a ^ Suc i)" using assms unfolding square_free_factorization_def
    by blast
  then have "degree p = degree (∏(a, i)∈set fs. a ^ Suc i)"
    using assms square_free_factorizationD(4) by fastforce
  also have "... = degree (prod_list (map (λ(f,c). f ^ Suc c) fs))"
    by (metis assms prod.distinct_set_conv_list square_free_factorizationD(5))
  also have "... = (∑(a, i)←fs. degree (a ^ Suc i))"
    apply (subst degree_prod_list_eq)
    apply (auto simp add: o_def)
    using assms degree_0 square_free_factorizationD(2) apply blast
    using assms degree_0 square_free_factorizationD(2) apply blast
    by (simp add: prod.case_distrib)
  ultimately show ?thesis
    by (smt Polynomial.degree_power_eq add.commute assms degree_0 map_eq_conv plus_1_eq_Suc split_cong square_free_factorizationD(2))  
qed

subsection "Auxiliary Polynomial Lemmas"
definition roots_of_coprime_r:: "real poly list ⇒ real set"
  where "roots_of_coprime_r qs = {x. poly (coprime_r qs) x = 0}"

lemma crb_lem_pos: 
  fixes x:: "real"
  fixes p:: "real poly"
  assumes x: "poly p x = 0" 
  assumes p: "p ≠ 0" 
  shows "x < crb p"
  using cauchy_root_bound[of p x] apply (auto)
  unfolding crb_def apply (auto)
  using p x 
  by linarith 

lemma crb_lem_neg: 
  fixes x:: "real"
  fixes p:: "real poly"
  assumes x: "poly p x = 0" 
  assumes p: "p ≠ 0" 
  shows "x > -crb p"
  using cauchy_root_bound[of p x] apply (auto)
  unfolding crb_def apply (auto)
  using p x by linarith 

(* Show that the product of the polynomial list is 0 at x iff there is a polynomial 
  in the list that is 0 at x *)
lemma prod_zero:
  shows "∀x . poly (prod_list (qs:: rat poly list)) x = 0 ⟷ (∃q ∈ set (qs). poly q x = 0)" 
  apply auto
  using poly_prod_list_zero_iff apply blast
  using poly_prod_list_zero_iff by blast

lemma coprime_r_zero1: "poly (coprime_r qs) (crb (prod_list qs)) = 0"
  by (simp add: coprime_r_def) 

lemma coprime_r_zero2: "poly (coprime_r qs) (-crb (prod_list qs)) = 0"
  by (simp add: coprime_r_def) 

lemma coprime_mult:
  fixes a:: "real poly"
  fixes b:: "real poly"
  fixes c:: "real poly"
  assumes "algebraic_semidom_class.coprime a b"
  assumes "algebraic_semidom_class.coprime a c"
  shows "algebraic_semidom_class.coprime a (b*c)"
  using assms(1) assms(2) by auto

(* Will be needed when we call the BKR roots on coprime_r *)
lemma coprime_r_coprime_prop:
  fixes ps:: "rat poly list"
  assumes "factorize_polys ps = (fs,data)"
  assumes "fss = map real_of_rat_poly fs"
  shows "⋀f. f ∈ set fss ⟹ coprime f (coprime_r fss)"
proof clarsimp 
  fix f:: "real poly"
  assume f_in: "f ∈ set fss"
  have nonz_prod: "prod_list fss ≠ 0" using factorize_polys_map_square_free apply (auto)
    using assms(1) assms(2) square_free_def by fastforce 
  have nonz_f: "f ≠ 0" using f_in factorize_polys_map_square_free apply (auto)
    using assms(1) assms(2) square_free_def by fastforce 
  have copr_pderiv: "algebraic_semidom_class.coprime f (pderiv (prod_list fss))" using factorize_polys_map_coprime_pderiv
    apply (auto)
    using f_in assms(1) assms(2) by auto 
  have z_iff: "∀x. poly f x = 0 ⟶ poly (prod_list fss) x = 0"
    using f_in apply (auto)
    using poly_prod_list_zero_iff by blast 
  let ?inf_p = "[:-(crb (prod_list fss)),1:]::real poly"
  have copr_inf: "algebraic_semidom_class.coprime f ([:-(crb (prod_list fss)),1:])"
  proof - 
    have zero_prop: "∀x. poly ?inf_p x = 0 ⟷ x = crb (prod_list fss)" 
      by auto
    have  "poly (prod_list fss) (crb (prod_list fss)) ≠ 0"
    proof -
      have h: "∀x. poly (prod_list fss) x = 0 ⟶ x < (crb (prod_list fss))"
        using nonz_prod crb_lem_pos[where p = "prod_list fss"] 
        by auto
      then  show ?thesis by auto
    qed
    then have nonzero: "poly f (crb (prod_list fss)) ≠ 0"
      using z_iff by auto
    then have "¬(∃x. poly f x = 0 ∧ poly ?inf_p x = 0)"
      by simp
    have is_unit_gcd: "is_unit (gcd ?inf_p f)"
      using prime_elem_imp_gcd_eq  prime_elem_iff_irreducible linear_irreducible_field
      apply (auto) using nonzero
    proof -
      have f1: "∀x0. - (x0::real) = - 1 * x0"
        by simp
      have "(1::real) ≠ 0"
        by auto
      then have "is_unit (gcd (pCons (- 1 * real_of_int (crb (prod_list fss))) 1) f)"
        using f1 by (metis (no_types) is_unit_gcd nonzero one_poly_eq_simps(1) poly_eq_0_iff_dvd prime_elem_imp_coprime prime_elem_linear_field_poly)
      then show "degree (gcd (pCons (- real_of_int (crb (prod_list fss))) 1) f) = 0"
        by simp
    qed 
    then show ?thesis
      using is_unit_gcd
      by (metis gcd.commute gcd_eq_1_imp_coprime is_unit_gcd_iff) 
  qed
  let ?ninf_p = "[:(crb (prod_list fss)),1:]::real poly"
  have copr_neg_inf: "algebraic_semidom_class.coprime f ([:(crb (prod_list fss)),1:])"
  proof - 
    have h: "∀x. poly f x = 0 ⟶ poly (prod_list fss) x = 0"
      using f_in apply (auto)
      using poly_prod_list_zero_iff by blast 
    have zero_prop: "∀x. poly ?ninf_p x = 0 ⟷ x = -crb (prod_list fss)" 
      by auto
    have  "poly (prod_list fss) (-crb (prod_list fss)) ≠ 0"
    proof  -
      have h: "∀x. poly (prod_list fss) x = 0 ⟶ x > (-crb (prod_list fss))"
        using nonz_prod crb_lem_neg[where p = "prod_list fss"] 
        by auto
      then  show ?thesis by auto
    qed
    then have nonzero: "poly f (-crb (prod_list fss)) ≠ 0"
      using z_iff by auto
    then have "¬(∃x. poly f x = 0 ∧ poly ?ninf_p x = 0)"
      using zero_prop by auto
    have is_unit_gcd: "is_unit (gcd ?ninf_p f)"
      using prime_elem_imp_gcd_eq  prime_elem_iff_irreducible linear_irreducible_field
      apply (auto) using nonzero
    proof -
      have f1: "(1::real) ≠ 0"
        by auto
      have "¬ pCons (real_of_int (crb (prod_list fss))) 1 dvd f"
        using nonzero by auto
      then show "degree (gcd (pCons (real_of_int (crb (prod_list fss))) 1) f) = 0"
        using f1 by (metis (no_types) Missing_Polynomial_Factorial.is_unit_field_poly coprime_imp_gcd_eq_1 is_unit_gcd_iff one_poly_eq_simps(1) prime_elem_imp_coprime prime_elem_linear_field_poly)
    qed
    then show ?thesis 
      using is_unit_gcd
      by (metis gcd.commute gcd_eq_1_imp_coprime is_unit_gcd_iff) 
  qed
  show "algebraic_semidom_class.coprime f (coprime_r fss)" 
    using copr_pderiv coprime_mult unfolding coprime_r_def
    using copr_inf copr_neg_inf by blast 
qed

lemma coprime_r_nonzero:
  fixes ps:: "rat poly list"
  assumes "factorize_polys ps = (fs,data)"
  assumes nonempty_fs: "fs ≠ []"
  assumes fss_is: "fss = map real_of_rat_poly fs"
  shows "(coprime_r fss) ≠ 0"
proof - 
  have nonempty_fss: "fss ≠ []" using nonempty_fs fss_is by auto
  have deg_f: "∀f ∈ set (fs). degree f > 0"
    using factorize_polys_irreducible_monic 
    apply (auto)
    using assms(1) irreducible_degree_field by blast
  then have deg_fss: "∀f ∈ set (fss). degree f > 0"
    using fss_is by simp 
  then have fss_nonz: "∀f ∈ set (fss). f ≠ 0"
    by auto
  have "fss ≠ [] ⟶ ((∀f ∈ set (fss). (degree f > 0 ∧ f ≠ 0)) ⟶ degree (prod_list fss) > 0)"
  proof (induct fss)
    case Nil
    then show ?case
      by blast 
  next
    case (Cons a fss)
    show ?case 
    proof clarsimp 
      assume z_lt: "0 < degree a"
      assume anonz: "a ≠ 0"
      assume fnonz: "∀f∈set fss. 0 < degree f ∧ f ≠ 0"
      have h: "degree (a * prod_list fss) = degree a + degree (prod_list fss) "
        using degree_mult_eq[where p = "a", where q = "prod_list fss"] anonz fnonz 
        by auto
      then show "0 < degree (a * prod_list fss)"
        using z_lt Cons.hyps by auto
    qed
  qed
  then have "degree (prod_list fss) > 0"
    using nonempty_fss deg_fss fss_nonz by auto
  then have pderiv_nonzero: "pderiv (prod_list fss) ≠ 0"
    by (simp add: pderiv_eq_0_iff)
  have "(([:-(crb (prod_list fss)),1:]) * ([:(crb (prod_list fss)),1:])) ≠ 0"
    by auto
  then show ?thesis using pderiv_nonzero
    unfolding coprime_r_def apply (auto)
    by (metis offset_poly_eq_0_lemma right_minus_eq synthetic_div_unique_lemma) 
qed

lemma Rolle_pderiv:
  fixes q:: "real poly"
  fixes x1 x2:: "real"
  shows "(x1 < x2 ∧ poly q x1 = 0 ∧ poly q x2 = 0) ⟶ (∃w. x1 < w ∧ w < x2 ∧ poly (pderiv q) w = 0)"
  using Rolle_deriv apply (auto)
  by (metis DERIV_unique Rolle continuous_at_imp_continuous_on poly_DERIV poly_differentiable poly_isCont) 

lemma coprime_r_roots_prop: 
  fixes qs:: "real poly list"
  assumes pairwise_rel_prime: "∀q1 q2. (q1 ≠ q2 ∧ (List.member qs q1) ∧ (List.member qs q2))⟶ coprime q1 q2"
  shows "∀x1. ∀x2. ((x1 < x2 ∧ (∃q1 ∈ set (qs). (poly q1 x1) = 0) ∧ (∃q2∈ set(qs). (poly q2 x2) = 0)) ⟶ (∃q. x1 < q ∧ q < x2 ∧ poly (coprime_r qs) q = 0))"
proof clarsimp
  fix x1:: "real"
  fix x2:: "real"
  fix q1:: "real poly"
  fix q2:: "real poly"
  assume "x1 < x2"
  assume q1_in: "q1 ∈ set qs"
  assume q1_0: "poly q1 x1 = 0"
  assume q2_in: "q2 ∈ set qs"
  assume q2_0: "poly q2 x2 = 0"
  have prod_z_x1: "poly (prod_list qs) x1 = 0" using q1_in q1_0
    using poly_prod_list_zero_iff by blast  
  have prod_z_x2: "poly (prod_list qs) x2 = 0" using q2_in q2_0 
    using poly_prod_list_zero_iff by blast 
  have "∃w>x1. w < x2 ∧ poly (pderiv (prod_list qs)) w = 0" 
    using Rolle_pderiv[where q = "prod_list qs"] prod_z_x1 prod_z_x2
    using ‹x1 < x2› by blast
  then obtain w where w_def: "w > x1 ∧w < x2 ∧ poly (pderiv (prod_list qs)) w = 0"
    by auto
  then have "poly (coprime_r qs) w = 0"
    unfolding coprime_r_def
    by simp 
  then show "∃q>x1. q < x2 ∧ poly (coprime_r qs) q = 0"
    using w_def by blast 
qed

subsection "Setting Up the Procedure: Lemmas"
definition has_no_zeros::"rat list ⇒ bool"
  where "has_no_zeros l = (0 ∉ set l)"

lemma hnz_prop: "has_no_zeros l ⟷ ¬(∃k < length l. l ! k = 0)"
  unfolding has_no_zeros_def
  by (simp add: in_set_conv_nth)

definition cast_rat_list:: "rat poly list ⇒ real poly list"
  where "cast_rat_list qs = map real_of_rat_poly qs"

definition consistent_sign_vectors_r::"real poly list ⇒ real set ⇒ rat list set"
  where "consistent_sign_vectors_r qs S = (signs_at qs) ` S"

lemma consistent_sign_vectors_consistent_sign_vectors_r:
  shows"consistent_sign_vectors_r (cast_rat_list qs) S = consistent_sign_vectors qs S"
  unfolding consistent_sign_vectors_r_def cast_rat_list_def consistent_sign_vectors_def
    sign_vec_def signs_at_def
  by auto

(* Relies on coprime_rat_poly_real_poly *)
lemma coprime_over_reals_coprime_over_rats:
  fixes qs:: "rat poly list"
  assumes csa_in: "csa ∈ (consistent_sign_vectors qs UNIV)"
  assumes p1p2: "p1≠p2 ∧ p1 < length csa ∧ p2 < length csa ∧ csa ! p1 = 0 ∧ csa ! p2 = 0"
  shows "¬ algebraic_semidom_class.coprime (nth qs p1) (nth qs p2) "
proof - 
  have isx: "∃(x::real). csa = (sign_vec qs x)" 
    using csa_in unfolding consistent_sign_vectors_def by auto
  then obtain x where havex: "csa = (sign_vec qs x)" by auto
  then have expolys: "poly (real_of_rat_poly (nth qs p1)) x = 0 ∧ poly (real_of_rat_poly (nth qs p2)) x = 0"
    using havex unfolding sign_vec_def squash_def
    by (smt class_field.neg_1_not_0 length_map map_map nth_map one_neq_zero p1p2)
  then have "¬ coprime (real_of_rat_poly (nth qs p1)) ((real_of_rat_poly (nth qs p2))::real poly)"
    apply (auto) using coprime_poly_0 
    by blast 
  then show ?thesis
    using coprime_rat_poly_real_poly by auto
qed

(* This and the following lemma are designed to show that if you have two sgas that aren't the same, 
   there's a 0 in between! The proof uses IVT. It hones in on the component that's changed sign 
   (either from 1 to {0, -1} or from -1 to {0, 1}) *)
lemma zero_above: 
  fixes qs:: "rat poly list"
  fixes x1:: "real"
  assumes hnz: "has_no_zeros (sign_vec qs x1)"
  shows "(∀ x2 > x1. ((sign_vec qs x1) ≠ (sign_vec qs x2)) ⟶ 
  (∃(r::real) > x1. (r ≤ x2 ∧ (∃ q ∈ set(qs). rpoly q r = 0))))"
proof clarsimp 
  fix x2
  assume x1_lt: "x1 < x2"
  assume diff_sign_vec: "sign_vec qs x1 ≠ sign_vec qs x2"
  then have "∃q ∈ set qs. squash (rpoly q x1) ≠ squash (rpoly q x2)"
    unfolding sign_vec_def
    by simp 
  then obtain q where q_prop: "q ∈ set qs ∧ squash (rpoly q x1) ≠ squash (rpoly q x2)"
    by auto
  then have q_in: "q ∈ set qs" by auto
  have poss1: "squash (rpoly q x1) = -1 ∧ squash (rpoly q x2) = 1 ⟶ (∃r>x1. r ≤ x2 ∧ (∃q∈set qs. rpoly q r = 0))"
    using poly_IVT_pos[of x1 x2] using x1_lt unfolding squash_def apply (auto)
    using q_prop by fastforce 
  have poss2: "squash (rpoly q x1) = 1 ∧ squash (rpoly q x2) = -1 ⟶ (∃r>x1. r ≤ x2 ∧ (∃q∈set qs. rpoly q r = 0))"
    using poly_IVT_neg[of x1 x2] using x1_lt unfolding squash_def apply (auto)
    using q_prop by fastforce 
  have poss3: "squash (rpoly q x2) = 0 ⟶ (∃r>x1. r ≤ x2 ∧ (∃q∈set qs. rpoly q r = 0))"
    using x1_lt unfolding squash_def apply (auto)
    using q_prop by blast
  have "(q ∈ set qs ∧ rpoly q x1 = 0) ⟶ ¬has_no_zeros (sign_vec qs x1)"
    unfolding has_no_zeros_def sign_vec_def apply auto
    by (smt image_iff squash_def)
  have not_poss4: "squash (rpoly q x1) ≠ 0" 
    using hnz q_in unfolding squash_def
    using ‹q ∈ set qs ∧ rpoly q x1 = 0 ⟶ ¬ has_no_zeros (sign_vec qs x1)› by auto 
  then show "∃r>x1. r ≤ x2 ∧ (∃q∈set qs. rpoly q r = 0)"
    using q_prop poss1 poss2 poss3 not_poss4 
    apply (auto)
     apply (meson squash_def)
      apply (metis squash_def)
       apply (metis squash_def) by (meson squash_def)
qed

lemma zero_below: 
  fixes qs:: "rat poly list"
  fixes x1:: "real"
  assumes hnz: "has_no_zeros (sign_vec qs x1)"
  shows "∀x2 < x1. ((sign_vec qs x1) ≠ (sign_vec qs x2)) ⟶ 
  (∃(r::real) < x1. (r ≥ x2 ∧ (∃ q ∈ set(qs). rpoly q r = 0)))"
proof clarsimp 
  fix x2
  assume x1_gt: "x2 < x1"
  assume diff_sign_vec: "sign_vec qs x1 ≠ sign_vec qs x2"
  then have "∃q ∈ set qs. squash (rpoly q x1) ≠ squash (rpoly q x2)"
    unfolding sign_vec_def
    by simp 
  then obtain q where q_prop: "q ∈ set qs ∧ squash (rpoly q x1) ≠ squash (rpoly q x2)"
    by auto
  then have q_in: "q ∈ set qs" by auto
  have poss1: "squash (rpoly q x1) = -1 ∧ squash (rpoly q x2) = 1 ⟶ (∃r<x1. (r ≥ x2 ∧ (∃ q ∈ set(qs). rpoly q r = 0)))"
    using poly_IVT_neg[of x2 x1] using x1_gt unfolding squash_def apply (auto)
    using q_prop by fastforce 
  have poss2: "squash (rpoly q x1) = 1 ∧ squash (rpoly q x2) = -1 ⟶ (∃r<x1. (r ≥ x2 ∧ (∃ q ∈ set(qs). rpoly q r = 0)))"
    using poly_IVT_pos[of x2 x1] using x1_gt unfolding squash_def apply (auto)
    using q_prop by fastforce 
  have poss3: "squash (rpoly q x2) = 0 ⟶ (∃r<x1. (r ≥ x2 ∧ (∃ q ∈ set(qs). rpoly q r = 0)))"
    using x1_gt unfolding squash_def apply (auto)
    using q_prop by blast
  have "(q ∈ set qs ∧ rpoly q x1 = 0) ⟶ ¬has_no_zeros (sign_vec qs x1)"
    unfolding has_no_zeros_def sign_vec_def apply auto
    using image_iff squash_def
    by (smt image_iff squash_def)
  have not_poss4: "squash (rpoly q x1) ≠ 0" 
    using hnz q_in unfolding squash_def
    using ‹q ∈ set qs ∧ rpoly q x1 = 0 ⟶ ¬ has_no_zeros (sign_vec qs x1)› by auto 
  then show "(∃r<x1. (r ≥ x2 ∧ (∃ q ∈ set(qs). rpoly q r = 0)))"
    using q_prop poss1 poss2 poss3 not_poss4 apply (auto)
    apply (meson squash_def)
     apply (metis squash_def)
      apply (metis squash_def)
       by (meson squash_def)
qed

lemma sorted_list_lemma:
  fixes l:: "real list"
  fixes a b:: "real"
  fixes n:: "nat"
  assumes "a < b"
  assumes "(n + 1) < length l"
  assumes strict_sort: "sorted_wrt (<) l"
  assumes lt_a: "(l ! n) < a"
  assumes b_lt: "b < (l ! (n+1))"
  shows "¬(∃(x::real). (List.member l x ∧ a ≤ x ∧ x ≤ b))"
proof -
  have sorted_hyp_var: "∀q1 < length l. ∀q2 < length l. (q1 < q2 ⟶
     (l ! q1) < (l ! q2))" 
    using strict_sort by (auto simp: sorted_wrt_iff_nth_less)
  then have sorted_hyp_var2:  "∀q1 < length l. ∀q2 < length l. ((l ! q1) < (l ! q2)) ⟶ q1 < q2"
    using linorder_neqE_nat
    by (metis Groups.add_ac(2) add_mono_thms_linordered_field(5) less_irrefl) 
  have "(∃(x::real). (List.member l x ∧ a ≤ x ∧ x ≤ b)) ⟹ False"
  proof -
    assume "(∃(x::real). (List.member l x ∧ a ≤ x ∧ x ≤ b))"
    then obtain x where x_prop: "List.member l x ∧ a ≤ x ∧ x ≤ b" by auto
    then have l_prop: "List.member l x ∧ (l ! n) < x ∧ x < (l ! (n+1))"
      using lt_a b_lt by auto
    have nth_l: "l ! n < x" using l_prop by auto
    have np1th_l: "x < l ! (n+1)" using l_prop by auto
    have "∃k. k < length l ∧ nth l k = x" using l_prop
      by (meson in_set_member index_of_lookup(1) index_of_lookup(2))
    then obtain k where k_prop: "k < length l ∧ nth l k = x" by auto
    have n_lt: "n < k"
      using nth_l sorted_hyp_var k_prop add_lessD1 assms(2) linorder_neqE_nat nat_SN.gt_trans
      by (meson sorted_hyp_var2)
    have k_gt: "k < n + 1" 
      using sorted_hyp_var np1th_l k_prop
      using assms(2) sorted_hyp_var2 by blast 
    show False 
      using n_lt k_gt by auto
  qed
  then show ?thesis by auto
qed

(* This lemma shows that any zero of coprime_r is either between two roots or it's smaller than the 
least root (neg inf) or it's greater than the biggest root (pos inf). *)
lemma roots_of_coprime_r_location_property: 
  fixes qs:: "rat poly list"
  fixes sga:: "rat list"
  fixes zer_list
  assumes pairwise_rel_prime: "pairwise_coprime_list qs"
  assumes all_squarefree: "⋀q. q ∈ set qs ⟹ rsquarefree q"
  assumes x1_prop: "sga = sign_vec qs x1"
  assumes hnz: "has_no_zeros sga"
  assumes zer_list_prop: "zer_list = sorted_list_of_set {(x::real). (∃q ∈ set(qs). (rpoly q x = 0))}"
  shows "zer_list ≠ [] ⟶ ((x1 < (zer_list ! 0)) ∨ (x1 > (zer_list ! (length zer_list - 1)) ∨
    (∃ n < (length zer_list - 1). x1 > (zer_list ! n) ∧ x1 < (zer_list ! (n+1)))))"
proof - 
  let ?zer_list = "sorted_list_of_set {(x::real). (∃q ∈ set(qs). (rpoly q x = 0))} :: real list"
  show ?thesis 
  proof -
    have "((∀q. (List.member qs q) ⟶ q ≠ 0) ∧ has_no_zeros (sign_vec qs x1)) ⟹ ¬ List.member ?zer_list x1"
    proof (induct qs)
      case Nil
      then show ?case  apply (auto)
        by (simp add: member_rec(2)) 
    next
      case (Cons a qs)
      then show ?case 
      proof clarsimp 
        assume imp: "((∀q. List.member qs q ⟶ q ≠ 0) ∧
     has_no_zeros (sign_vec qs x1) ⟹
     ¬ List.member (sorted_list_of_set {x. ∃q∈set qs. rpoly q x = 0})
         x1)"
        assume nonz: "∀q. List.member (a # qs) q ⟶ q ≠ 0"
        assume hnz: " has_no_zeros (sign_vec (a # qs) x1)"
        assume mem_list: "List.member
     (sorted_list_of_set {x. rpoly a x = 0 ∨ (∃q∈set qs. rpoly q x = 0)})
     x1"
        have "has_no_zeros (sign_vec (a # qs) x1) ⟹ has_no_zeros (sign_vec qs x1)"
        proof - 
          assume hnz: "has_no_zeros (sign_vec (a # qs) x1)"
          have same_vec: "sign_vec (a # qs) x1 = ((if rpoly a x1 > 0 then 1 else if rpoly a x1 = 0 then 0 else -1) # sign_vec (qs) x1)"
            unfolding sign_vec_def squash_def by auto
          have "has_no_zeros ((if rpoly a x1 > 0 then 1 else if rpoly a x1 = 0 then 0 else -1) # sign_vec (qs) x1)
            ⟹  has_no_zeros (sign_vec (qs) x1)"
            by (simp add: has_no_zeros_def)
          then show "has_no_zeros (sign_vec qs x1)" using hnz same_vec by auto
        qed
        then have nmem: "¬ List.member (sorted_list_of_set {x. ∃q∈set qs. rpoly q x = 0}) x1"
          using hnz nonz imp apply (auto)
          by (simp add: member_rec(1))
        have "∀q ∈set qs. q ≠ 0"
          using nonz using in_set_member apply (auto) by fastforce 
        then have "∀q ∈set qs. finite {x. rpoly q x = 0}"
          by (simp add: poly_roots_finite)
        then have fin_set: "finite {x. ∃q∈set qs. rpoly q x = 0}"
          by auto
        have not_in: "x1 ∉ {x. ∃q∈set qs. rpoly q x = 0}" using fin_set nmem set_sorted_list_of_set
            all_squarefree
          apply (auto)
          by (simp add: List.member_def ‹finite {x. ∃q∈set qs. rpoly q x = 0}›)
        have x1_in: "x1 ∈ {x. rpoly a x = 0 ∨ (∃q∈set qs. rpoly q x = 0)}"
          using mem_list sorted_list_of_set 
        proof -
          have f1: "∀r R. ((r::real) ∈ R ∨ ¬ List.member (sorted_list_of_set R) r) ∨ infinite R"
            by (metis in_set_member set_sorted_list_of_set)
          have "finite {r. rpoly a (r::real) = 0}"
            by (metis (full_types) List.finite_set member_rec(1) nonz real_roots_of_rat_poly(1))
          then show ?thesis
            using f1 ‹finite {x. ∃q∈set qs. rpoly q x = 0}› mem_list by fastforce
        qed
        have "rpoly a x1 ≠ 0" using hnz 
          unfolding has_no_zeros_def sign_vec_def squash_def by auto
        then show "False" using not_in x1_in 
          by auto
      qed
    qed
    then have non_mem: "¬ List.member ?zer_list x1"
      using all_squarefree unfolding rsquarefree_def hnz apply (auto)
      using hnz x1_prop
      by (simp add: in_set_member) 
    have "?zer_list ≠ [] ⟹ ((x1 ≥ (?zer_list ! 0)) ∧ (x1 ≤ (?zer_list ! (length ?zer_list - 1))))
⟹ (∃ n < (length ?zer_list - 1). x1 > (?zer_list ! n) ∧ x1 < (?zer_list ! (n+1)))"
    proof - 
      assume "?zer_list ≠ []"
      assume x1_asm: "(x1 ≥ (?zer_list ! 0)) ∧ (x1 ≤ (?zer_list ! (length ?zer_list - 1)))"
      have nm1: "x1 ≠ ?zer_list ! 0" using non_mem
        using ‹sorted_list_of_set {x. ∃q∈set qs. rpoly q x = 0} ≠ []› in_set_member
        by (metis (no_types, lifting) in_set_conv_nth length_greater_0_conv)
      have nm2: "x1 ≠ ?zer_list ! (length ?zer_list -1)" using non_mem
        by (metis (no_types, lifting) One_nat_def ‹sorted_list_of_set {x. ∃q∈set qs. rpoly q x = 0} ≠ []› diff_Suc_less in_set_member length_greater_0_conv nth_mem) 
      then have x_asm_var: "x1 > (?zer_list ! 0) ∧ x1 < ?zer_list ! (length ?zer_list -1)"
        using x1_asm nm1 nm2 by auto
      have "(∀n. (n < (length ?zer_list - 1) ∧  x1 ≥ (?zer_list ! n) ⟶ x1 ≥ (?zer_list ! (n+1)))) ⟹ False"
      proof - 
        assume assump: "(∀n. (n < (length ?zer_list - 1) ∧  x1 ≥ (?zer_list ! n) ⟶ x1 ≥ (?zer_list ! (n+1))))"
        have zer_case: "x1 ≥ ?zer_list ! 0" using x_asm_var by auto