Session Count_Complex_Roots

Theory Extended_Sturm

(*
    Author:     Wenda Li <wl302@cam.ac.uk / liwenda1990@hotmail.com>
*)

section ‹An alternative Sturm sequences›

theory Extended_Sturm imports 
  "Sturm_Tarski.Sturm_Tarski" 
  "Winding_Number_Eval.Cauchy_Index_Theorem"
begin 
  
text ‹The main purpose of this theory is to provide an effective way to compute 
  @{term "cindexE a b f"} when @{term f} is a rational function. The idea is similar to and based
  on the evaluation of @{term cindex_poly} through @{thm cindex_poly_changes_itv_mods}.›   
  
text ‹
This alternative version of remainder sequences is inspired by the paper 
  "The Fundamental Theorem of Algebra made effective: 
  an elementary real-algebraic proof via Sturm chains" 
by Michael Eisermann.
›  
  
hide_const Permutations.sign  
  
subsection ‹Misc›
   
lemma is_unit_pCons_ex_iff:
  fixes p::"'a::field poly"
  shows "is_unit p  (a. a0  p=[:a:])"
  using is_unit_poly_iff is_unit_triv by auto

lemma poly_gcd_iff: 
  "poly (gcd p q) x=0  poly p x=0  poly q x=0 "
  by (simp add: poly_eq_0_iff_dvd)
  
lemma eventually_poly_nz_at_within:
  fixes x :: "'a::{idom,euclidean_space} "
  assumes "p0" 
  shows "eventually (λx. poly p x0) (at x within S)"     
proof -
  have "eventually (λx. poly p x0) (at x within S) 
      = (F x in (at x within S). yproots p. x  y)"
    apply (rule eventually_subst,rule eventuallyI)
    by (auto simp add:proots_def)
  also have "... = (yproots p. F x in (at x within S). x  y)"
    apply (subst eventually_ball_finite_distrib)
    using p0 by auto
  also have "..."
    unfolding eventually_at
    by (metis gt_ex not_less_iff_gr_or_eq zero_less_dist_iff) 
  finally show ?thesis .
qed
    
lemma sgn_power:
  fixes x::"'a::linordered_idom"
  shows "sgn (x^n) = (if n=0 then 1 else if even n then ¦sgn x¦ else sgn x)"
  apply (induct n)
  by (auto simp add:sgn_mult)

lemma poly_divide_filterlim_at_top: 
  fixes p q::"real poly"
  defines "ll( if degree q<degree p then 
                    at 0 
                else if degree q=degree p then 
                    nhds (lead_coeff q / lead_coeff p)
                else if sgn_pos_inf q * sgn_pos_inf p > 0 then 
                    at_top
                else 
                    at_bot)"
  assumes "p0" "q0"
  shows "filterlim (λx. poly q x / poly p x) ll at_top"
proof -
  define pp where "pp=(λx. poly p x / x^(degree p))"    
  define qq where "qq=(λx. poly q x / x^(degree q))"
  define dd where "dd=(λx::real. if degree p>degree q then 1/x^(degree p - degree q) else 
                                x^(degree q - degree p))"
  have divide_cong:"Fx in at_top. poly q x / poly p x = qq x / pp x * dd x"
  proof (rule eventually_at_top_linorderI[of 1])
    fix x assume "(x::real)1"
    then have "x0" by auto
    then show "poly q x / poly p x = qq x / pp x * dd x"
      unfolding qq_def pp_def dd_def using assms 
      by (auto simp add:field_simps power_diff) 
  qed
  have qqpp_tendsto:"((λx. qq x / pp x)  lead_coeff q / lead_coeff p) at_top"
  proof -
    have "(qq  lead_coeff q) at_top"
      unfolding qq_def using poly_divide_tendsto_aux[of q]  
      by (auto elim!:filterlim_mono simp:at_top_le_at_infinity)
    moreover have "(pp  lead_coeff p) at_top"
      unfolding pp_def using poly_divide_tendsto_aux[of p]  
      by (auto elim!:filterlim_mono simp:at_top_le_at_infinity)
    ultimately show ?thesis using p0 by (auto intro!:tendsto_eq_intros)
  qed
  
  have ?thesis when "degree q<degree p"
  proof -
    have "filterlim (λx. poly q x / poly p x) (at 0) at_top"
    proof (rule filterlim_atI)
      show "((λx. poly q x / poly p x)  0) at_top"
        using poly_divide_tendsto_0_at_infinity[OF that]
        by (auto elim:filterlim_mono simp:at_top_le_at_infinity)
      have "F x in at_top. poly q x 0" "F x in at_top. poly p x 0"
        using poly_eventually_not_zero[OF q0] poly_eventually_not_zero[OF p0]
              filter_leD[OF at_top_le_at_infinity]
        by auto
      then show "F x in at_top. poly q x / poly p x  0"
        apply eventually_elim
        by auto
    qed
    then show ?thesis unfolding ll_def using that by auto
  qed
  moreover have ?thesis when "degree q=degree p"
  proof -
    have "((λx. poly q x / poly p x)  lead_coeff q / lead_coeff p) at_top"
      using divide_cong qqpp_tendsto that unfolding dd_def
      by (auto dest:tendsto_cong)
    then show ?thesis unfolding ll_def using that by auto
  qed
  moreover have ?thesis when "degree q>degree p" "sgn_pos_inf q * sgn_pos_inf p > 0"
  proof -
    have "filterlim (λx. (qq x / pp x) * dd x) at_top at_top"
    proof (subst filterlim_tendsto_pos_mult_at_top_iff[OF qqpp_tendsto])
      show "0 < lead_coeff q / lead_coeff p" using that(2) unfolding sgn_pos_inf_def
        by (simp add: zero_less_divide_iff zero_less_mult_iff)
      show "filterlim dd at_top at_top"
        unfolding dd_def using that(1) 
        by (auto intro!:filterlim_pow_at_top simp:filterlim_ident)
    qed
    then have "LIM x at_top. poly q x / poly p x :> at_top" 
      using filterlim_cong[OF _ _ divide_cong] by blast
    then show ?thesis unfolding ll_def using that by auto
  qed
  moreover have ?thesis  when "degree q>degree p" "¬ sgn_pos_inf q * sgn_pos_inf p > 0"
  proof -
    have "filterlim (λx. (qq x / pp x) * dd x) at_bot at_top"
    proof (subst filterlim_tendsto_neg_mult_at_bot_iff[OF qqpp_tendsto])
      show "lead_coeff q / lead_coeff p < 0" 
        using that(2) p0 q0 unfolding sgn_pos_inf_def
        by (metis divide_eq_0_iff divide_sgn leading_coeff_0_iff 
            linorder_neqE_linordered_idom sgn_divide sgn_greater)
      show "filterlim dd at_top at_top"
        unfolding dd_def using that(1) 
        by (auto intro!:filterlim_pow_at_top simp:filterlim_ident)
    qed
    then have "LIM x at_top. poly q x / poly p x :> at_bot" 
      using filterlim_cong[OF _ _ divide_cong] by blast
    then show ?thesis unfolding ll_def using that by auto
  qed
  ultimately show ?thesis by linarith
qed

lemma poly_divide_filterlim_at_bot: 
  fixes p q::"real poly"
  defines "ll( if degree q<degree p then 
                    at 0 
                else if degree q=degree p then 
                    nhds (lead_coeff q / lead_coeff p)
                else if sgn_neg_inf q * sgn_neg_inf p > 0 then 
                    at_top
                else 
                    at_bot)"
  assumes "p0" "q0"
  shows "filterlim (λx. poly q x / poly p x) ll at_bot"
proof -
  define pp where "pp=(λx. poly p x / x^(degree p))"    
  define qq where "qq=(λx. poly q x / x^(degree q))"
  define dd where "dd=(λx::real. if degree p>degree q then 1/x^(degree p - degree q) else 
                                x^(degree q - degree p))"
  have divide_cong:"Fx in at_bot. poly q x / poly p x = qq x / pp x * dd x"
  proof (rule eventually_at_bot_linorderI[of "-1"])
    fix x assume "(x::real)-1"
    then have "x0" by auto
    then show "poly q x / poly p x = qq x / pp x * dd x"
      unfolding qq_def pp_def dd_def using assms 
      by (auto simp add:field_simps power_diff) 
  qed
  have qqpp_tendsto:"((λx. qq x / pp x)  lead_coeff q / lead_coeff p) at_bot"
  proof -
    have "(qq  lead_coeff q) at_bot"
      unfolding qq_def using poly_divide_tendsto_aux[of q]  
      by (auto elim!:filterlim_mono simp:at_bot_le_at_infinity)
    moreover have "(pp  lead_coeff p) at_bot"
      unfolding pp_def using poly_divide_tendsto_aux[of p]  
      by (auto elim!:filterlim_mono simp:at_bot_le_at_infinity)
    ultimately show ?thesis using p0 by (auto intro!:tendsto_eq_intros)
  qed
  
  have ?thesis when "degree q<degree p"
  proof -
    have "filterlim (λx. poly q x / poly p x) (at 0) at_bot"
    proof (rule filterlim_atI)
      show "((λx. poly q x / poly p x)  0) at_bot"
        using poly_divide_tendsto_0_at_infinity[OF that]
        by (auto elim:filterlim_mono simp:at_bot_le_at_infinity)
      have "F x in at_bot. poly q x 0" "F x in at_bot. poly p x 0"
        using poly_eventually_not_zero[OF q0] poly_eventually_not_zero[OF p0]
              filter_leD[OF at_bot_le_at_infinity]
        by auto
      then show "F x in at_bot. poly q x / poly p x  0"
        by eventually_elim auto
    qed
    then show ?thesis unfolding ll_def using that by auto
  qed
  moreover have ?thesis when "degree q=degree p"
  proof -
    have "((λx. poly q x / poly p x)  lead_coeff q / lead_coeff p) at_bot"
      using divide_cong qqpp_tendsto that unfolding dd_def
      by (auto dest:tendsto_cong)
    then show ?thesis unfolding ll_def using that by auto
  qed
  moreover have ?thesis when "degree q>degree p" "sgn_neg_inf q * sgn_neg_inf p > 0"
  proof -
    define cc where "cc=lead_coeff q / lead_coeff p"
    have "(cc > 0  even (degree q - degree p))  (cc<0  odd (degree q - degree p))"
    proof -
      have "even (degree q - degree p)  
            (even (degree q)  even (degree p))  (odd (degree q)  odd (degree p))"
        using ‹degree q>degree p by auto
      then show ?thesis
        using that p0 q0 unfolding sgn_neg_inf_def cc_def zero_less_mult_iff 
          divide_less_0_iff zero_less_divide_iff 
        apply (simp add:if_split[of "(<) 0"] if_split[of "(>) 0"])
        by argo
    qed
    moreover have "filterlim (λx. (qq x / pp x) * dd x) at_top at_bot"
      when "cc>0" "even (degree q - degree p)"
    proof (subst filterlim_tendsto_pos_mult_at_top_iff[OF qqpp_tendsto])
      show "0 < lead_coeff q / lead_coeff p" using cc>0 unfolding cc_def by auto
      show "filterlim dd at_top at_bot"
        unfolding dd_def using ‹degree q>degree p that(2)
        by (auto intro!:filterlim_pow_at_bot_even simp:filterlim_ident)
    qed
    moreover have "filterlim (λx. (qq x / pp x) * dd x) at_top at_bot"
      when "cc<0" "odd (degree q - degree p)"
    proof (subst filterlim_tendsto_neg_mult_at_top_iff[OF qqpp_tendsto])
      show "0 > lead_coeff q / lead_coeff p" using cc<0 unfolding cc_def by auto
      show "filterlim dd at_bot at_bot"
        unfolding dd_def using ‹degree q>degree p that(2)
        by (auto intro!:filterlim_pow_at_bot_odd simp:filterlim_ident)
    qed
    ultimately have "filterlim (λx. (qq x / pp x) * dd x) at_top at_bot"
      by blast
    then have "LIM x at_bot. poly q x / poly p x :> at_top"
      using filterlim_cong[OF _ _ divide_cong] by blast
    then show ?thesis unfolding ll_def using that by auto
  qed
  moreover have ?thesis  when "degree q>degree p" "¬ sgn_neg_inf q * sgn_neg_inf p > 0"
  proof -
    define cc where "cc=lead_coeff q / lead_coeff p"
    have "(cc < 0  even (degree q - degree p))  (cc > 0  odd (degree q - degree p))"
    proof -
      have "even (degree q - degree p)  
            (even (degree q)  even (degree p))  (odd (degree q)  odd (degree p))"
        using ‹degree q>degree p by auto
      then show ?thesis
        using that p0 q0 unfolding sgn_neg_inf_def cc_def zero_less_mult_iff 
          divide_less_0_iff zero_less_divide_iff
        apply (simp add:if_split[of "(<) 0"] if_split[of "(>) 0"])
        by (metis leading_coeff_0_iff linorder_neqE_linordered_idom)
    qed
    moreover have "filterlim (λx. (qq x / pp x) * dd x) at_bot at_bot"
      when "cc<0" "even (degree q - degree p)"
    proof (subst filterlim_tendsto_neg_mult_at_bot_iff[OF qqpp_tendsto])
      show "0 > lead_coeff q / lead_coeff p" using cc<0 unfolding cc_def by auto
      show "filterlim dd at_top at_bot"
        unfolding dd_def using ‹degree q>degree p that(2)
        by (auto intro!:filterlim_pow_at_bot_even simp:filterlim_ident)
    qed
    moreover have "filterlim (λx. (qq x / pp x) * dd x) at_bot at_bot"
      when "cc>0" "odd (degree q - degree p)"
    proof (subst filterlim_tendsto_pos_mult_at_bot_iff[OF qqpp_tendsto])
      show "0 < lead_coeff q / lead_coeff p" using cc>0 unfolding cc_def by auto
      show "filterlim dd at_bot at_bot"
        unfolding dd_def using ‹degree q>degree p that(2)
        by (auto intro!:filterlim_pow_at_bot_odd simp:filterlim_ident)
    qed
    ultimately have "filterlim (λx. (qq x / pp x) * dd x) at_bot at_bot"
      by blast
    then have "LIM x at_bot. poly q x / poly p x :> at_bot"
      using filterlim_cong[OF _ _ divide_cong] by blast
    then show ?thesis unfolding ll_def using that by auto
  qed
  ultimately show ?thesis by linarith
qed
        
subsection ‹Alternative definition of cross›
 
definition cross_alt :: "real poly real poly  real  real  int" where
  "cross_alt p q a b=¦sign (poly p a) - sign (poly q a)¦ - ¦sign (poly p b) - sign(poly q b)¦"  

lemma cross_alt_coprime_0:
  assumes "coprime p q" "p=0q=0"
  shows "cross_alt p q a b=0" 
proof -
  have ?thesis when "p=0" 
  proof -
    have "is_unit q" using that ‹coprime p q 
      by simp
    then obtain a where "a0" "q=[:a:]" using is_unit_pCons_ex_iff by blast
    then show ?thesis using that unfolding cross_alt_def by auto
  qed
  moreover have ?thesis when "q=0" 
  proof -
    have "is_unit p" using that ‹coprime p q 
      by simp
    then obtain a where "a0" "p=[:a:]" using is_unit_pCons_ex_iff by blast
    then show ?thesis using that unfolding cross_alt_def by auto
  qed 
  ultimately show ?thesis using p=0q=0 by auto
qed  
  
lemma cross_alt_0[simp]: "cross_alt 0 0 a b=0" unfolding cross_alt_def by auto 
  
lemma cross_alt_poly_commute:
  "cross_alt p q a b = cross_alt q p a b"
  unfolding cross_alt_def by auto
    
lemma cross_alt_clear_n:
  assumes "coprime p q"
  shows "cross_alt p q a b = cross_alt 1 (p*q) a b"
proof -
  have "¦sign (poly p a) - sign (poly q a)¦  = ¦1 - sign (poly p a) * sign (poly q a)¦"
  proof (cases "poly p a=0  poly q a=0")
    case True
    then have False using assms using coprime_poly_0 by blast
    then show ?thesis by auto
  next
    case False
    then show ?thesis 
      unfolding Sturm_Tarski.sign_def
      by force
  qed
  moreover have "¦sign (poly p b) - sign (poly q b)¦  = ¦1 - sign (poly p b) * sign (poly q b)¦"
  proof (cases "poly p b=0  poly q b=0")
    case True
    then have False using assms using coprime_poly_0 by blast
    then show ?thesis by auto
  next
    case False
    then show ?thesis 
      unfolding Sturm_Tarski.sign_def
      by force
  qed  
  ultimately show ?thesis
    by (simp add: cross_alt_def sign_times)
qed   
  
subsection ‹Alternative sign variation sequencse›  
  
fun changes_alt:: "('a ::linordered_idom) list  int" where
  "changes_alt [] = 0"|
  "changes_alt [_] = 0" |
  "changes_alt (x1#x2#xs) = abs(sign x1 - sign x2) + changes_alt (x2#xs)"  
  
definition changes_alt_poly_at::"('a ::linordered_idom) poly list  'a  int" where
  "changes_alt_poly_at ps a= changes_alt (map (λp. poly p a) ps)"

definition changes_alt_itv_smods:: "real  real real poly  real poly   int" where
  "changes_alt_itv_smods a b p q= (let ps= smods p q 
                                    in changes_alt_poly_at ps a - changes_alt_poly_at ps b)"  
 
lemma changes_alt_itv_smods_rec:
  assumes "a<b" "coprime p q" 
  shows "changes_alt_itv_smods a b p q  = cross_alt p q a b + changes_alt_itv_smods a b q (-(p mod q))"
proof (cases "p = 0  q = 0  q dvd p")
  case True
  moreover have "p=0  q=0  ?thesis"
    using cross_alt_coprime_0[OF ‹coprime p q] 
    unfolding changes_alt_itv_smods_def changes_alt_poly_at_def by fastforce
  moreover have "p0;q0;p mod q = 0  ?thesis"  
    unfolding changes_alt_itv_smods_def changes_alt_poly_at_def cross_alt_def
    by (simp add:sign_times)
  ultimately show ?thesis
    by auto (auto elim: dvdE)
next
  case False
  hence "p0" "q0" "p mod q0" by auto
  then obtain ps where ps:"smods p q=p#q#-(p mod q)#ps" "smods q (-(p mod q)) = q#-(p mod q)#ps"
    by auto
  define changes_diff where "changes_diffλx. changes_alt_poly_at (p#q#-(p mod q)#ps) x 
    - changes_alt_poly_at (q#-(p mod q)#ps) x"
  have "changes_diff a - changes_diff b=cross_alt p q a b" 
    unfolding changes_diff_def changes_alt_poly_at_def cross_alt_def by simp
  thus ?thesis unfolding changes_alt_itv_smods_def changes_diff_def changes_alt_poly_at_def ps 
    by force
qed  
  
subsection ‹jumpF on polynomials›

definition jumpF_polyR:: "real poly  real poly  real  real" where
  "jumpF_polyR q p a = jumpF (λx. poly q x / poly p x) (at_right a)"
  
definition jumpF_polyL:: "real poly  real poly  real  real" where
  "jumpF_polyL q p a = jumpF (λx. poly q x / poly p x) (at_left a)"

definition jumpF_poly_top:: "real poly  real poly  real" where
  "jumpF_poly_top q p = jumpF (λx. poly q x / poly p x) at_top"

definition jumpF_poly_bot:: "real poly  real poly  real" where
  "jumpF_poly_bot q p = jumpF (λx. poly q x / poly p x) at_bot"

  
lemma jumpF_polyR_0[simp]: "jumpF_polyR 0 p a = 0" "jumpF_polyR q 0 a = 0" 
  unfolding jumpF_polyR_def by auto
    
lemma jumpF_polyL_0[simp]: "jumpF_polyL 0 p a = 0" "jumpF_polyL q 0 a = 0" 
  unfolding jumpF_polyL_def by auto
    
lemma jumpF_polyR_mult_cancel:
  assumes "p'0"
  shows "jumpF_polyR (p' * q) (p' * p) a = jumpF_polyR q p a"
unfolding jumpF_polyR_def
proof (rule jumpF_cong)
  obtain ub where "a < ub" "z. a < z  z  ub  poly p' z  0"
    using next_non_root_interval[OF p'0,of a] by auto
  then show "F x in at_right a. poly (p' * q) x / poly (p' * p) x = poly q x / poly p x"
    apply (unfold eventually_at_right)
    apply (intro exI[where x=ub])
    by auto
qed simp
  
lemma jumpF_polyL_mult_cancel:
  assumes "p'0"
  shows "jumpF_polyL (p' * q) (p' * p) a = jumpF_polyL q p a"
unfolding jumpF_polyL_def
proof (rule jumpF_cong)
  obtain lb where "lb < a" "z. lb  z  z < a  poly p' z  0 "
    using last_non_root_interval[OF p'0,of a] by auto
  then show "F x in at_left a. poly (p' * q) x / poly (p' * p) x = poly q x / poly p x"
    apply (unfold eventually_at_left)
    apply (intro exI[where x=lb])
    by auto
qed simp  
      
lemma jumpF_poly_noroot: 
  assumes "poly p a0"
  shows "jumpF_polyL q p a = 0" "jumpF_polyR q p a = 0" 
  subgoal unfolding jumpF_polyL_def using assms
    apply (intro jumpF_not_infinity)
    by (auto intro!:continuous_intros)  
  subgoal unfolding jumpF_polyR_def using assms
    apply (intro jumpF_not_infinity)
    by (auto intro!:continuous_intros)
  done
  

lemma jumpF_polyR_coprime:
  assumes "coprime p q"
  shows "jumpF_polyR q p x = (if p  0  q  0  poly p x=0 then 
                                if sign_r_pos p x  poly q x>0 then 1/2 else - 1/2 else 0)"
proof (cases "p=0  q=0  poly p x0")
  case True
  then show ?thesis using jumpF_poly_noroot by fastforce
next
  case False
  then have asm:"p0" "q0" "poly p x=0" by auto  
  then have "poly q x0" using assms using coprime_poly_0 by blast
  have ?thesis when "sign_r_pos p x  poly q x>0"
  proof -
    have "(poly p has_sgnx sgn (poly q x)) (at_right x)"
      by (metis False ‹poly q x  0 add.inverse_neutral has_sgnx_imp_sgnx less_not_sym 
          neg_less_iff_less poly_has_sgnx_values(2) sgn_if sign_r_pos_sgnx_iff that 
          trivial_limit_at_right_real zero_less_one)
    then have "LIM x at_right x. poly q x / poly p x :> at_top"    
      apply (subst filterlim_divide_at_bot_at_top_iff[of _ "poly q x"])
      apply (auto simp add:‹poly q x0)
      by (metis asm(3) poly_tendsto(3))
    then have "jumpF_polyR q p x = 1/2"
      unfolding jumpF_polyR_def jumpF_def by auto
    then show ?thesis using that False by auto
  qed
  moreover have ?thesis when "¬ (sign_r_pos p x  poly q x>0)"
  proof -
    have "(poly p has_sgnx - sgn (poly q x)) (at_right x)"
    proof -
      have "(0::real) < 1  ¬ (1::real) < 0  sign_r_pos p x 
           (poly p has_sgnx - sgn (poly q x)) (at_right x)"
        by simp
      then show ?thesis
        by (metis (no_types) False ‹poly q x  0 add.inverse_inverse has_sgnx_imp_sgnx 
            neg_less_0_iff_less poly_has_sgnx_values(2) sgn_if sgn_less sign_r_pos_sgnx_iff 
            that trivial_limit_at_right_real)
    qed
    then have "LIM x at_right x. poly q x / poly p x :> at_bot"    
      apply (subst filterlim_divide_at_bot_at_top_iff[of _ "poly q x"])
      apply (auto simp add:‹poly q x0)
      by (metis asm(3) poly_tendsto(3)) 
    then have "jumpF_polyR q p x = - 1/2"
      unfolding jumpF_polyR_def jumpF_def by auto
    then show ?thesis using that False by auto  
  qed
  ultimately show ?thesis by auto
qed

lemma jumpF_polyL_coprime:
  assumes "coprime p q"
  shows "jumpF_polyL q p x = (if p  0  q  0  poly p x=0 then 
                if even (order x p)  sign_r_pos p x  poly q x>0 then 1/2 else - 1/2 else 0)"  
proof (cases "p=0  q=0  poly p x0")
  case True
  then show ?thesis using jumpF_poly_noroot by fastforce
next  
  case False
  then have asm:"p0" "q0" "poly p x=0" by auto  
  then have "poly q x0" using assms using coprime_poly_0 by blast
  have ?thesis when "even (order x p)  sign_r_pos p x  poly q x>0"
  proof -
    consider (lt) "poly q x>0" | (gt) " poly q x<0" using ‹poly q x0 by linarith
    then have "sgnx (poly p) (at_left x) = sgn (poly q x)"
      apply cases
      subgoal using that sign_r_pos_sgnx_iff poly_sgnx_values[OF p0,of x]
        apply (subst poly_sgnx_left_right[OF p0])
        by auto
      subgoal using that sign_r_pos_sgnx_iff poly_sgnx_values[OF p0,of x]
        apply (subst poly_sgnx_left_right[OF p0])
        by auto
      done
    then have "(poly p has_sgnx sgn (poly q x)) (at_left x)"
      by (metis sgnx_able_poly(2) sgnx_able_sgnx)
    then have "LIM x at_left x. poly q x / poly p x :> at_top"    
      apply (subst filterlim_divide_at_bot_at_top_iff[of _ "poly q x"])
      apply (auto simp add:‹poly q x0)
      by (metis asm(3) poly_tendsto(2))   
    then have "jumpF_polyL q p x = 1/2"
      unfolding jumpF_polyL_def jumpF_def by auto
    then show ?thesis using that False by auto
  qed
  moreover have ?thesis when "¬ (even (order x p)  sign_r_pos p x  poly q x>0)"
  proof -
    consider (lt) "poly q x>0" | (gt) " poly q x<0" using ‹poly q x0 by linarith
    then have "sgnx (poly p) (at_left x) = - sgn (poly q x)"
      apply cases
      subgoal using that sign_r_pos_sgnx_iff poly_sgnx_values[OF p0,of x]
        apply (subst poly_sgnx_left_right[OF p0])
        by auto
      subgoal using that sign_r_pos_sgnx_iff poly_sgnx_values[OF p0,of x]
        apply (subst poly_sgnx_left_right[OF p0])
        by auto
      done
    then have "(poly p has_sgnx - sgn (poly q x)) (at_left x)"
      by (metis sgnx_able_poly(2) sgnx_able_sgnx)
    then have "LIM x at_left x. poly q x / poly p x :> at_bot"    
      apply (subst filterlim_divide_at_bot_at_top_iff[of _ "poly q x"])
      apply (auto simp add:‹poly q x0)
      by (metis asm(3) poly_tendsto(2))   
    then have "jumpF_polyL q p x = - 1/2"
      unfolding jumpF_polyL_def jumpF_def by auto
    then show ?thesis using that False by auto 
  qed
  ultimately show ?thesis by auto
qed    
    
lemma jumpF_times:
  assumes tendsto:"(f  c) F" and "c0" "Fbot"
  shows "jumpF (λx. f x * g x) F = sgn c * jumpF g F"  
proof -
  have "c>0  c<0" using c0 by auto
  moreover have ?thesis when "c>0"
  proof -
    note filterlim_tendsto_pos_mult_at_top_iff[OF tendsto c>0,of g]
    moreover note filterlim_tendsto_pos_mult_at_bot_iff[OF tendsto c>0,of g]
    moreover have "sgn c = 1" using c>0 by auto
    ultimately show ?thesis unfolding jumpF_def by auto
  qed
  moreover have ?thesis when "c<0"
  proof -
    define atbot where "atbot = filterlim g at_bot F"
    define attop where "attop = filterlim g at_top F"    
    have "jumpF (λx. f x * g x) F = (if atbot then 1 / 2 else if attop then - 1 / 2 else 0)"
    proof -
      note filterlim_tendsto_neg_mult_at_top_iff[OF tendsto c<0,of g]
      moreover note filterlim_tendsto_neg_mult_at_bot_iff[OF tendsto c<0,of g]
      ultimately show ?thesis unfolding jumpF_def atbot_def attop_def by auto
    qed
    also have "... = - (if attop then 1 / 2 else if atbot then - 1 / 2 else 0)"
    proof -
      have False when atbot attop
        using filterlim_at_top_at_bot[OF _ _ Fbot›] that unfolding atbot_def attop_def by auto
      then show ?thesis by fastforce    
    qed
    also have "... = sgn c * jumpF g F"
      using c<0 unfolding jumpF_def attop_def atbot_def by auto
    finally show ?thesis .
  qed
  ultimately show ?thesis by auto
qed

  
lemma jumpF_polyR_inverse_add:
  assumes "coprime p q"
  shows "jumpF_polyR q p x + jumpF_polyR p q x = jumpF_polyR 1 (q*p) x"
proof (cases "p=0  q=0")
  case True
  then show ?thesis by auto
next
  case False
  have jumpF_add:
    "jumpF_polyR q p x= jumpF_polyR 1 (q*p) x" when "poly p x=0" "coprime p q" for p q
  proof (cases "p=0")
    case True
    then show ?thesis by auto
  next
    case False
    have "poly q x0" using that coprime_poly_0 by blast
    then have "q0" by auto
    moreover have "sign_r_pos p x = (0 < poly q x)  sign_r_pos (q * p) x"
      using sign_r_pos_mult[OF q0 p0] sign_r_pos_rec[OF q0] ‹poly q x0
      by auto
    ultimately show ?thesis using ‹poly p x=0  
      unfolding jumpF_polyR_coprime[OF ‹coprime p q,of x] jumpF_polyR_coprime[of "q*p" 1 x,simplified]
      by auto
  qed
  have False when "poly p x=0" "poly q x=0" 
    using ‹coprime p q that coprime_poly_0 by blast
  moreover have ?thesis when "poly p x=0" "poly q x0" 
  proof -
    have "jumpF_polyR p q x = 0" using jumpF_poly_noroot[OF ‹poly q x0] by auto
    then show ?thesis using jumpF_add[OF ‹poly p x=0 ‹coprime p q] by auto
  qed
  moreover have ?thesis when "poly p x0" "poly q x=0" 
  proof -
    have "jumpF_polyR q p x = 0" using jumpF_poly_noroot[OF ‹poly p x0] by auto
    then show ?thesis using jumpF_add[OF ‹poly q x=0,of p] ‹coprime p q 
      by (simp add: ac_simps)
  qed  
  moreover have ?thesis when "poly p x0" "poly q x0" 
    by (simp add: jumpF_poly_noroot(2) that(1) that(2))
  ultimately show ?thesis by auto
qed

lemma jumpF_polyL_inverse_add:
  assumes "coprime p q"
  shows "jumpF_polyL q p x + jumpF_polyL p q x = jumpF_polyL 1 (q*p) x"  
proof (cases "p=0  q=0")
  case True
  then show ?thesis by auto
next
  case False
  have jumpF_add:
    "jumpF_polyL q p x= jumpF_polyL 1 (q*p) x" when "poly p x=0" "coprime p q" for p q
  proof (cases "p=0")
    case True
    then show ?thesis by auto
  next
    case False
    have "poly q x0" using that coprime_poly_0 by blast
    then have "q0" by auto
    moreover have "sign_r_pos p x = (0 < poly q x)  sign_r_pos (q * p) x"
      using sign_r_pos_mult[OF q0 p0] sign_r_pos_rec[OF q0] ‹poly q x0
      by auto
    moreover have "order x p = order x (q * p)" 
      by (metis ‹poly q x  0 add_cancel_right_left divisors_zero order_mult order_root)
    ultimately show ?thesis using ‹poly p x=0  
      unfolding jumpF_polyL_coprime[OF ‹coprime p q,of x] jumpF_polyL_coprime[of "q*p" 1 x,simplified]
      by auto
  qed
  have False when "poly p x=0" "poly q x=0" 
    using ‹coprime p q that coprime_poly_0 by blast
  moreover have ?thesis when "poly p x=0" "poly q x0" 
  proof -
    have "jumpF_polyL p q x = 0" using jumpF_poly_noroot[OF ‹poly q x0] by auto
    then show ?thesis using jumpF_add[OF ‹poly p x=0 ‹coprime p q] by auto
  qed
  moreover have ?thesis when "poly p x0" "poly q x=0" 
  proof -
    have "jumpF_polyL q p x = 0" using jumpF_poly_noroot[OF ‹poly p x0] by auto
    then show ?thesis using jumpF_add[OF ‹poly q x=0,of p] ‹coprime p q 
      by (simp add: ac_simps)
  qed  
  moreover have ?thesis when "poly p x0" "poly q x0" 
    by (simp add: jumpF_poly_noroot that(1) that(2))
  ultimately show ?thesis by auto
qed    
  

lemma jumpF_polyL_smult_1:
  "jumpF_polyL (smult c q) p x = sgn c * jumpF_polyL q p x"
proof (cases "c=0")
  case True
  then show ?thesis by auto
next
  case False
  then show ?thesis 
    unfolding jumpF_polyL_def 
    apply (subst jumpF_times[of "λ_. c",symmetric])
    by auto
qed  

lemma jumpF_polyR_smult_1:
  "jumpF_polyR (smult c q) p x = sgn c * jumpF_polyR q p x"
proof (cases "c=0")
  case True
  then show ?thesis by auto
next
  case False
  then show ?thesis
    unfolding jumpF_polyR_def using False 
    apply (subst jumpF_times[of "λ_. c",symmetric])
    by auto
qed  
  

lemma
  shows jumpF_polyR_mod:"jumpF_polyR q p x = jumpF_polyR (q mod p) p x" and
        jumpF_polyL_mod:"jumpF_polyL q p x = jumpF_polyL (q mod p) p x"
proof -
  define f where "f=(λx. poly (q div p) x)"
  define g where "g=(λx. poly (q mod p) x / poly p x)"
  have jumpF_eq:"jumpF (λx. poly q x / poly p x) (at y within S) = jumpF g (at y within S)" 
    when "p0" for y S
  proof -
    let ?F = "at y within S"
    have "F x in at y within S. poly p x  0" 
      using eventually_poly_nz_at_within[OF p0,of y S] .
    then have "eventually (λx. (poly q x / poly p x) = (f x+ g x)) ?F"
    proof (rule eventually_mono)
      fix x
      assume P: "poly p x  0"
      have "poly q x = poly (q div p * p + q mod p) x"
        by simp
      also have " = f x * poly p x + poly (q mod p) x"
        by (simp only: poly_add poly_mult f_def g_def)
      moreover have "poly (q mod p) x = g x * poly p x"
        using P by (simp add: g_def)
      ultimately show "poly q x / poly p x = f x + g x"
        using P by simp
    qed
    then have "jumpF (λx. poly q x / poly p x) ?F = jumpF (λx. f x+ g x) ?F"
      by (intro jumpF_cong,auto)
    also have "... = jumpF g ?F"  
    proof -
      have "(f  f y) (at y within S)"
        unfolding f_def by (intro tendsto_intros)  
      from filterlim_tendsto_add_at_bot_iff[OF this,of g] filterlim_tendsto_add_at_top_iff[OF this,of g] 
      show ?thesis unfolding jumpF_def by auto
    qed
    finally show ?thesis .
  qed
  show "jumpF_polyR q p x = jumpF_polyR (q mod p) p x"
    apply (cases "p=0")
    subgoal by auto
    subgoal using jumpF_eq unfolding g_def jumpF_polyR_def by auto
    done
  show "jumpF_polyL q p x = jumpF_polyL (q mod p) p x"
    apply (cases "p=0")
    subgoal by auto
    subgoal using jumpF_eq unfolding g_def jumpF_polyL_def by auto
    done  
qed 


lemma jumpF_poly_top_0[simp]: "jumpF_poly_top 0 p = 0" "jumpF_poly_top q 0 = 0"
  unfolding jumpF_poly_top_def by auto

lemma jumpF_poly_bot_0[simp]: "jumpF_poly_bot 0 p = 0" "jumpF_poly_bot q 0 = 0"
  unfolding jumpF_poly_bot_def by auto

lemma jumpF_poly_top_code:
  "jumpF_poly_top q p = (if p0  q0  degree q>degree p then 
          if sgn_pos_inf q * sgn_pos_inf p > 0 then 1/2 else -1/2 else 0)"
proof (cases "p0  q0  degree q>degree p")
  case True
  have ?thesis when "sgn_pos_inf q * sgn_pos_inf p > 0"
  proof -
    have "LIM x at_top. poly q x / poly p x :> at_top"
      using poly_divide_filterlim_at_top[of p q] True that by auto
    then have "jumpF (λx. poly q x / poly p x) at_top = 1/2"
      unfolding jumpF_def by auto
    then show ?thesis unfolding jumpF_poly_top_def using that True by auto
  qed
  moreover have ?thesis when "¬ sgn_pos_inf q * sgn_pos_inf p > 0"
  proof -
    have "LIM x at_top. poly q x / poly p x :> at_bot"
      using poly_divide_filterlim_at_top[of p q] True that by auto
    then have "jumpF (λx. poly q x / poly p x) at_top = - 1/2"
      unfolding jumpF_def by auto
    then show ?thesis unfolding jumpF_poly_top_def using that True by auto
  qed
  ultimately show ?thesis by auto
next
  case False
  define P where "P= (¬ (LIM x at_top. poly q x / poly p x :> at_bot) 
                       ¬ (LIM x at_top. poly q x / poly p x :> at_top))"
  have P when "p=0  q=0"
    unfolding P_def using that 
    by (auto elim!:filterlim_at_bot_nhds filterlim_at_top_nhds)
  moreover have P when "p0" "q0" "degree p> degree q"
  proof -
    have "LIM x at_top. poly q x / poly p x :> at 0"
      using poly_divide_filterlim_at_top[OF that(1,2)] that(3) by auto
    then show ?thesis unfolding P_def 
      by (auto elim!:filterlim_at_bot_nhds filterlim_at_top_nhds simp:filterlim_at)
  qed 
  moreover have P when "p0" "q0" "degree p = degree q"
  proof -
    have "((λx. poly q x / poly p x)  lead_coeff q / lead_coeff p) at_top"
      using poly_divide_filterlim_at_top[OF that(1,2)] using that by auto
    then show ?thesis unfolding P_def 
      by (auto elim!:filterlim_at_bot_nhds filterlim_at_top_nhds)
  qed
  ultimately have P using False by fastforce
  then have "jumpF (λx. poly q x / poly p x) at_top = 0"
    unfolding jumpF_def P_def by auto
  then show ?thesis unfolding jumpF_poly_top_def using False by presburger
qed

lemma jumpF_poly_bot_code:
  "jumpF_poly_bot q p = (if p0  q0  degree q>degree p then 
          if sgn_neg_inf q * sgn_neg_inf p > 0 then 1/2 else -1/2 else 0)"
proof (cases "p0  q0  degree q>degree p")
  case True
  have ?thesis when "sgn_neg_inf q * sgn_neg_inf p > 0"
  proof -
    have "LIM x at_bot. poly q x / poly p x :> at_top"
      using poly_divide_filterlim_at_bot[of p q] True that by auto
    then have "jumpF (λx. poly q x / poly p x) at_bot = 1/2"
      unfolding jumpF_def by auto
    then show ?thesis unfolding jumpF_poly_bot_def using that True by auto
  qed
  moreover have ?thesis when "¬ sgn_neg_inf q * sgn_neg_inf p > 0"
  proof -
    have "LIM x at_bot. poly q x / poly p x :> at_bot"
      using poly_divide_filterlim_at_bot[of p q] True that by auto
    then have "jumpF (λx. poly q x / poly p x) at_bot = - 1/2"
      unfolding jumpF_def by auto
    then show ?thesis unfolding jumpF_poly_bot_def using that True by auto
  qed
  ultimately show ?thesis by auto
next
  case False
  define P where "P= (¬ (LIM x at_bot. poly q x / poly p x :> at_bot) 
                       ¬ (LIM x at_bot. poly q x / poly p x :> at_top))"
  have P when "p=0  q=0"
    unfolding P_def using that 
    by (auto elim!:filterlim_at_bot_nhds filterlim_at_top_nhds)
  moreover have P when "p0" "q0" "degree p> degree q"
  proof -
    have "LIM x at_bot. poly q x / poly p x :> at 0"
      using poly_divide_filterlim_at_bot[OF that(1,2)] that(3) by auto
    then show ?thesis unfolding P_def 
      by (auto elim!:filterlim_at_bot_nhds filterlim_at_top_nhds simp:filterlim_at)
  qed 
  moreover have P when "p0" "q0" "degree p = degree q"
  proof -
    have "((λx. poly q x / poly p x)  lead_coeff q / lead_coeff p) at_bot"
      using poly_divide_filterlim_at_bot[OF that(1,2)] using that by auto
    then show ?thesis unfolding P_def 
      by (auto elim!:filterlim_at_bot_nhds filterlim_at_top_nhds)
  qed
  ultimately have P using False by fastforce
  then have "jumpF (λx. poly q x / poly p x) at_bot = 0"
    unfolding jumpF_def P_def by auto
  then show ?thesis unfolding jumpF_poly_bot_def using False by presburger
qed
  
subsection ‹The extended Cauchy index on polynomials›

definition cindex_polyE:: "real  real  real poly  real poly  real" where
  "cindex_polyE a b q p = jumpF_polyR q p a + cindex_poly a b q p - jumpF_polyL q p b"
  
definition cindex_poly_ubd::"real poly  real poly  int" where
  "cindex_poly_ubd q p = (THE l. (F r in at_top. cindexE (-r) r (λx. poly q x/poly p x) = of_int l))"
   
lemma cindex_polyE_0[simp]: "cindex_polyE a b 0 p = 0" "cindex_polyE a b q 0 = 0"
  unfolding cindex_polyE_def by auto
  
lemma cindex_polyE_mult_cancel:
  fixes p q p'::"real poly"
  assumes "p' 0"  
  shows "cindex_polyE a b (p' * q ) (p' * p) = cindex_polyE a b q p"
  unfolding cindex_polyE_def
  using cindex_poly_mult[OF p'0] jumpF_polyL_mult_cancel[OF p'0] 
    jumpF_polyR_mult_cancel[OF p'0] 
  by simp
  
lemma cindexE_eq_cindex_polyE: 
  assumes "a<b"
  shows "cindexE a b (λx. poly q x/poly p x) = cindex_polyE a b q p"
proof (cases "p=0  q=0")
  case True
  then show ?thesis by (auto simp add: cindexE_constI)
next
  case False
  then have "p0" "q0" by auto
  define g where "g=gcd p q"
  define p' q' where "p'=p div g" and "q' = q div g"
  define f' where "f'=(λx. poly q' x / poly p' x)"
  have "g0" using False g_def by auto  
  have pq_f:"p=g*p'" "q=g*q'" and "coprime p' q'"
    unfolding g_def p'_def q'_def 
    apply simp_all 
    using False div_gcd_coprime by blast    
  have "cindexE a b (λx. poly q x/poly p x) = cindexE a b (λx. poly q' x/poly p' x)" 
  proof -
    define f where "f=(λx. poly q x / poly p x)"
    define f' where "f'=(λx. poly q' x / poly p' x)"
    have "jumpF f (at_right x) = jumpF f' (at_right x)" for x
    proof (rule jumpF_cong)
      obtain ub where "x < ub" "z. x < z  z  ub  poly g z  0"
        using next_non_root_interval[OF g0,of x] by auto
      then show "F x in at_right x. f x = f' x" 
        unfolding eventually_at_right f_def f'_def pq_f
        apply (intro exI[where x=ub])
        by auto
    qed simp
    moreover have "jumpF f (at_left x) = jumpF f' (at_left x)" for x 
    proof (rule jumpF_cong)
      obtain lb where "lb < x" "z. lb  z  z < x  poly g z  0 "
        using last_non_root_interval[OF g0,of x] by auto
      then show "F x in at_left x. f x = f' x" 
        unfolding eventually_at_left f_def f'_def pq_f
        apply (intro exI[where x=lb])
        by auto    
    qed simp    
    ultimately show ?thesis unfolding cindexE_def
      apply (fold f_def f'_def)
      by auto
  qed
  also have "... = jumpF f' (at_right a) + real_of_int (cindex a b f') - jumpF f' (at_left b)" 
    unfolding f'_def 
    apply (rule cindex_eq_cindexE_divide)
    subgoal using a<b .
    subgoal using False poly_roots_finite pq_f(1) pq_f(2) by fastforce
    subgoal using ‹coprime p' q' poly_gcd_iff by force
    subgoal by (auto intro:continuous_intros)
    subgoal by (auto intro:continuous_intros)
    done
  also have "... = cindex_polyE a b q' p'"
    using cindex_eq_cindex_poly unfolding cindex_polyE_def jumpF_polyR_def jumpF_polyL_def f'_def
    by auto
  also have "... = cindex_polyE a b q p"
    using cindex_polyE_mult_cancel[OF g0] unfolding pq_f by auto
  finally show ?thesis .
qed
 
lemma cindex_polyE_cross:
  fixes p::"real poly" and a b::real
  assumes "a<b" 
  shows "cindex_polyE a b 1 p = cross_alt 1 p a b / 2" 
proof (induct "degree p" arbitrary:p rule:nat_less_induct) 
  case induct:1
  have ?case when "p=0" 
    using that unfolding cross_alt_def by auto
  moreover have ?case when "p0" and noroot:"{x.  a< x x< b  poly p x=0 } = {}"
  proof -
    have "cindex_polyE a b 1 p = jumpF_polyR 1 p a  - jumpF_polyL 1 p b" 
    proof -
      have "cindex_poly a b 1 p = 0" unfolding cindex_poly_def
        apply (rule sum.neutral)
        using that by auto
      then show ?thesis unfolding cindex_polyE_def by auto
    qed
    also have "... = cross_alt 1 p a b / 2"  
    proof -
      define f where "f = (λx. 1 / poly p x)"
      define ja where "ja = jumpF f (at_right a)"  
      define jb where "jb = jumpF f (at_left b)"
      define right where "right = (λR. R ja (0::real)  (continuous (at_right a) f  R (poly p a) 0))"
      define left where "left = (λR. R jb (0::real)  (continuous (at_left b) f  R (poly p b) 0))" 
        
      note ja_alt=jumpF_polyR_coprime[of p 1 a,unfolded jumpF_polyR_def,simplified,folded f_def ja_def]
      note jb_alt=jumpF_polyL_coprime[of p 1 b,unfolded jumpF_polyL_def,simplified,folded f_def jb_def]  
      
      have [simp]:"0 < ja  jumpF_polyR 1 p a = 1/2" "0 > ja  jumpF_polyR 1 p a = -1/2"
           "0 < jb  jumpF_polyL 1 p b = 1/2" "0 > jb  jumpF_polyL 1 p b = -1/2"
        unfolding ja_def jb_def jumpF_polyR_def jumpF_polyL_def f_def jumpF_def
        by auto           
      have [simp]: 
          "poly p a 0  continuous (at_right a) f"
          "poly p b 0  continuous (at_left b) f"  
        unfolding f_def by (auto intro!: continuous_intros )  
      have not_right_left: False when "(right greater  left less  right less  left greater)"
      proof - 
        have [simp]: "f a > 0  poly p a > 0" "f a < 0  poly p a < 0"
             "f b > 0  poly p b > 0" "f b < 0  poly p b < 0" 
           unfolding f_def by auto  
        have "continuous_on {a<..<b} f" 
          unfolding f_def using noroot by (auto intro!: continuous_intros)
        then have "x>a. x < b  f x = 0" 
          apply (elim jumpF_IVT[OF a<b,of f])
          using that unfolding right_def left_def by (fold ja_def jb_def,auto)
        then show False using noroot using f_def by auto
      qed
      have ?thesis when "poly p a>0  poly p b>0  poly p a<0  poly p b<0" 
        using that jumpF_poly_noroot unfolding cross_alt_def by auto
      moreover have False when "poly p a>0  poly p b<0  poly p a<0  poly p b>0" 
        apply (rule not_right_left)
        unfolding right_def left_def using that by auto
      moreover have ?thesis when "poly p a=0" "poly p b>0  poly p b <0" 
      proof -
        have "ja>0  ja < 0" using ja_alt p0 ‹poly p a=0 by argo
        moreover have False when "ja > 0  poly p b<0  ja < 0  poly p b>0" 
          apply (rule not_right_left)
          unfolding right_def left_def using that by fastforce
        moreover have ?thesis when "ja >0  poly p b>0  ja < 0  poly p b<0"
          using that jumpF_poly_noroot ‹poly p a=0 unfolding cross_alt_def by auto 
        ultimately show ?thesis using that jumpF_poly_noroot unfolding cross_alt_def by auto 
      qed
      moreover have ?thesis when "poly p b=0" "poly p a>0  poly p a <0" 
      proof -
        have "jb>0  jb < 0" using jb_alt p0 ‹poly p b=0 by argo
        moreover have False when "jb > 0  poly p a<0  jb < 0  poly p a>0" 
          apply (rule not_right_left)
          unfolding right_def left_def using that by fastforce
        moreover have ?thesis when "jb >0  poly p a>0  jb < 0  poly p a<0"
          using that jumpF_poly_noroot ‹poly p b=0 unfolding cross_alt_def by auto 
        ultimately show ?thesis using that jumpF_poly_noroot unfolding cross_alt_def by auto 
      qed  
      moreover have ?thesis when "poly p a=0" "poly p b=0"
      proof -
        have "jb>0  jb < 0" using jb_alt p0 ‹poly p b=0 by argo
        moreover have "ja>0  ja < 0" using ja_alt p0 ‹poly p a=0 by argo
        moreover have False when "ja>0  jb<0  ja<0  jb>0"
          apply (rule not_right_left)
          unfolding right_def left_def using that by fastforce
        moreover have ?thesis when "ja>0  jb>0  ja<0  jb<0"
          using that jumpF_poly_noroot ‹poly p b=0 ‹poly p a=0 
          unfolding cross_alt_def by auto
        ultimately show ?thesis by blast
      qed
      ultimately show ?thesis by argo
    qed
    finally show ?thesis .
  qed    
  moreover have ?case when "p0" and no_empty:"{x.  a< x x< b  poly p x=0 }  {}"
  proof -
    define roots where "roots{x.  a< x x< b  poly p x=0 }"
    have "finite roots" unfolding roots_def using poly_roots_finite[OF p0] by auto
    define max_r where "max_rMax roots"
    hence "poly p max_r=0" and "a<max_r" and "max_r<b" 
      using Max_in[OF ‹finite roots] no_empty  unfolding roots_def by auto
    define max_rp where "max_rp[:-max_r,1:]^order max_r p"
    then obtain p' where p'_def:"p=p'*max_rp" and "¬ [:-max_r,1:] dvd p'"  
      by (metis p0 mult.commute order_decomp)
    hence "p'0" and "max_rp0" and max_r_nz:"poly p' max_r0"(*and "poly p' a≠0" and "poly p' b≠0" *)
      (*and  "poly max_rp a≠0" and "poly max_rp b≠0"*) 
      using p0 by (auto simp add: dvd_iff_poly_eq_0)
        
    define max_r_sign where "max_r_signif odd(order max_r p) then -1 else 1::int"
    define roots' where "roots'{x.  a< x x< b  poly p' x=0}"
  
    have "cindex_polyE a b 1 p = jumpF_polyR 1 p a + (xroots. jump_poly 1 p x) - jumpF_polyL 1 p b"  
      unfolding cindex_polyE_def cindex_poly_def roots_def by (simp,meson)
    also have "... = max_r_sign * cindex_poly a b 1 p' + jump_poly 1 p max_r 
        + max_r_sign * jumpF_polyR 1 p' a - jumpF_polyL 1 p' b"
    proof -
      have "(xroots. jump_poly 1 p x) = max_r_sign * cindex_poly a b 1 p' + jump_poly 1 p max_r"  
      proof -
        have "(xroots. jump_poly 1 p x)= (xroots'. jump_poly 1 p x)+ jump_poly 1 p max_r"
        proof -
          have "roots = insert max_r roots'" 
            unfolding roots_def roots'_def p'_def 
            using ‹poly p max_r=0 a<max_r max_r<b p0 order_root
            apply (subst max_rp_def)
            by auto
          moreover have "finite roots'" 
            unfolding roots'_def using poly_roots_finite[OF p'0] by auto 
          moreover have "max_r  roots'" 
            unfolding roots'_def using max_r_nz
            by auto
          ultimately show ?thesis using sum.insert[of roots' max_r] by auto   
        qed
        moreover have "(xroots'. jump_poly 1 p x) = max_r_sign * cindex_poly a b 1 p'"
        proof -
          have "(xroots'. jump_poly 1 p x) = (xroots'. max_r_sign * jump_poly 1 p' x)"
          proof (rule sum.cong,rule refl)
            fix x assume "x  roots'" 
            hence "xmax_r" using max_r_nz unfolding roots'_def 
              by auto
            hence "poly max_rp x0" using poly_power_n_eq unfolding max_rp_def by auto
            hence "order x max_rp=0"  by (metis order_root)
            moreover have "jump_poly 1 max_rp x=0" 
              using ‹poly max_rp x0 by (metis jump_poly_not_root)
            moreover have "xroots"
              using x  roots' unfolding roots_def roots'_def p'_def by auto
            hence "x<max_r" 
              using Max_ge[OF ‹finite roots,of x] xmax_r by (fold max_r_def,auto)
            hence "sign (poly max_rp x) = max_r_sign" 
              using ‹poly max_rp x  0 unfolding max_r_sign_def max_rp_def sign_def
              by (subst poly_power,simp add:linorder_class.not_less zero_less_power_eq)
            ultimately show "jump_poly 1 p x = max_r_sign * jump_poly 1 p' x" 
              using jump_poly_1_mult[of p' x max_rp]  unfolding p'_def 
              by (simp add: ‹poly max_rp x  0)
          qed
          also have "... = max_r_sign * (xroots'. jump_poly 1 p' x)" 
            by (simp add: sum_distrib_left) 
          also have "... = max_r_sign * cindex_poly a b 1 p'"
            unfolding cindex_poly_def roots'_def by meson
          finally show ?thesis .
        qed
        ultimately show ?thesis by simp
      qed
      moreover have "jumpF_polyR 1 p a = max_r_sign * jumpF_polyR 1 p' a"
      proof -
        define f where "f = (λx. 1 / poly max_rp x)"
        define g where "g = (λx. 1 / poly p' x)"
        have "jumpF_polyR 1 p a = jumpF (λx. f x * g x) (at_right a)"  
          unfolding jumpF_polyR_def f_def g_def p'_def 
          by (auto simp add:field_simps)
        also have "... = sgn (f a) * jumpF g (at_right a)" 
        proof (rule jumpF_times) 
          have [simp]: "poly max_rp a 0" 
            unfolding max_rp_def using max_r>a by auto  
          show "(f  f a) (at_right a)" "f a  0"
            unfolding f_def by (auto intro:tendsto_intros)
        qed auto      
        also have "... = max_r_sign * jumpF_polyR 1 p' a"
        proof -
          have "sgn (f a) = max_r_sign" 
            unfolding max_r_sign_def f_def max_rp_def using a<max_r
            by (auto simp add:sgn_power)
          then show ?thesis unfolding jumpF_polyR_def g_def by auto
        qed
        finally show ?thesis .
      qed
      moreover have "jumpF_polyL 1 p b =  jumpF_polyL 1 p' b"
      proof -
        define f where "f = (λx. 1 / poly max_rp x)"
        define g where "g = (λx. 1 / poly p' x)"
        have "jumpF_polyL 1 p b = jumpF (λx. f x * g x) (at_left b)"  
          unfolding jumpF_polyL_def f_def g_def p'_def 
          by (auto simp add:field_simps)
        also have "... = sgn (f b) * jumpF g (at_left b)" 
        proof (rule jumpF_times) 
          have [simp]: "poly max_rp b 0" 
            unfolding max_rp_def using max_r<b by auto  
          show "(f  f b) (at_left b)" "f b  0"
            unfolding f_def by (auto intro:tendsto_intros)
        qed auto      
        also have "... = jumpF_polyL 1 p' b"
        proof -
          have "sgn (f b) = 1" 
            unfolding max_r_sign_def f_def max_rp_def using b>max_r
            by (auto simp add:sgn_power)
          then show ?thesis unfolding jumpF_polyL_def g_def by auto
        qed
        finally show ?thesis .
      qed 
      ultimately show ?thesis by auto
    qed
    also have "... = max_r_sign * cindex_polyE a b 1 p' + jump_poly 1 p max_r 
        + (max_r_sign - 1) * jumpF_polyL 1 p' b"
      unfolding cindex_polyE_def roots'_def by (auto simp add:algebra_simps)
    also have "... = max_r_sign * cross_alt 1 p' a b / 2 + jump_poly 1 p max_r 
        + (max_r_sign - 1) * jumpF_polyL 1 p' b"
    proof -
      have "degree max_rp>0" unfolding max_rp_def degree_linear_power
        using ‹poly p max_r=0 order_root p0 by blast
      then have "degree p'<degree p" unfolding p'_def 
        using degree_mult_eq[OF p'0 max_rp0] by auto
      from induct[rule_format, OF this] 
      have "cindex_polyE a b 1 p' = real_of_int (cross_alt 1 p' a b) / 2" by auto
      then show ?thesis by auto
    qed
    also have "... = real_of_int (cross_alt 1 p a b) / 2"
    proof -
      have sjump_p:"jump_poly 1 p max_r = (if odd (order max_r p) then sign (poly p' max_r) else 0)"
      proof -
        note max_r_nz 
        moreover then have "poly max_rp max_r=0" 
          using ‹poly p max_r = 0 p'_def by auto
        ultimately have "jump_poly 1 p max_r = sign (poly p' max_r) * jump_poly 1 max_rp max_r"
          unfolding p'_def using jump_poly_1_mult[of p' max_r max_rp] 
          by auto
        also have "... = (if odd (order max_r max_rp) then sign (poly p' max_r) else 0)"  
        proof -
          have "sign_r_pos max_rp max_r"
            unfolding max_rp_def using sign_r_pos_power by auto
          then show ?thesis using max_rp0 unfolding jump_poly_def by auto
        qed
        also have "... = (if odd (order max_r p) then sign (poly p' max_r) else 0)"
        proof -
          have "order max_r p'=0" by (simp add: ‹poly p' max_r  0 order_0I)
          then have "order max_r max_rp = order max_r p" 
            unfolding p'_def using p'0 max_rp0
            apply (subst order_mult)
            by auto  
          then show ?thesis by auto
        qed
        finally show ?thesis .
      qed
      have ?thesis when "even (order max_r p)"
      proof -
        have "sign (poly p x) =  sign (poly p' x)" when "xmax_r" for x
        proof -
          have "sign (poly max_rp x) = 1"
            unfolding max_rp_def using ‹even (order max_r p) that 
            apply (simp add:sign_power )
            by (simp add: Sturm_Tarski.sign_def)
          then show ?thesis unfolding p'_def by (simp add:sign_times)
        qed    
        from this[of a] this[of b] a<max_r max_r<b 
        have "cross_alt 1 p' a b = cross_alt 1 p a b" 
          unfolding cross_alt_def by auto 
        then show ?thesis using that unfolding max_r_sign_def sjump_p by auto
      qed
      moreover have ?thesis when "odd (order max_r p)" 
      proof -
        let ?thesis2 = "sign (poly p' max_r) * 2 - cross_alt 1 p' a b - 4 * jumpF_polyL 1 p' b 
              = cross_alt 1 p a b"    
        have ?thesis2 when "poly p' b=0"
        proof -
          have "jumpF_polyL 1 p' b = 1/2  jumpF_polyL 1 p' b=-1/2"  
            using jumpF_polyL_coprime[of p' 1 b,simplified] p'0 ‹poly p' b=0 by auto
          moreover have "poly p' max_r>0  poly p' max_r<0" 
            using max_r_nz by auto
          moreover have False when "poly p' max_r>0  jumpF_polyL 1 p' b=-1/2 
                 poly p' max_r<0  jumpF_polyL 1 p' b=1/2"
          proof -
            define f where "f= (λx. 1/ poly p' x)"
            have noroots:"poly p' x0" when "x{max_r<..<b}" for x
            proof (rule ccontr)
              assume " ¬ poly p' x  0"
              then have "poly p x =0" unfolding p'_def by auto
              then have "xroots" unfolding roots_def using that a<max_r by auto
              then have "xmax_r" using Max_ge[OF ‹finite roots] unfolding max_r_def by auto
              moreover have "x>max_r" using that by auto
              ultimately show False by auto
            qed  
            have "continuous_on {max_r<..<b} f"
              unfolding f_def using noroots by (auto intro!:continuous_intros)
            moreover have "continuous (at_right max_r) f" 
              unfolding f_def using max_r_nz
              by (auto intro!:continuous_intros)
            moreover have "f max_r>0  jumpF f (at_left b)<0 
                 f max_r<0  jumpF f (at_left b)>0"
              using that unfolding f_def jumpF_polyL_def by auto  
            ultimately have "x>max_r. x < b  f x = 0"
              apply (intro jumpF_IVT[OF max_r<b])
              by auto
            then show False using noroots unfolding f_def by auto
          qed
          moreover have ?thesis when "poly p' max_r>0  jumpF_polyL 1 p' b=1/2
               poly p' max_r<0  jumpF_polyL 1 p' b=-1/2"
          proof -
            have "poly max_rp a < 0" "poly max_rp b>0"
              unfolding max_rp_def using ‹odd (order max_r p) a<max_r max_r<b
              by (simp_all add:zero_less_power_eq)
            then have "cross_alt 1 p a b = - cross_alt 1 p' a b" 
              unfolding cross_alt_def p'_def using ‹poly p' b=0 
              apply (simp add:sign_times)
              by (simp add: Sturm_Tarski.sign_def)
            with that show ?thesis by auto
          qed
          ultimately show ?thesis by blast  
        qed
        moreover have ?thesis2 when "poly p' b0"
        proof -
          have [simp]:"jumpF_polyL 1 p' b = 0" 
            using jumpF_polyL_coprime[of p' 1 b,simplified] ‹poly p' b0 by auto 
          have [simp]:"poly max_rp a < 0" "poly max_rp b>0"
            unfolding max_rp_def using ‹odd (order max_r p) a<max_r max_r<b
            by (simp_all add:zero_less_power_eq)
          have "poly p' b>0  poly p' b<0" 
            using ‹poly p' b0 by auto
          moreover have "poly p' max_r>0  poly p' max_r<0" 
            using max_r_nz by auto
          moreover have ?thesis when "poly p' b>0  poly p' max_r>0 "  
            using that unfolding cross_alt_def p'_def
            apply (simp add:sign_times)
            by (simp add: Sturm_Tarski.sign_def)
          moreover have ?thesis when "poly p' b<0  poly p' max_r<0"      
            using that unfolding cross_alt_def p'_def
            apply (simp add:sign_times)
            by (simp add: Sturm_Tarski.sign_def)  
          moreover have False when "poly p' b>0  poly p' max_r<0  poly p' b<0  poly p' max_r>0"
          proof -
            have "x>max_r. x < b  poly p' x = 0"
              apply (rule poly_IVT[OF max_r<b,of p'])
              using that mult_less_0_iff by blast
            then obtain x where "max_r<x" "x<b" "poly p x=0" unfolding p'_def by auto
            then have "xroots" using a<max_r unfolding roots_def by auto
            then have "xmax_r" unfolding max_r_def using Max_ge[OF ‹finite roots] by auto    
            then show False using max_r<x by auto
          qed
          ultimately show ?thesis by blast
        qed
        ultimately have ?thesis2 by auto 
        then show ?thesis unfolding max_r_sign_def sjump_p using that by simp
      qed
      ultimately show ?thesis by auto
    qed
    finally show ?thesis . 
  qed
  ultimately show ?case by fast
qed      
     
lemma cindex_polyE_inverse_add:
  fixes p q::"real poly" 
  assumes cp:"coprime p q"
  shows "cindex_polyE a b q p + cindex_polyE a b p q=cindex_polyE a b 1 (q*p)"
  unfolding cindex_polyE_def
  using cindex_poly_inverse_add[OF cp,symmetric] jumpF_polyR_inverse_add[OF cp,symmetric] 
    jumpF_polyL_inverse_add[OF cp,symmetric]
  by auto    

lemma cindex_polyE_inverse_add_cross:
  fixes p q::"real poly"
  assumes "a < b" "coprime p q" 
  shows "cindex_polyE a b q p  + cindex_polyE a b p q = cross_alt p q a b / 2" 
  apply (subst cindex_polyE_inverse_add[OF ‹coprime p q])
  apply (subst cindex_polyE_cross[OF a<b])
  apply (subst mult.commute)  
  apply (subst cross_alt_clear_n[OF ‹coprime p q])
  by simp
      
lemma cindex_polyE_smult_1: 
  fixes p q::"real poly" and c::real
  shows "cindex_polyE a b (smult c q) p =  (sgn c) * cindex_polyE a b q p"
  unfolding cindex_polyE_def jumpF_polyL_smult_1 jumpF_polyR_smult_1 cindex_poly_smult_1 
  by (auto simp add:sgn_sign_eq[symmetric] algebra_simps)    

lemma cindex_polyE_mod:
  fixes p q::"real poly" 
  shows "cindex_polyE a b q p =  cindex_polyE a b (q mod p) p"
  unfolding cindex_polyE_def
  apply (subst cindex_poly_mod)
  apply (subst jumpF_polyR_mod)
  apply (subst jumpF_polyL_mod)
  by simp

lemma cindex_polyE_rec:
  fixes p q::"real poly"
  assumes "a < b" "coprime p q"
  shows "cindex_polyE a b q p  = cross_alt q p a b/2  +  cindex_polyE a b (- (p mod q)) q"  
proof -
  note cindex_polyE_inverse_add_cross[OF assms]
  moreover have "cindex_polyE a b (- (p mod q)) q = - cindex_polyE a b p q"
    using cindex_polyE_mod cindex_polyE_smult_1[of a b "-1"]
    by auto
  ultimately show ?thesis by (auto simp add:field_simps cross_alt_poly_commute)
qed    
  
lemma cindex_polyE_changes_alt_itv_mods: 
  assumes "a<b" "coprime p q"
  shows "cindex_polyE a b q p = changes_alt_itv_smods a b p q / 2" using ‹coprime p q
proof (induct "smods p q" arbitrary:p q)
  case Nil
  then have "p=0" by (metis smods_nil_eq)
  then show ?case by (simp add:changes_alt_itv_smods_def changes_alt_poly_at_def) 
next
  case (Cons x xs)
  then have "p0" by auto
  have ?case when "q=0" 
    using that by (simp add:changes_alt_itv_smods_def changes_alt_poly_at_def)
  moreover have ?case when "q0"
  proof -
    define r where "r- (p mod q)"
    obtain ps where ps:"smods p q=p#q#ps" "smods q r=q#ps" and "xs=q#ps"
      unfolding r_def using q0 p0 x # xs = smods p q 
      by (metis list.inject smods.simps)
    from Cons.prems q  0 have "coprime q r" 
      by (simp add: r_def ac_simps)
    then have "cindex_polyE a b r q = real_of_int (changes_alt_itv_smods a b q r) / 2"
      apply (rule_tac Cons.hyps(1))
      using ps xs=q#ps by simp_all  
    moreover have "changes_alt_itv_smods a b p q = cross_alt p q a b + changes_alt_itv_smods a b q r" 
      using changes_alt_itv_smods_rec[OF a<b ‹coprime p q,folded r_def] .
    moreover have "cindex_polyE a b q p = real_of_int (cross_alt q p a b) / 2 + cindex_polyE a b r q"
      using cindex_polyE_rec[OF a<b ‹coprime p q,folded r_def] .
    ultimately show ?case 
      by (auto simp add:field_simps cross_alt_poly_commute)
  qed
  ultimately show ?case by blast
qed
  
lemma cindex_poly_ubd_eventually:
  shows "F r in at_top. cindexE (-r) r (λx. poly q x/poly p x) = of_int (cindex_poly_ubd q p)" 
proof -
  define f where "f=(λx. poly q x/poly p x)"
  obtain R where R_def: "R>0" "proots p  {-R<..<R}"
    if "p0" 
  proof (cases "p=0")
    case True
    then show ?thesis using that[of 1] by auto
  next
    case False
    then have "finite (proots p)" by auto
    from finite_ball_include[OF this,of 0] 
    obtain r where "r>0" and r_ball:"proots p  ball 0 r"
      by auto
    have "proots p  {-r<..<r}"
    proof 
      fix x assume "x  proots p"
      then have "xball 0 r" using r_ball by auto
      then have "abs x<r" using mem_ball_0 by auto
      then show "x  {- r<..<r}" using r>0 by auto
    qed
    then show ?thesis using that[of r] False r>0 by auto
  qed
  define l where "l=(if p=0  then 0 else cindex_poly (-R) R q p)"  
  define P where "P=(λl. (F r in at_top. cindexE (-r) r f = of_int l))"
  have "P l " 
  proof (cases "p=0 ")
    case True
    then show ?thesis
      unfolding P_def f_def l_def using True
      by (auto intro!: eventuallyI cindexE_constI)
  next
    case False
    have "P l" unfolding P_def
    proof (rule eventually_at_top_linorderI[of R])  
      fix r assume "R  r" 
      then have "cindexE (- r) r f =  cindex_polyE (-r) r q p"
        unfolding f_def using R_def[OF p0] by (auto intro: cindexE_eq_cindex_polyE)
      also have "... = of_int (cindex_poly (- r) r q p)"
      proof -
        have "jumpF_polyR q p (- r) = 0"
          apply (rule jumpF_poly_noroot)
          using Rr R_def[OF p0] by auto
        moreover have "jumpF_polyL q p r = 0"
          apply (rule jumpF_poly_noroot)
          using Rr R_def[OF p0] by auto
        ultimately show ?thesis unfolding cindex_polyE_def by auto
      qed  
      also have "... = of_int (cindex_poly (- R) R q p)"
      proof -
        define rs where "rs={x. poly p x = 0  - r < x  x < r}"
        define Rs where "Rs={x. poly p x = 0  - R < x  x < R}"
        have "rs=Rs" 
          using R_def[OF p0] Rr unfolding rs_def Rs_def by force    
        then show ?thesis
          unfolding cindex_poly_def by (fold rs_def Rs_def,auto)
      qed
      also have "... = of_int l" unfolding l_def using False by auto
      finally show "cindexE (- r) r f = real_of_int l" .
    qed
    then show ?thesis unfolding P_def by auto
  qed
  moreover have "x=l" when "P x" for x 
  proof -
    have "F r in at_top. cindexE (- r) r f = real_of_int x"
         "F r in at_top. cindexE (- r) r f = real_of_int l"
      using P x P l unfolding P_def by auto
    from eventually_conj[OF this] 
    have "F r::real in at_top. real_of_int x = real_of_int l"
      by (elim eventually_mono,auto)
    then have "real_of_int x = real_of_int l" by auto
    then show ?thesis by simp
  qed
  ultimately have "P (THE x. P x)" using theI[of P l] by blast
  then show ?thesis unfolding P_def f_def cindex_poly_ubd_def by auto 
qed
  
lemma cindex_poly_ubd_0:
  assumes "p=0  q=0"
  shows "cindex_poly_ubd q p = 0"  
proof -
  have "F r in at_top. cindexE (-r) r (λx. poly q x/poly p x) = 0"
    apply (rule eventuallyI)
    using assms by (auto intro:cindexE_constI)
  from eventually_conj[OF this cindex_poly_ubd_eventually[of q p]]
  have "F r::real in at_top.  (cindex_poly_ubd q p) = (0::int)"
    apply (elim eventually_mono)
    by auto
  then show ?thesis by auto
qed
  
lemma cindex_poly_ubd_code:
  shows "cindex_poly_ubd q p = changes_R_smods p q"
proof (cases "p=0")
  case True
  then show ?thesis using cindex_poly_ubd_0 by auto
next
  case False
  define ps where "pssmods p q"
  have "pset ps" using ps_def p0 by auto
  obtain lb where lb:"pset ps. x. poly p x=0  x>lb"
      and lb_sgn:"xlb. pset ps. sgn (poly p x) = sgn_neg_inf p"
      and "lb<0"
    using root_list_lb[OF no_0_in_smods,of p q,folded ps_def] 
    by auto
  obtain ub where ub:"pset ps. x. poly p x=0  x<ub"
      and ub_sgn:"xub. pset ps. sgn (poly p x) = sgn_pos_inf p"
      and "ub>0"
    using root_list_ub[OF no_0_in_smods,of p q,folded ps_def] 
    by auto
  define f where "f=(λt. poly q t/ poly p t)"
  define P where "P=(λl. (F r in at_top. cindexE (-r) r f = of_int l))"
  have "P (changes_R_smods p q)" unfolding P_def
  proof (rule eventually_at_top_linorderI[of "max ¦lb¦ ¦ub¦ + 1"])
    fix r assume r_asm:"rmax ¦lb¦ ¦ub¦ + 1"
    have "cindexE (- r) r f =  cindex_polyE (-r) r q p"
      unfolding f_def using r_asm by (auto intro: cindexE_eq_cindex_polyE)
    also have "... = of_int (cindex_poly (- r) r q p)"
    proof -
      have "jumpF_polyR q p (- r) = 0"
        apply (rule jumpF_poly_noroot)
        using r_asm lb[rule_format,OF pset ps,of "-r"] by linarith
      moreover have "jumpF_polyL q p r = 0"
        apply (rule jumpF_poly_noroot)
        using r_asm ub[rule_format,OF pset ps,of "r"] by linarith
      ultimately show ?thesis unfolding cindex_polyE_def by auto
    qed
    also have "... = of_int (changes_itv_smods (- r) r p q)"
      apply (rule cindex_poly_changes_itv_mods[THEN arg_cong])
      using r_asm lb[rule_format,OF pset ps,of "-r"] ub[rule_format,OF pset ps,of "r"]
      by linarith+
    also have "... = of_int (changes_R_smods p q)"
    proof -
      have "map (sgn  (λp. poly p (-r))) ps = map sgn_neg_inf ps"
          and "map (sgn  (λp. poly p r)) ps = map sgn_pos_inf ps"
        using lb_sgn[THEN spec,of "-r",simplified] ub_sgn[THEN spec,of r,simplified] r_asm 
        by auto
      hence "changes_poly_at ps (-r)=changes_poly_neg_inf ps
           changes_poly_at ps r=changes_poly_pos_inf ps"
        unfolding changes_poly_neg_inf_def changes_poly_at_def changes_poly_pos_inf_def
        by (subst (1 3)  changes_map_sgn_eq,metis map_map)
      thus ?thesis unfolding changes_R_smods_def changes_itv_smods_def ps_def
        by metis
    qed
    finally show "cindexE (- r) r f =  of_int (changes_R_smods p q)" .
  qed
  moreover have "x = changes_R_smods p q" when "P x" for x 
  proof -
    have "F r in at_top. cindexE (- r) r f = real_of_int (changes_R_smods p q)" 
        "F r in at_top. cindexE (- r) r f = real_of_int x"
      using P (changes_R_smods p q) P x unfolding P_def by auto
    from eventually_conj[OF this] 
    have "F (r::real) in at_top. of_int x = of_int (changes_R_smods p q)"
      by (elim eventually_mono,auto)
    then have "of_int x = of_int (changes_R_smods p q)" 
      using eventually_const_iff by auto
    then show ?thesis using of_int_eq_iff by blast
  qed
  ultimately have "(THE x. P x) = changes_R_smods p q" 
    using the_equality[of P "changes_R_smods p q"] by blast
  then show ?thesis unfolding cindex_poly_ubd_def P_def f_def by auto
qed  


lemma cindexE_ubd_poly: "cindexE_ubd (λx. poly q x/poly p x) = cindex_poly_ubd q p"
proof (cases "p=0")
  case True
  then show ?thesis using cindex_poly_ubd_0 unfolding cindexE_ubd_def 
    by auto
next
  case False
  define mx mn where "mx = Max {x. poly p x = 0}" and "mn = Min {x. poly p x=0}"
  define rr where "rr= 1+ (max ¦mx¦ ¦mn¦)"
  have rr:"-rr < x  x< rr" when "poly p x = 0 " for x
  proof -
    have "finite {x. poly p x = 0}" using p0 poly_roots_finite by blast
    then have "mn  x" "xmx"
      using Max_ge Min_le that unfolding mn_def mx_def by simp_all
    then show ?thesis unfolding rr_def by auto
  qed
  define f where "f=(λx. poly q x / poly p x)"
  have "F r in at_top. cindexE (- r) r f = cindexE_ubd f"
  proof (rule eventually_at_top_linorderI[of rr])
    fix r assume "rrr"
    define R1 R2 where "R1={x. jumpF f (at_right x)  0  - r  x  x < r}"
                       and "R2 = {x. jumpF f (at_right x)  0}"
    define L1 L2 where "L1={x. jumpF f (at_left x)  0  - r < x  x  r}"
                       and "L2={x. jumpF f (at_left x)  0}"
    have "R1=R2"
    proof -
      have "jumpF f (at_right x) = 0" when "¬ (- r  x  x < r)" for x 
      proof -
        have "jumpF f (at_right x) = jumpF_polyR q p x"
          unfolding f_def jumpF_polyR_def by simp
        also have "... = 0"
          apply (rule jumpF_poly_noroot)
          using  that rrr by (auto dest:rr)
        finally show ?thesis .
      qed
      then show ?thesis unfolding R1_def R2_def by blast
    qed
    moreover have "L1=L2"
    proof -
      have "jumpF f (at_left x) = 0" when "¬ (- r < x  x  r)" for x 
      proof -
        have "jumpF f (at_left x) = jumpF_polyL q p x"
          unfolding f_def jumpF_polyL_def by simp
        also have "... = 0"
          apply (rule jumpF_poly_noroot)
          using that rrr by (auto dest:rr)
        finally show ?thesis .
      qed
      then show ?thesis unfolding L1_def L2_def by blast
    qed
    ultimately show "cindexE (- r) r f = cindexE_ubd f" 
      unfolding cindexE_def cindexE_ubd_def
      apply (fold R1_def R2_def L1_def L2_def)
      by auto
  qed
  moreover have "F r in at_top. cindexE (- r) r f = cindex_poly_ubd q p"
    using cindex_poly_ubd_eventually unfolding f_def by auto
  ultimately have "F r in at_top. cindexE (- r) r f = cindexE_ubd f 
                           cindexE (- r) r f = cindex_poly_ubd q p"
    using eventually_conj by auto
  then have "F (r::real) in at_top. cindexE_ubd f = cindex_poly_ubd q p"
    by (elim eventually_mono) auto
  then show ?thesis unfolding f_def by auto
qed
  
end

Theory More_Polynomials

(*
    Author:     Wenda Li <wl302@cam.ac.uk / liwenda1990@hotmail.com>
*)
section ‹More useful lemmas related polynomials›

theory More_Polynomials imports 
  Winding_Number_Eval.Missing_Algebraic
  Winding_Number_Eval.Missing_Transcendental
  Sturm_Tarski.PolyMisc
  Budan_Fourier.BF_Misc
begin

subsection ‹More about @{term order}

lemma order_normalize[simp]:"order x (normalize p) = order x p"      
  by (metis dvd_normalize_iff normalize_eq_0_iff order_1 order_2 order_unique_lemma)

lemma order_gcd:
  assumes "p0" "q0"
  shows "order x (gcd p q) = min (order x p) (order x q)"
proof -
  define xx op oq where "xx=[:- x, 1:]" and "op = order x p" and "oq = order x q"
  obtain pp where pp:"p = xx ^ op * pp" "¬ xx dvd pp"
    using order_decomp[OF p0,of x,folded xx_def op_def] by auto
  obtain qq where qq:"q = xx ^ oq * qq" "¬ xx dvd qq"
    using order_decomp[OF q0,of x,folded xx_def oq_def] by auto
  define pq where "pq = gcd pp qq"

  have p_unfold:"p = (pq * xx ^ (min op oq)) * ((pp div pq) * xx ^ (op - min op oq))"
        and [simp]:"coprime xx (pp div pq)" and "pp0"
  proof -
    have "xx ^ op = xx ^ (min op oq) * xx ^ (op - min op oq)" 
      by (simp flip:power_add)
    moreover have "pp = pq * (pp div pq)" 
      unfolding pq_def by simp
    ultimately show "p = (pq * xx ^ (min op oq)) * ((pp div pq) * xx ^ (op - min op oq))"
      unfolding pq_def pp by(auto simp:algebra_simps)
    show "coprime xx (pp div pq)" 
      apply (rule prime_elem_imp_coprime[OF 
                    prime_elem_linear_poly[of 1 "-x",simplified],folded xx_def])
      using pp = pq * (pp div pq) pp(2) by auto
  qed (use pp p0 in auto)
  have q_unfold:"q = (pq * xx ^ (min op oq)) * ((qq div pq) * xx ^ (oq - min op oq))"
         and [simp]:"coprime xx (qq div pq)" 
  proof -
    have "xx ^ oq = xx ^ (min op oq) * xx ^ (oq - min op oq)" 
      by (simp flip:power_add)
    moreover have "qq = pq * (qq div pq)" 
      unfolding pq_def by simp
    ultimately show "q = (pq * xx ^ (min op oq)) * ((qq div pq) * xx ^ (oq - min op oq))"
      unfolding pq_def qq by(auto simp:algebra_simps)
    show "coprime xx (qq div pq)" 
      apply (rule prime_elem_imp_coprime[OF 
                    prime_elem_linear_poly[of 1 "-x",simplified],folded xx_def])
      using qq = pq * (qq div pq) qq(2) by auto
  qed

  have "gcd p q=normalize (pq * xx ^ (min op oq))"
  proof -
    have "coprime (pp div pq * xx ^ (op - min op oq)) (qq div pq * xx ^ (oq - min op oq))"
    proof (cases "op>oq")
      case True
      then have "oq - min op oq = 0" by auto
      moreover have "coprime (xx ^ (op - min op oq)) (qq div pq)" by auto
      moreover have "coprime (pp div pq) (qq div pq)"
        apply (rule div_gcd_coprime[of pp qq,folded pq_def])
        using pp0 by auto
      ultimately show ?thesis by auto
    next
      case False
      then have "op - min op oq = 0" by auto
      moreover have "coprime (pp div pq) (xx ^ (oq - min op oq))" 
        by (auto simp:coprime_commute)
      moreover have "coprime (pp div pq) (qq div pq)"
        apply (rule div_gcd_coprime[of pp qq,folded pq_def])
        using pp0 by auto
      ultimately show ?thesis by auto
    qed 
    then show ?thesis unfolding p_unfold q_unfold
      apply (subst gcd_mult_left)
      by auto
  qed
  then have "order x (gcd p q) = order x pq + order x (xx ^ (min op oq))"
    apply simp
    apply (subst order_mult)
    using assms(1) p_unfold by auto
  also have "... = order x (xx ^ (min op oq))"
    using pp(2) qq(2) unfolding pq_def xx_def 
    by (auto simp add: order_0I poly_eq_0_iff_dvd)
  also have "... = min op oq"
    unfolding xx_def by (rule order_power_n_n)
  also have "... = min (order x p) (order x q)" unfolding op_def oq_def by simp
  finally show ?thesis .
qed

lemma pderiv_power: "pderiv (p ^ n) = smult (of_nat n) (p ^ (n-1)) * pderiv p"
  apply (cases n)
  using pderiv_power_Suc by auto

(*TODO: to replace the one (with the same name) in the library, as this version does 
  not require 'a to be a field?*)
lemma order_pderiv:
  fixes p::"'a::{idom,semiring_char_0} poly"
  assumes "p0" "poly p x=0"
  shows "order x p = Suc (order x (pderiv p))" using assms
proof -
  define xx op where "xx=[:- x, 1:]" and "op = order x p"
  have "op 0" unfolding op_def using assms order_root by blast
  obtain pp where pp:"p = xx ^ op * pp" "¬ xx dvd pp"
    using order_decomp[OF p0,of x,folded xx_def op_def] by auto
  have p_der:"pderiv p = smult (of_nat op) (xx^(op -1)) * pp + xx^op*pderiv pp"
    unfolding pp(1) by (auto simp:pderiv_mult pderiv_power xx_def algebra_simps pderiv_pCons)
  have "xx^(op -1) dvd (pderiv p)"
    unfolding p_der 
    by (metis One_nat_def Suc_pred assms(1) assms(2) dvd_add dvd_mult_right dvd_triv_left 
        neq0_conv op_def order_root power_Suc smult_dvd_cancel)
  moreover have "¬ xx^op dvd (pderiv p)"
  proof 
    assume "xx ^ op dvd pderiv p"
    then have "xx ^ op dvd smult (of_nat op) (xx^(op -1) * pp)"
      unfolding p_der by (simp add: dvd_add_left_iff)
    then have "xx ^ op dvd (xx^(op -1)) * pp"
      apply (elim dvd_monic[rotated])
      using op0 by (auto simp:lead_coeff_power xx_def)
    then have "xx ^ (op-1) * xx dvd (xx^(op -1))"
      using ¬ xx dvd pp by (simp add: op  0 mult.commute power_eq_if)
    then have "xx dvd 1" 
      using assms(1) pp(1) by auto
    then show False unfolding xx_def by (meson assms(1) dvd_trans one_dvd order_decomp)
  qed
  ultimately have "op - 1 = order x (pderiv p)"
    using order_unique_lemma[of x "op-1" "pderiv p",folded xx_def] op0 
    by auto
  then show ?thesis using op0 unfolding op_def by auto
qed

subsection ‹More about @{term rsquarefree}

lemma rsquarefree_0[simp]: "¬ rsquarefree 0"
  unfolding rsquarefree_def by simp

lemma rsquarefree_times:
  assumes "rsquarefree (p*q)"
  shows "rsquarefree q" using assms
proof (induct p rule:poly_root_induct_alt)
  case 0
  then show ?case by simp
next
  case (no_proots p)
  then have [simp]:"p0" "q0" "a. order a p = 0" 
    using order_0I by auto
  have "order a (p * q) = 0  order a q = 0"
       "order a (p * q) = 1  order a q = 1"
       for a
    subgoal by (subst order_mult) auto
    subgoal by (subst order_mult) auto
    done
  then show ?case using ‹rsquarefree (p * q)
    unfolding rsquarefree_def by simp
next
  case (root a p)
  define pq aa where "pq = p * q" and "aa = [:- a, 1:]"
  have [simp]:"pq0" "aa0" "order a aa=1"
    subgoal using pq_def root.prems by auto
    subgoal by (simp add: aa_def)
    subgoal by (metis aa_def order_power_n_n power_one_right)
    done
  have "rsquarefree (aa * pq)"
    unfolding aa_def pq_def using root(2) by (simp add:algebra_simps)
  then have "rsquarefree pq"
    unfolding rsquarefree_def by (auto simp add:order_mult)
  from root(1)[OF this[unfolded pq_def]] show ?case .
qed

lemma rsquarefree_smult_iff:
  assumes "s0"
  shows "rsquarefree (smult s p)  rsquarefree p"
  unfolding rsquarefree_def using assms by (auto simp add:order_smult)

lemma card_proots_within_rsquarefree:
  assumes "rsquarefree p"
  shows "proots_count p s = card (proots_within p s)" using assms
proof (induct rule:poly_root_induct[of _ "λx. xs"])
  case 0
  then have False by simp
  then show ?case by simp
next
  case (no_roots p)
  then show ?case 
    by (metis all_not_in_conv card.empty proots_count_def proots_within_iff sum.empty)
next
  case (root a p)
  have "proots_count ([:a, - 1:] * p) s = 1 + proots_count p s"
    apply (subst proots_count_times)
    subgoal using root.prems rsquarefree_def by blast
    subgoal by (metis (no_types, hide_lams) add.inverse_inverse add.inverse_neutral 
                  minus_pCons proots_count_pCons_1_iff proots_count_uminus root.hyps(1))  
    done
  also have "... = 1 + card (proots_within p s)"
  proof -
    have "rsquarefree p" using ‹rsquarefree ([:a, - 1:] * p) 
      by (elim rsquarefree_times)
    from root(2)[OF this] show ?thesis by simp
  qed
  also have "... = card (proots_within ([:a, - 1:] * p) s)" unfolding proots_within_times 
  proof (subst card_Un_disjoint)
    have [simp]:"p0" using root.prems by auto
    show "finite (proots_within [:a, - 1:] s)" "finite (proots_within p s)"
      by auto
    show " 1 + card (proots_within p s) = card (proots_within [:a, - 1:] s)
               + card (proots_within p s)"
      using a  s
      apply (subst proots_within_pCons_1_iff)
      by simp
    have "poly p a0" 
    proof (rule ccontr)
      assume "¬ poly p a  0"
      then have "order a p >0" by (simp add: order_root)
      moreover have "order a [:a,-1:] = 1"
        by (metis (no_types, hide_lams) add.inverse_inverse add.inverse_neutral minus_pCons 
            order_power_n_n order_uminus power_one_right)
      ultimately have "order a  ([:a, - 1:] * p) > 1"
        apply (subst order_mult)
        subgoal using root.prems by auto
        subgoal by auto
        done
      then show False using ‹rsquarefree ([:a, - 1:] * p) 
        unfolding rsquarefree_def using gr_implies_not0 less_not_refl2 by blast
    qed
    then show " proots_within [:a, - 1:] s  proots_within p s = {}"
      using proots_within_pCons_1_iff(2) by auto
  qed
  finally show ?case .
qed

lemma rsquarefree_gcd_pderiv:
  fixes p::"'a::{factorial_ring_gcd,semiring_gcd_mult_normalize,semiring_char_0} poly"
  assumes "p0"
  shows "rsquarefree (p div (gcd p (pderiv p)))"
proof (cases "pderiv p = 0")
  case True
  have "poly (unit_factor p) x 0" for x 
    using unit_factor_is_unit[OF p0] 
    by (meson assms dvd_trans order_decomp poly_eq_0_iff_dvd unit_factor_dvd)
  then have "order x (unit_factor p) = 0" for x
    using order_0I by blast
  then show ?thesis using True p0 unfolding rsquarefree_def by simp
next
  case False
  define q where "q = p div (gcd p (pderiv p))"
  have "q0" unfolding q_def by (simp add: assms dvd_div_eq_0_iff)

  have order_pq:"order x p = order x q + min (order x p) (order x (pderiv p))"
    for x
  proof -
    have *:"p = q * gcd p (pderiv p)"
      unfolding q_def by simp
    show ?thesis
      apply (subst *)
      using q0 p0 ‹pderiv p0 by (simp add:order_mult order_gcd)
  qed
  have "order x q = 0  order x q=1" for x
  proof (cases "poly p x=0")
    case True
    from order_pderiv[OF p0 this] 
    have "order x p = order x (pderiv p) + 1" by simp
    then show ?thesis using order_pq[of x] by auto
  next
    case False
    then have "order x p = 0" by (simp add: order_0I)
    then have "order x q = 0" using order_pq[of x] by simp
    then show ?thesis by simp
  qed 
  then show ?thesis using q0 unfolding rsquarefree_def q_def
    by auto
qed

lemma poly_gcd_pderiv_iff:
  fixes p::"'a::{semiring_char_0,factorial_ring_gcd,semiring_gcd_mult_normalize} poly"
  shows "poly (p div (gcd p (pderiv p))) x =0  poly p x=0"
proof (cases "pderiv p=0")
  case True
  then obtain a where "p=[:a:]" using pderiv_iszero by auto
  then show ?thesis by (auto simp add: unit_factor_poly_def)
next
  case False
  then have "p0" using pderiv_0 by blast
  define q where "q = p div (gcd p (pderiv p))"
  have "q0" unfolding q_def by (simp add: p0 dvd_div_eq_0_iff)

  have order_pq:"order x p = order x q + min (order x p) (order x (pderiv p))" for x
  proof -
    have *:"p = q * gcd p (pderiv p)"
      unfolding q_def by simp
    show ?thesis
      apply (subst *)
      using q0 p0 ‹pderiv p0 by (simp add:order_mult order_gcd)
  qed

  have "order x q =0  order x p = 0" 
  proof (cases "poly p x=0")
    case True
    from order_pderiv[OF p0 this] 
    have "order x p = order x (pderiv p) + 1" by simp
    then show ?thesis using order_pq[of x] by auto
  next
    case False
    then have "order x p = 0" by (simp add: order_0I)
    then have "order x q = 0" using order_pq[of x] by simp
    then show ?thesis using ‹order x p = 0 by simp
  qed 
  then show ?thesis 
    apply (fold q_def)
    unfolding order_root using p0 q0 by auto
qed

subsection ‹Composition of a polynomial and a circular path›

lemma poly_circlepath_tan_eq:
  fixes z0::complex and r::real and p::"complex poly"
  defines "q1 fcompose p [:(z0+r)*𝗂,z0-r:] [:𝗂,1:]" and "q2  [:𝗂,1:] ^ degree p"
  assumes "0t" "t1" "t1/2"
  shows "poly p (circlepath z0 r t) = poly q1 (tan (pi*t)) / poly q2 (tan (pi*t))" 
    (is "?L = ?R")
proof -
  have "?L = poly p (z0+ r*exp (2 * of_real pi * 𝗂  * t))" 
    unfolding circlepath by simp
  also have "... = ?R"
  proof -
    define f where "f = (poly p  (λx::real. z0 + r * exp (𝗂 * x)))"
    have f_eq:"f t = ((λx::real. poly q1 x / poly q2 x) o  (λx. tan (x/2)) ) t" 
      when "cos (t / 2)  0" for t
    proof -
      have "f t = poly p (z0 + r * (cos t + 𝗂 * sin t)) " 
        unfolding f_def exp_Euler  by (auto simp add:cos_of_real sin_of_real)
      also have "... = poly p ((λx. ((z0-r)*x+(z0+r)*𝗂) / (𝗂+x)) (tan (t/2)))"
      proof -
        define tt where "tt=complex_of_real (tan (t / 2))"
        define rr where "rr = complex_of_real r"
        have "cos t = (1-tt*tt) / (1 + tt * tt)" 
             "sin t = 2*tt  / (1 + tt * tt)"
          unfolding sin_tan_half[of "t/2",simplified] cos_tan_half[of "t/2",OF that, simplified] tt_def
          by (auto simp add:power2_eq_square)
        moreover have "1 + tt * tt  0" unfolding tt_def 
          apply (fold of_real_mult)
          by (metis (no_types, hide_lams) mult_numeral_1 numeral_One of_real_add of_real_eq_0_iff
              of_real_numeral sum_squares_eq_zero_iff zero_neq_one)
        ultimately have "z0 +  r * ( (cos t) + 𝗂 * (sin t))
            =(z0*(1+tt*tt)+rr*(1-tt*tt)+𝗂*rr*2*tt ) / (1 + tt * tt) "
          apply (fold rr_def,simp add:add_divide_distrib)
          by (simp add:algebra_simps)
        also have "... = ((z0-rr)*tt+z0*𝗂+rr*𝗂) / (tt + 𝗂)"
        proof -
          have "tt + 𝗂  0" 
            using 1 + tt * tt  0 
            by (metis i_squared neg_eq_iff_add_eq_0 square_eq_iff)
          then show ?thesis 
            using 1 + tt * tt  0 by (auto simp add:divide_simps algebra_simps)
        qed
        finally have "z0 +  r * ( (cos t) + 𝗂 * (sin t)) = ((z0-rr)*tt+z0*𝗂+rr*𝗂) / (tt + 𝗂)" .
        then show ?thesis unfolding tt_def rr_def 
          by (auto simp add:algebra_simps power2_eq_square)
      qed
      also have "... = (poly p o ((λx. ((z0-r)*x+(z0+r)*𝗂) / (𝗂+x)) o (λx. tan (x/2)) )) t"
        unfolding comp_def by (auto simp:tan_of_real)
      also have "... = ((λx::real. poly q1 x / poly q2 x) o  (λx. tan (x/2)) ) t"
        unfolding q2_def q1_def
        apply (subst fcompose_poly[symmetric])
        subgoal for x
          apply simp
          by (metis Re_complex_of_real add_cancel_right_left complex_i_not_zero imaginary_unit.sel(1) plus_complex.sel(1) rcis_zero_arg rcis_zero_mod)
        subgoal by (auto simp:tan_of_real algebra_simps)
        done
      finally show ?thesis .
    qed
    
    have "cos (pi * t) 0" unfolding cos_zero_iff_int2
    proof 
      assume "i. pi * t = real_of_int i * pi + pi / 2"
      then obtain i where "pi * t = real_of_int i * pi + pi / 2" by auto
      then have "pi * t=pi * (real_of_int i + 1 / 2)" by (simp add:algebra_simps)
      then have "t=real_of_int i + 1 / 2" by auto
      then show False using 0t t1 t1/2 by auto
    qed
    from f_eq[of "2*pi*t",simplified,OF this] 
    show "?thesis"
      unfolding f_def comp_def by (auto simp add:algebra_simps)
  qed
  finally show ?thesis .
qed

end

Theory Count_Complex_Roots

(*
    Author:     Wenda Li <wl302@cam.ac.uk / liwenda1990@hotmail.com>
*)

section ‹Procedures to count the number of complex roots›

theory Count_Complex_Roots imports 
  Winding_Number_Eval.Winding_Number_Eval 
  Extended_Sturm
  More_Polynomials
  Budan_Fourier.Sturm_Multiple_Roots
begin

subsection ‹Misc›
    
(*refined version of the library one with the same name by dropping unnecessary assumptions*) 
corollary path_image_part_circlepath_subset:
  assumes "r0"
  shows "path_image(part_circlepath z r st tt)  sphere z r"
proof (cases "sttt")
  case True
  then show ?thesis 
    by (auto simp: assms path_image_part_circlepath sphere_def dist_norm algebra_simps norm_mult) 
next
  case False
  then have "path_image(part_circlepath z r tt st)  sphere z r"
    by (auto simp: assms path_image_part_circlepath sphere_def dist_norm algebra_simps norm_mult)
  moreover have "path_image(part_circlepath z r tt st) = path_image(part_circlepath z r st tt)"
    using path_image_reversepath by fastforce
  ultimately show ?thesis by auto
qed

(*refined version of the library one with the same name by dropping unnecessary assumptions*)
proposition in_path_image_part_circlepath:
  assumes "w  path_image(part_circlepath z r st tt)" "0  r"
  shows "norm(w - z) = r"
proof -
  have "w  {c. dist z c = r}"
    by (metis (no_types) path_image_part_circlepath_subset sphere_def subset_eq assms)
  thus ?thesis
    by (simp add: dist_norm norm_minus_commute)
qed  

lemma infinite_ball:
  fixes a :: "'a::euclidean_space"
  assumes "r > 0" 
  shows "infinite (ball a r)"
  using uncountable_ball[OF assms, THEN uncountable_infinite] .

lemma infinite_cball:
  fixes a :: "'a::euclidean_space"
  assumes "r > 0" 
  shows "infinite (cball a r)"
  using uncountable_cball[OF assms, THEN uncountable_infinite,of a] .

(*FIXME: to generalise*)
lemma infinite_sphere:
  fixes a :: complex
  assumes "r > 0" 
  shows "infinite (sphere a r)" 
proof -
  have "uncountable (path_image (circlepath a r))"
    apply (rule simple_path_image_uncountable)
    using simple_path_circlepath assms by simp
  then have "uncountable (sphere a r)"
    using assms by simp
  from uncountable_infinite[OF this] show ?thesis .
qed

lemma infinite_halfspace_Im_gt: "infinite {x. Im x > b}"
  apply (rule connected_uncountable[THEN uncountable_infinite,of _ "(b+1)* 𝗂" "(b+2)*𝗂"])
  by (auto intro!:convex_connected simp add: convex_halfspace_Im_gt)
       
lemma (in ring_1) Ints_minus2: "- a    a  "
  using Ints_minus[of "-a"] by auto

lemma dvd_divide_Ints_iff:
  "b dvd a  b=0  of_int a / of_int b  ( :: 'a :: {field,ring_char_0} set)"  
proof
  assume asm:"b dvd a  b=0"
  let ?thesis = "of_int a / of_int b  ( :: 'a :: {field,ring_char_0} set)"
  have ?thesis when "b dvd a"
  proof -
    obtain c where "a=b * c" using b dvd a unfolding dvd_def by auto
    then show ?thesis by (auto simp add:field_simps)
  qed
  moreover have ?thesis when "b=0" 
    using that by auto
  ultimately show ?thesis using asm by auto
next
  assume "of_int a / of_int b  ( :: 'a :: {field,ring_char_0} set)"
  from Ints_cases[OF this] obtain c where *:"(of_int::_  'a) c= of_int a / of_int b"
    by metis 
  have "b dvd a" when "b0"
  proof -
    have "(of_int::_  'a) a = of_int b * of_int c" using that * by auto
    then have "a = b * c" using of_int_eq_iff by fastforce
    then show ?thesis unfolding dvd_def by auto
  qed
  then show " b dvd a  b = 0" by auto
qed
 
lemma of_int_div_field:
  assumes "d dvd n"
  shows "(of_int::_'a::field_char_0) (n div d) = of_int n / of_int d" 
  apply (subst (2) dvd_mult_div_cancel[OF assms,symmetric])
  by (auto simp add:field_simps)

lemma powr_eq_1_iff:
  assumes "a>0"
  shows "(a::real) powr b =1  a=1  b=0"
proof 
  assume "a powr b = 1"
  have "b * ln a = 0"
    using a powr b = 1 ln_powr[of a b] assms by auto
  then have "b=0  ln a =0" by auto
  then show "a = 1  b = 0" using assms by auto
qed (insert assms, auto)

lemma tan_inj_pi:
  "- (pi/2) < x  x < pi/2  - (pi/2) < y  y < pi/2  tan x = tan y  x = y"
  by (metis arctan_tan)

(*TODO: can we avoid fcompose in the proof?*)
lemma finite_ReZ_segments_poly_circlepath:
          "finite_ReZ_segments (poly p  circlepath z0 r) 0"
proof (cases "t({0..1} - {1/2}). Re ((poly p  circlepath z0 r) t) = 0")
  case True
  have "isCont (Re  poly p  circlepath z0 r) (1/2)"
    by (auto intro!:continuous_intros simp:circlepath)
  moreover have "(Re  poly p  circlepath z0 r) 1/2  0"
  proof -
    have "F x in at (1 / 2). (Re  poly p  circlepath z0 r) x = 0" 
      unfolding eventually_at_le 
      apply (rule exI[where x="1/2"])
      unfolding dist_real_def abs_diff_le_iff
      by (auto intro!:True[rule_format, unfolded comp_def])
    then show ?thesis by (rule tendsto_eventually)
  qed
  ultimately have "Re ((poly p  circlepath z0 r) (1/2)) = 0"
    unfolding comp_def by (simp add: LIM_unique continuous_within)
  then have "t{0..1}. Re ((poly p  circlepath z0 r) t) = 0"
    using True by blast
  then show ?thesis 
    apply (rule_tac finite_ReZ_segments_constI[THEN finite_ReZ_segments_congE])
    by auto
next
  case False
   define q1 q2 where "q1=fcompose p [:(z0+r)*𝗂,z0-r:] [:𝗂,1:]" and 
                      "q2=([:𝗂, 1:] ^ degree p)"
  define q1R q1I where "q1R=map_poly Re q1" and "q1I=map_poly Im q1"
  define q2R q2I where "q2R=map_poly Re q2" and "q2I=map_poly Im q2"
  define qq where "qq=q1R*q2R + q1I*q2I"
  
  have poly_eq:"Re ((poly p  circlepath z0 r) t) = 0  poly qq (tan (pi * t)) = 0"
    when "0t" "t1" "t1/2" for t 
  proof -
    define tt where "tt=tan (pi * t)"
    have "Re ((poly p  circlepath z0 r) t) = 0  Re (poly q1 tt / poly q2 tt) = 0"
      unfolding comp_def
      apply (subst poly_circlepath_tan_eq[of t p z0 r,folded q1_def q2_def tt_def])
      using that by simp_all
    also have "...  poly q1R tt * poly q2R tt + poly q1I tt * poly q2I tt = 0"
      unfolding q1I_def q1R_def q2R_def q2I_def
      by (simp add:Re_complex_div_eq_0 Re_poly_of_real Im_poly_of_real)
    also have "...  poly qq tt = 0"
      unfolding qq_def by simp
    finally show ?thesis unfolding tt_def .
  qed

  have "finite {t. Re ((poly p  circlepath z0 r) t) = 0  0  t  t  1}"
  proof - 
    define P where "P=(λt. Re ((poly p  circlepath z0 r) t) = 0)"
    define A where "A= ({0..1}::real set)"
    define S where "S={tA-{1,1/2}. P t}"
    have "finite {t. poly qq (tan (pi * t)) = 0  0  t  t < 1  t1/2}"
    proof -
      define A where "A={t::real. 0  t  t < 1  t  1 / 2}"
      have "finite ((λt. tan (pi * t))  -`  {x. poly qq x=0}  A)"
      proof (rule finite_vimage_IntI)
        have "x = y" when "tan (pi * x) = tan (pi * y)" "xA" "yA" for x y
        proof -
          define x' where "x'=(if x<1/2 then x else x-1)"
          define y' where "y'=(if y<1/2 then y else y-1)"
          have "x'*pi = y'*pi" 
          proof (rule tan_inj_pi)
            have *:"- 1 / 2 < x'" "x' < 1 / 2" "- 1 / 2 < y'" "y' < 1 / 2" 
              using that(2,3) unfolding x'_def y'_def A_def by simp_all
            show "- (pi / 2) < x'  * pi" "x'  * pi < pi / 2" "- (pi / 2) < y'  * pi" 
                  "y'*pi < pi / 2"
              using mult_strict_right_mono[OF *(1),of pi] 
                    mult_strict_right_mono[OF *(2),of pi] 
                    mult_strict_right_mono[OF *(3),of pi] 
                    mult_strict_right_mono[OF *(4),of pi] 
              by auto
          next
            have "tan (x' * pi) = tan (x * pi)"
              unfolding x'_def using tan_periodic_int[of _ "- 1",simplified]
              by (auto simp add:algebra_simps)
            also have "... = tan (y * pi)"
              using ‹tan (pi * x) = tan (pi * y) by (auto simp:algebra_simps)
            also have "... = tan (y' * pi)"
              unfolding y'_def using tan_periodic_int[of _ "- 1",simplified]
              by (auto simp add:algebra_simps)
            finally show "tan (x' * pi) = tan (y' * pi)" .
          qed
          then have "x'=y'" by auto
          then show ?thesis 
            using that(2,3) unfolding x'_def y'_def A_def by (auto split:if_splits)
        qed
        then show "inj_on (λt. tan (pi * t)) A"
          unfolding inj_on_def by blast
      next
        have "qq0"
        proof (rule ccontr)
          assume "¬ qq  0"
          then have "Re ((poly p  circlepath z0 r) t) = 0" when "t{0..1} - {1/2}" for t
            apply (subst poly_eq)
            using that by auto
          then show False using False by blast
        qed
        then show "finite {x. poly qq x = 0}" by (simp add: poly_roots_finite)
      qed
      then show ?thesis by (elim rev_finite_subset) (auto simp:A_def)
    qed
    moreover have "{t. poly qq (tan (pi * t)) = 0  0  t  t < 1  t1/2} = S"
      unfolding S_def P_def A_def using poly_eq by force
    ultimately have "finite S" by blast
    then have "finite (S  (if P 1 then {1} else {})  (if P (1/2) then {1/2} else {}))"
      by auto
    moreover have "(S  (if P 1 then {1} else {})  (if P (1/2) then {1/2} else {}))
                      = {t. P t  0  t  t  1}" 
    proof -
      have "1A" "1/2 A" unfolding A_def by auto
      then have "(S  (if P 1 then {1} else {})  (if P (1/2) then {1/2} else {}))
                      = {tA. P t}"
        unfolding S_def
        apply auto
        by (metis eq_divide_eq_numeral1(1) zero_neq_numeral)+
      also have "... = {t. P t  0  t  t  1}"
        unfolding A_def by auto
      finally show ?thesis .
    qed
    ultimately have "finite {t. P t  0  t  t  1}" by auto
    then show ?thesis unfolding P_def by simp
  qed
  then show ?thesis 
    apply (rule_tac finite_imp_finite_ReZ_segments)
    by auto
qed

subsection ‹Some useful conformal/@{term bij_betw} properties›

lemma bij_betw_plane_ball:"bij_betw (λx. (𝗂-x)/(𝗂+x)) {x. Im x>0} (ball 0 1)"
proof (rule bij_betw_imageI)
  have neq:"𝗂 + x 0" when "Im x>0" for x 
    using that 
    by (metis add_less_same_cancel2 add_uminus_conv_diff diff_0 diff_add_cancel 
        imaginary_unit.simps(2) not_one_less_zero uminus_complex.sel(2))
  then show "inj_on (λx. (𝗂 - x) / (𝗂 + x)) {x. 0 < Im x}"
    unfolding inj_on_def by (auto simp add:divide_simps algebra_simps)
  have "cmod ((𝗂 - x) / (𝗂 + x)) < 1" when "0 < Im x " for x
  proof -
    have "cmod (𝗂 - x) < cmod (𝗂 + x)" 
      unfolding norm_lt inner_complex_def using that 
      by (auto simp add:algebra_simps)
    then show ?thesis
      unfolding norm_divide using neq[OF that] by auto
  qed
  moreover have "x  (λx. (𝗂 - x) / (𝗂 + x)) ` {x. 0 < Im x}" when "cmod x < 1" for x 
  proof (rule rev_image_eqI[of "𝗂*(1-x)/(1+x)"])
    have "1 + x0" "𝗂 * 2 + 𝗂 * (x * 2) 0" 
      subgoal using that by (metis complex_mod_triangle_sub norm_one norm_zero not_le pth_7(1))
      subgoal using that by (metis 1 + x  0 complex_i_not_zero div_mult_self4 mult_2 
            mult_zero_right nonzero_mult_div_cancel_left nonzero_mult_div_cancel_right 
            one_add_one zero_neq_numeral)
      done        
    then show "x = (𝗂 - 𝗂 * (1 - x) / (1 + x)) / (𝗂 + 𝗂 * (1 - x) / (1 + x))"
      by (auto simp add:field_simps)
    show " 𝗂 * (1 - x) / (1 + x)  {x. 0 < Im x}"
      apply (auto simp:Im_complex_div_gt_0 algebra_simps)
      using that unfolding cmod_def by (auto simp:power2_eq_square)
  qed
  ultimately show "(λx. (𝗂 - x) / (𝗂 + x)) ` {x. 0 < Im x} = ball 0 1"
    by auto
qed
    
lemma bij_betw_axis_sphere:"bij_betw (λx. (𝗂-x)/(𝗂+x)) {x. Im x=0} (sphere 0 1 - {-1})"
proof (rule bij_betw_imageI)
  have neq:"𝗂 + x 0" when "Im x=0" for x 
    using that 
    by (metis add_diff_cancel_left' imaginary_unit.simps(2) minus_complex.simps(2) 
         right_minus_eq zero_complex.simps(2) zero_neq_one)
  then show "inj_on (λx. (𝗂 - x) / (𝗂 + x)) {x. Im x = 0}"
    unfolding inj_on_def by (auto simp add:divide_simps algebra_simps)
  have "cmod ((𝗂 - x) / (𝗂 + x)) = 1" "(𝗂 - x) / (𝗂 + x)  - 1" when "Im x = 0" for x 
  proof -
    have "cmod (𝗂 + x) = cmod (𝗂 - x)" 
      using that unfolding cmod_def by auto
    then show "cmod ((𝗂 - x) / (𝗂 + x)) = 1"
      unfolding norm_divide using neq[OF that] by auto
    show "(𝗂 - x) / (𝗂 + x)  - 1" using neq[OF that] by (auto simp add:divide_simps)
  qed
  moreover have "x  (λx. (𝗂 - x) / (𝗂 + x)) ` {x. Im x = 0}" 
    when "cmod x = 1" "x-1" for x 
  proof (rule rev_image_eqI[of "𝗂*(1-x)/(1+x)"])
    have "1 + x0" "𝗂 * 2 + 𝗂 * (x * 2) 0" 
      subgoal using that(2) by algebra 
      subgoal using that by (metis 1 + x  0 complex_i_not_zero div_mult_self4 mult_2 
            mult_zero_right nonzero_mult_div_cancel_left nonzero_mult_div_cancel_right 
            one_add_one zero_neq_numeral)
      done        
    then show "x = (𝗂 - 𝗂 * (1 - x) / (1 + x)) / (𝗂 + 𝗂 * (1 - x) / (1 + x))"
      by (auto simp add:field_simps)
    show " 𝗂 * (1 - x) / (1 + x)  {x. Im x = 0}"
      apply (auto simp:algebra_simps Im_complex_div_eq_0)
      using that(1) unfolding cmod_def by (auto simp:power2_eq_square)
  qed
  ultimately show "(λx. (𝗂 - x) / (𝗂 + x)) ` {x. Im x = 0} = sphere 0 1 - {- 1}"
    by force
qed

lemma bij_betw_ball_uball:
  assumes "r>0"
  shows "bij_betw (λx. complex_of_real r*x + z0) (ball 0 1) (ball z0 r)"
proof (rule bij_betw_imageI)
  show "inj_on (λx. complex_of_real r * x + z0) (ball 0 1)"
    unfolding inj_on_def using assms by simp
  have "dist z0 (complex_of_real r * x + z0) < r" when "cmod x<1" for x 
    using that assms by (auto simp:dist_norm norm_mult abs_of_pos)
  moreover have "x  (λx. complex_of_real r * x + z0) ` ball 0 1" when "dist z0 x < r" for x 
    apply (rule rev_image_eqI[of "(x-z0)/r"])
    using that assms by (auto simp add: dist_norm norm_divide norm_minus_commute)
  ultimately show "(λx. complex_of_real r  * x + z0) ` ball 0 1 = ball z0 r" 
    by auto
qed

lemma bij_betw_sphere_usphere:
  assumes "r>0"
  shows "bij_betw (λx. complex_of_real r*x + z0) (sphere 0 1) (sphere z0 r)"
proof (rule bij_betw_imageI)
  show "inj_on (λx. complex_of_real r * x + z0) (sphere 0 1)"
    unfolding inj_on_def using assms by simp
  have "dist z0 (complex_of_real r * x + z0) = r" when "cmod x=1" for x 
    using that assms by (auto simp:dist_norm norm_mult abs_of_pos)
  moreover have "x  (λx. complex_of_real r * x + z0) ` sphere 0 1" when "dist z0 x = r" for x 
    apply (rule rev_image_eqI[of "(x-z0)/r"])
    using that assms by (auto simp add: dist_norm norm_divide norm_minus_commute)
  ultimately show "(λx. complex_of_real r  * x + z0) ` sphere 0 1 = sphere z0 r" 
    by auto
qed

lemma proots_ball_plane_eq:
  defines "q1[:𝗂,-1:]" and "q2[:𝗂,1:]"
  assumes "p0"
  shows "proots_count p (ball 0 1) = proots_count (fcompose p q1 q2) {x. 0 < Im x}"
  unfolding q1_def q2_def 
proof (rule proots_fcompose_bij_eq[OF _ p0])
  show "x{x. 0 < Im x}. poly [:𝗂, 1:] x  0" 
    apply simp 
    by (metis add_less_same_cancel2 imaginary_unit.simps(2) not_one_less_zero 
          plus_complex.simps(2) zero_complex.simps(2))
  show "infinite (UNIV::complex set)" by (simp add: infinite_UNIV_char_0)
qed (use bij_betw_plane_ball in auto)

lemma proots_sphere_axis_eq:
  defines "q1[:𝗂,-1:]" and "q2[:𝗂,1:]"
  assumes "p0"
  shows "proots_count p (sphere 0 1 - {- 1}) = proots_count (fcompose p q1 q2) {x. 0 = Im x}"
unfolding q1_def q2_def 
proof (rule proots_fcompose_bij_eq[OF _ p0])
  show "x{x. 0 = Im x}. poly [:𝗂, 1:] x  0" by (simp add: Complex_eq_0 plus_complex.code)
  show "infinite (UNIV::complex set)" by (simp add: infinite_UNIV_char_0)
qed (use bij_betw_axis_sphere in auto)

lemma proots_card_ball_plane_eq:
  defines "q1[:𝗂,-1:]" and "q2[:𝗂,1:]"
  assumes "p0"
  shows "card (proots_within p (ball 0 1)) = card (proots_within (fcompose p q1 q2) {x. 0 < Im x})"
unfolding q1_def q2_def
proof (rule proots_card_fcompose_bij_eq[OF _ p0])
  show "x{x. 0 < Im x}. poly [:𝗂, 1:] x  0" 
    apply simp 
    by (metis add_less_same_cancel2 imaginary_unit.simps(2) not_one_less_zero 
          plus_complex.simps(2) zero_complex.simps(2))
qed (use bij_betw_plane_ball infinite_UNIV_char_0 in auto)

lemma proots_card_sphere_axis_eq:
  defines "q1[:𝗂,-1:]" and "q2[:𝗂,1:]"
  assumes "p0"
  shows "card (proots_within p (sphere 0 1 - {- 1})) 
            = card (proots_within (fcompose p q1 q2) {x. 0 = Im x})"
unfolding q1_def q2_def
proof (rule proots_card_fcompose_bij_eq[OF _ p0])
  show "x{x. 0 = Im x}. poly [:𝗂, 1:] x  0" by (simp add: Complex_eq_0 plus_complex.code)
qed (use bij_betw_axis_sphere infinite_UNIV_char_0 in auto)

lemma proots_uball_eq:
  fixes z0::complex and r::real
  defines "q[:z0, of_real r:]"
  assumes "p0" and "r>0"
  shows "proots_count p (ball z0 r) = proots_count (p p q) (ball 0 1)"
proof -
  show ?thesis
    apply (rule proots_pcompose_bij_eq[OF _ p0])
    subgoal unfolding q_def using bij_betw_ball_uball[OF r>0,of z0] by (auto simp:algebra_simps)
    subgoal unfolding q_def using r>0 by auto
    done
qed

lemma proots_card_uball_eq:
  fixes z0::complex and r::real
  defines "q[:z0, of_real r:]"
  assumes "r>0"
  shows "card (proots_within p (ball z0 r)) = card (proots_within (p p q) (ball 0 1))"
proof -
  have ?thesis 
    when "p=0"
  proof -
    have "card (ball z0 r) = 0" "card (ball (0::complex) 1) = 0"
      using infinite_ball[OF r>0,of z0] infinite_ball[of 1 "0::complex"] by auto 
    then show ?thesis using that by auto
  qed
  moreover have ?thesis 
    when "p0"
    apply (rule proots_card_pcompose_bij_eq[OF _ p0])
    subgoal unfolding q_def using bij_betw_ball_uball[OF r>0,of z0] by (auto simp:algebra_simps)
    subgoal unfolding q_def using r>0 by auto
    done
  ultimately show ?thesis 
    by blast
qed

lemma proots_card_usphere_eq:
  fixes z0::complex and r::real
  defines "q[:z0, of_real r:]"
  assumes "r>0"
  shows "card (proots_within p (sphere z0 r)) = card (proots_within (p p q) (sphere 0 1))"
proof -
  have ?thesis 
    when "p=0"
  proof -
    have "card (sphere z0 r) = 0" "card (sphere (0::complex) 1) = 0"
      using infinite_sphere[OF r>0,of z0] infinite_sphere[of 1 "0::complex"] by auto 
    then show ?thesis using that by auto
  qed
  moreover have ?thesis
    when "p0"
    apply (rule proots_card_pcompose_bij_eq[OF _ p0])
    subgoal unfolding q_def using bij_betw_sphere_usphere[OF r>0,of z0] 
      by (auto simp:algebra_simps)
    subgoal unfolding q_def using r>0 by auto
    done
  ultimately show "card (proots_within p (sphere z0 r)) = card (proots_within (p p q) (sphere 0 1))" 
    by blast
qed
  
subsection ‹Combining two real polynomials into a complex one›  

definition cpoly_of:: "real poly  real poly  complex poly" where
  "cpoly_of pR pI = map_poly of_real pR + smult 𝗂 (map_poly of_real pI)"

lemma cpoly_of_eq_0_iff[iff]:
  "cpoly_of pR pI = 0  pR = 0  pI = 0"
proof -
  have "pR = 0  pI = 0" when "cpoly_of pR pI = 0"
  proof -
    have "complex_of_real (coeff pR n) + 𝗂 * complex_of_real (coeff pI n) = 0" for n
      using that unfolding poly_eq_iff cpoly_of_def by (auto simp:coeff_map_poly)
    then have "coeff pR n = 0  coeff pI n = 0" for n
      by (metis Complex_eq Im_complex_of_real Re_complex_of_real complex.sel(1) complex.sel(2) 
          of_real_0)
    then show ?thesis unfolding poly_eq_iff by auto 
  qed
  then show ?thesis by (auto simp:cpoly_of_def)
qed

lemma cpoly_of_decompose:
    "p = cpoly_of (map_poly Re p) (map_poly Im p)"
  unfolding cpoly_of_def 
  apply (induct p)
  by (auto simp add:map_poly_pCons map_poly_map_poly complex_eq)

lemma cpoly_of_dist_right:
    "cpoly_of (pR*g) (pI*g) = cpoly_of pR pI * (map_poly of_real g)"
  unfolding cpoly_of_def by (simp add: distrib_right)
  
lemma poly_cpoly_of_real:
    "poly (cpoly_of pR pI) (of_real x) = Complex (poly pR x) (poly pI x)"
  unfolding cpoly_of_def by (simp add: Complex_eq of_real_poly_map_poly)
    
lemma poly_cpoly_of_real_iff:
  shows "poly (cpoly_of pR pI) (of_real t) =0  poly pR t = 0  poly pI t=0 "  
  unfolding  poly_cpoly_of_real using Complex_eq_0 by blast  

lemma order_cpoly_gcd_eq:
  assumes "pR0  pI0"
  shows "order t (cpoly_of pR pI) = order t (gcd pR pI)"
proof -
  define g where "g = gcd pR pI"
  have [simp]:"g0" unfolding g_def using assms by auto
  obtain pr pi where pri: "pR = pr * g" "pI = pi * g" "coprime pr pi"
    unfolding g_def using assms(1) gcd_coprime_exists g  0 g_def by blast
  then have "pr 0  pi 0" using assms mult_zero_left by blast

  have "order t (cpoly_of pR pI) = order t (cpoly_of pr pi * (map_poly of_real g))"
    unfolding pri cpoly_of_dist_right by simp
  also have "... = order t (cpoly_of pr pi) + order t g" 
    apply (subst order_mult)
    using pr 0  pi 0 by (auto simp:map_poly_order_of_real)
  also have "... = order t g"
  proof -
    have "poly (cpoly_of pr pi) t 0" unfolding poly_cpoly_of_real_iff
      using ‹coprime pr pi coprime_poly_0 by blast
    then have "order t (cpoly_of pr pi) = 0" by (simp add: order_0I)
    then show ?thesis by auto
  qed
  finally show ?thesis unfolding g_def .
qed

subsection ‹Number of roots on a (bounded or unbounded) segment›

― ‹1 dimensional hyperplane›
definition unbounded_line::"'a::real_vector  'a  'a set" where 
   "unbounded_line a b = ({x. u::real. x= (1 - u) *R a + u *R b})"

definition proots_line_card:: "complex poly  complex  complex  nat" where
  "proots_line_card p st tt = card (proots_within p (open_segment st tt))"

definition proots_unbounded_line_card:: "complex poly  complex  complex  nat" where
  "proots_unbounded_line_card p st tt = card (proots_within p (unbounded_line st tt))"

definition proots_unbounded_line :: "complex poly  complex  complex  nat" where
  "proots_unbounded_line p st tt = proots_count p (unbounded_line st tt)"

lemma card_proots_open_segments:
  assumes "poly p st 0" "poly p tt  0"
  shows "card (proots_within p (open_segment st tt)) = 
                (let pc = pcompose p [:st, tt - st:];
                     pR = map_poly Re pc;
                     pI = map_poly Im pc;
                     g  = gcd pR pI
                 in changes_itv_smods 0 1 g (pderiv g))" (is "?L = ?R")
proof -
  define pc pR pI g where 
      "pc = pcompose p [:st, tt-st:]" and
      "pR = map_poly Re pc" and
      "pI = map_poly Im pc" and
      "g  = gcd pR pI"
  have poly_iff:"poly g t=0  poly pc t =0" for t
  proof -
    have "poly g t = 0  poly pR t =0  poly pI t =0" 
      unfolding g_def using poly_gcd_iff by auto
    also have "...  poly pc t =0"
    proof -
      have "cpoly_of pR pI = pc"
        unfolding pc_def pR_def pI_def using cpoly_of_decompose by auto
      then show ?thesis using poly_cpoly_of_real_iff by blast
    qed
    finally show ?thesis by auto
  qed      

  have "?R = changes_itv_smods 0 1 g (pderiv g)"
    unfolding pc_def g_def pI_def pR_def by (auto simp add:Let_def)
  also have "... = card {t. poly g t = 0  0 < t  t < 1}"
  proof -
    have "poly g 0  0" 
      using poly_iff[of 0] assms unfolding pc_def by (auto simp add:poly_pcompose)
    moreover have "poly g 1  0" 
      using poly_iff[of 1] assms unfolding pc_def by (auto simp add:poly_pcompose)
    ultimately show ?thesis using sturm_interval[of 0 1 g] by auto
  qed
  also have "... = card {t::real. poly pc t = 0  0 < t  t < 1}"
    unfolding poly_iff by simp
  also have "... = ?L"
  proof (cases "st=tt")
    case True
    then show ?thesis unfolding pc_def poly_pcompose using ‹poly p tt  0
      by auto
  next
    case False
    define ff where "ff = (λt::real. st + t*(tt-st))"
    define ll where "ll = {t. poly pc (complex_of_real t) = 0  0 < t  t < 1}"
    have "ff ` ll = proots_within p (open_segment st tt)"
    proof (rule equalityI)
      show "ff ` ll  proots_within p (open_segment st tt)"
        unfolding ll_def ff_def pc_def poly_pcompose 
        by (auto simp add:in_segment False scaleR_conv_of_real algebra_simps)
    next
      show "proots_within p (open_segment st tt)  ff ` ll"
      proof clarify
        fix x assume asm:"x  proots_within p (open_segment st tt)" 
        then obtain u where "0<u" and "u < 1" and u:"x = (1 - u) *R st + u *R tt"
          by (auto simp add:in_segment)
        then have "poly p ((1 - u) *R st + u *R tt) = 0" using asm by simp
        then have "u  ll"
          unfolding ll_def pc_def poly_pcompose 
          by (simp add:scaleR_conv_of_real algebra_simps 0<u u<1)
        moreover have "x = ff u"
          unfolding ff_def using u by (auto simp add:algebra_simps scaleR_conv_of_real)
        ultimately show "x  ff ` ll" by (rule rev_image_eqI[of "u"])
      qed
    qed
    moreover have "inj_on ff ll"
      unfolding ff_def using False inj_on_def by fastforce
    ultimately show ?thesis unfolding ll_def 
      using card_image[of ff] by fastforce
  qed
  finally show ?thesis by simp
qed

lemma unbounded_line_closed_segment: "closed_segment a b  unbounded_line a b"
  unfolding unbounded_line_def closed_segment_def by auto

lemma card_proots_unbounded_line:
  assumes "sttt"
  shows "card (proots_within p (unbounded_line st tt)) = 
                (let pc = pcompose p [:st, tt - st:];
                     pR = map_poly Re pc;
                     pI = map_poly Im pc;
                     g  = gcd pR pI
                 in nat (changes_R_smods g (pderiv g)))" (is "?L = ?R")
proof -
  define pc pR pI g where 
      "pc = pcompose p [:st, tt-st:]" and
      "pR = map_poly Re pc" and
      "pI = map_poly Im pc" and
      "g  = gcd pR pI"
  have poly_iff:"poly g t=0  poly pc t =0" for t
  proof -
    have "poly g t = 0  poly pR t =0  poly pI t =0" 
      unfolding g_def using poly_gcd_iff by auto
    also have "...  poly pc t =0"
    proof -
      have "cpoly_of pR pI = pc"
        unfolding pc_def pR_def pI_def using cpoly_of_decompose by auto
      then show ?thesis using poly_cpoly_of_real_iff by blast
    qed
    finally show ?thesis by auto
  qed      

  have "?R = nat (changes_R_smods g (pderiv g))"
    unfolding pc_def g_def pI_def pR_def by (auto simp add:Let_def)
  also have "... = card {t. poly g t = 0}"
    using sturm_R[of g] by simp
  also have "... = card {t::real. poly pc t = 0}"
    unfolding poly_iff by simp
  also have "... = ?L"
  proof (cases "st=tt")
    case True
    then show ?thesis unfolding pc_def poly_pcompose unbounded_line_def using assms
      by (auto simp add:proots_within_def)
  next
    case False
    define ff where "ff = (λt::real. st + t*(tt-st))"
    define ll where "ll = {t. poly pc (complex_of_real t) = 0}"
    have "ff ` ll = proots_within p (unbounded_line st tt)"
    proof (rule equalityI)
      show "ff ` ll  proots_within p (unbounded_line st tt)"
        unfolding ll_def ff_def pc_def poly_pcompose 
        by (auto simp add:unbounded_line_def False scaleR_conv_of_real algebra_simps)
    next
      show "proots_within p (unbounded_line st tt)  ff ` ll"
      proof clarify
        fix x assume asm:"x  proots_within p (unbounded_line st tt)" 
        then obtain u where u:"x = (1 - u) *R st + u *R tt"
          by (auto simp add:unbounded_line_def)
        then have "poly p ((1 - u) *R st + u *R tt) = 0" using asm by simp
        then have "u  ll"
          unfolding ll_def pc_def poly_pcompose 
          by (simp add:scaleR_conv_of_real algebra_simps unbounded_line_def)
        moreover have "x = ff u"
          unfolding ff_def using u by (auto simp add:algebra_simps scaleR_conv_of_real)
        ultimately show "x  ff ` ll" by (rule rev_image_eqI[of "u"])
      qed
    qed
    moreover have "inj_on ff ll"
      unfolding ff_def using False inj_on_def by fastforce
    ultimately show ?thesis unfolding ll_def 
      using card_image[of ff] by metis
  qed  
  finally show ?thesis by simp
qed

lemma proots_unbounded_line:
  assumes "sttt" "p0"
  shows "(proots_count p (unbounded_line st tt)) = 
                (let pc = pcompose p [:st, tt - st:];
                     pR = map_poly Re pc;
                     pI = map_poly Im pc;
                     g  = gcd pR pI
                 in nat (changes_R_smods_ext g (pderiv g)))" (is "?L = ?R")
proof -
  define pc pR pI g where 
      "pc = pcompose p [:st, tt-st:]" and
      "pR = map_poly Re pc" and
      "pI = map_poly Im pc" and
      "g  = gcd pR pI"
  have [simp]: "g0" "pc0"
  proof -
    show "pc0" using assms(1) assms(2) pc_def pcompose_eq_0 by fastforce
    then have "pR0  pI0" unfolding pR_def pI_def by (metis cpoly_of_decompose map_poly_0)
    then show "g0" unfolding g_def by simp
  qed
  have order_eq:"order t g = order t pc" for t
    apply (subst order_cpoly_gcd_eq[of pR pI,folded g_def,symmetric])
    subgoal using g0 unfolding g_def by simp
    subgoal unfolding pR_def pI_def by (simp add:cpoly_of_decompose[symmetric])
    done

  have "?R = nat (changes_R_smods_ext g (pderiv g))"
    unfolding pc_def g_def pI_def pR_def by (auto simp add:Let_def)
  also have "... = proots_count g UNIV"
    using sturm_ext_R[OF g0] by auto
  also have "... = proots_count (map_poly complex_of_real g) (of_real ` UNIV)"
    apply (subst proots_count_of_real)
    by auto
  also have "... = proots_count (map_poly complex_of_real g) {x. Im x = 0}"
    apply (rule arg_cong2[where f=proots_count])
    using Reals_def complex_is_Real_iff by auto
  also have "... = proots_count pc {x. Im x = 0}"
    apply (rule proots_count_cong)
    apply (metis (mono_tags) Im_complex_of_real Re_complex_of_real g  0 complex_surj 
            map_poly_order_of_real mem_Collect_eq order_eq)
    by auto
  also have "... = proots_count p (unbounded_line st tt)"
  proof -
    have "poly [:st, tt - st:] ` {x. Im x = 0} = unbounded_line st tt"
      unfolding unbounded_line_def 
      apply safe
      subgoal for _ x 
        apply (rule_tac x="Re x" in exI) 
        apply (simp add:algebra_simps)
        by (simp add: mult.commute scaleR_complex.code times_complex.code)
      subgoal for _ u
        apply (rule rev_image_eqI[of "of_real u"])
        by (auto simp:scaleR_conv_of_real algebra_simps)
      done
    then show ?thesis 
      unfolding pc_def 
      apply (subst proots_pcompose)
      using p0 sttt by auto
  qed
  finally show ?thesis by simp
qed

lemma proots_unbounded_line_card_code[code]:
  "proots_unbounded_line_card p st tt = 
              (if sttt then 
                (let pc = pcompose p [:st, tt - st:];
                     pR = map_poly Re pc;
                     pI = map_poly Im pc;
                     g  = gcd pR pI
                 in nat (changes_R_smods g (pderiv g))) 
              else 
                  Code.abort (STR ''proots_unbounded_line_card fails due to invalid hyperplanes.'') 
                      (λ_. proots_unbounded_line_card p st tt))"
  unfolding proots_unbounded_line_card_def using card_proots_unbounded_line[of st tt p] by auto

lemma proots_unbounded_line_code[code]:
  "proots_unbounded_line p st tt = 
              ( if sttt then 
                if p0 then 
                  (let pc = pcompose p [:st, tt - st:];
                     pR = map_poly Re pc;
                     pI = map_poly Im pc;
                     g  = gcd pR pI
                  in nat (changes_R_smods_ext g (pderiv g)))
                else 
                  Code.abort (STR ''proots_unbounded_line fails due to p=0'') 
                      (λ_. proots_unbounded_line p st tt)
                else 
                  Code.abort (STR ''proots_unbounded_line fails due to invalid hyperplanes.'') 
                      (λ_. proots_unbounded_line p st tt) )"
  unfolding proots_unbounded_line_def using proots_unbounded_line by auto
  
subsection ‹Checking if there a polynomial root on a closed segment›    
    
definition no_proots_line::"complex poly  complex  complex  bool" where
  "no_proots_line p st tt = (proots_within p (closed_segment st tt) = {})"

(*TODO: the proof can probably be simplified using Lemma card_proots_open_segments*)
lemma no_proots_line_code[code]: "no_proots_line p st tt = (if poly p st 0  poly p tt  0 then 
                (let pc = pcompose p [:st, tt - st:];
                     pR = map_poly Re pc;
                     pI = map_poly Im pc;
                     g  = gcd pR pI
                 in if changes_itv_smods 0 1 g (pderiv g) = 0 then True else False) else False)"
            (is "?L = ?R")
proof (cases "poly p st 0  poly p tt  0")
  case False
  thus ?thesis unfolding no_proots_line_def by auto
next
  case True
  then have "poly p st  0" "poly p tt  0" by auto
  define pc pR pI g where 
      "pc = pcompose p [:st, tt-st:]" and
      "pR = map_poly Re pc" and
      "pI = map_poly Im pc" and
      "g  = gcd pR pI"
  have poly_iff:"poly g t=0  poly pc t =0" for t
  proof -
    have "poly g t = 0  poly pR t =0  poly pI t =0" 
      unfolding g_def using poly_gcd_iff by auto
    also have "...  poly pc t =0"
    proof -
      have "cpoly_of pR pI = pc"
        unfolding pc_def pR_def pI_def using cpoly_of_decompose by auto
      then show ?thesis using poly_cpoly_of_real_iff by blast
    qed
    finally show ?thesis by auto
  qed      
  have "?R = (changes_itv_smods 0 1 g (pderiv g) = 0)" 
    using True unfolding pc_def g_def pI_def pR_def
    by (auto simp add:Let_def)      
  also have "... = (card {x. poly g x = 0  0 < x  x < 1} = 0)"  
  proof -
    have "poly g 0  0" 
      using poly_iff[of 0] True unfolding pc_def by (auto simp add:poly_pcompose)
    moreover have "poly g 1  0" 
      using poly_iff[of 1] True unfolding pc_def by (auto simp add:poly_pcompose)
    ultimately show ?thesis using sturm_interval[of 0 1 g] by auto
  qed
  also have "... = ({x. poly g x = 0  0 < x  x < 1} = {})"
  proof -
    have "g0"
    proof (rule ccontr)
      assume "¬ g  0"
      then have "poly pc 0 =0"
        using poly_iff[of 0] by auto
      then show False using True unfolding pc_def by (auto simp add:poly_pcompose)
    qed
    from poly_roots_finite[OF this] have "finite {x. poly g x = 0  0 < x  x < 1}"
      by auto
    then show ?thesis using card_eq_0_iff by auto
  qed
  also have "... = ?L"
  proof -
    have "(t. poly g t = 0  0 < t  t < 1)  (t::real. poly pc t = 0  0 < t  t < 1)"
      using poly_iff by auto
    also have "...  (x. x  closed_segment st tt  poly p x = 0)" 
    proof 
      assume "t. poly pc (complex_of_real t) = 0  0 < t  t < 1"
      then obtain t where *:"poly pc (of_real t) = 0" and "0 < t" "t < 1" by auto
      define x where "x=poly [:st, tt - st:] t"
      have "x  closed_segment st tt" using 0<t t<1 unfolding x_def in_segment
        by (intro exI[where x=t],auto simp add: algebra_simps scaleR_conv_of_real)
      moreover have "poly p x=0" using * unfolding pc_def x_def
        by (auto simp add:poly_pcompose)
      ultimately show "x. x  closed_segment st tt  poly p x = 0" by auto
    next
      assume "x. x  closed_segment st tt  poly p x = 0"
      then obtain x where "x  closed_segment st tt" "poly p x = 0" by auto
      then obtain t::real where *:"x = (1 - t) *R st + t *R tt" and "0t" "t1"  
        unfolding in_segment by auto
      then have "x=poly [:st, tt - st:] t" by (auto simp add: algebra_simps scaleR_conv_of_real)
      then have "poly pc (complex_of_real t) = 0" 
        using ‹poly p x=0 unfolding pc_def by (auto simp add:poly_pcompose)
      moreover have "t0" "t1" using True *  ‹poly p x=0 by auto
      then have "0<t" "t<1" using 0t t1 by auto
      ultimately show "t. poly pc (complex_of_real t) = 0  0 < t  t < 1" by auto
    qed        
    finally show ?thesis
      unfolding no_proots_line_def proots_within_def 
      by blast
  qed
  finally show ?thesis by simp
qed
   
subsection ‹Counting roots in a rectangle›  
  
definition proots_rectangle ::"complex poly  complex  complex  nat" where
  "proots_rectangle p lb ub = proots_count p (box lb ub)"  
  
lemma closed_segment_imp_Re_Im:
  fixes x::complex
  assumes "xclosed_segment lb ub" 
  shows "Re lb  Re ub  Re lb  Re x  Re x  Re ub" 
        "Im lb  Im ub  Im lb  Im x  Im x  Im ub"
proof -
  obtain u where x_u:"x=(1 - u) *R lb + u *R ub " and "0  u" "u  1"
    using assms unfolding closed_segment_def by auto
  have "Re lb  Re x" when "Re lb  Re ub"
  proof -
    have "Re x = Re ((1 - u) *R lb + u *R ub)"
      using x_u by blast
    also have "... = Re (lb + u *R (ub - lb))" by (auto simp add:algebra_simps)
    also have "... = Re lb + u * (Re ub - Re lb)" by auto
    also have "...  Re lb" using u0 ‹Re lb  Re ub by auto
    finally show ?thesis .
  qed
  moreover have "Im lb  Im x" when "Im lb  Im ub"
  proof -
    have "Im x = Im ((1 - u) *R lb + u *R ub)"
      using x_u by blast
    also have "... = Im (lb + u *R (ub - lb))" by (auto simp add:algebra_simps)
    also have "... = Im lb + u * (Im ub - Im lb)" by auto
    also have "...  Im lb" using u0 ‹Im lb  Im ub by auto
    finally show ?thesis .
  qed
  moreover have "Re x  Re ub" when "Re lb  Re ub"
  proof -
    have "Re x = Re ((1 - u) *R lb + u *R ub)"
      using x_u by blast
    also have "... = (1 - u) * Re lb + u * Re ub" by auto
    also have "...  (1 - u) * Re ub + u * Re ub"
      using u1 ‹Re lb  Re ub by (auto simp add: mult_left_mono)
    also have "... = Re ub" by (auto simp add:algebra_simps)
    finally show ?thesis .
  qed
  moreover have "Im x  Im ub" when "Im lb  Im ub"
  proof -
    have "Im x = Im ((1 - u) *R lb + u *R ub)"
      using x_u by blast
    also have "... = (1 - u) * Im lb + u * Im ub" by auto
    also have "...  (1 - u) * Im ub + u * Im ub"
      using u1 ‹Im lb  Im ub by (auto simp add: mult_left_mono)
    also have "... = Im ub" by (auto simp add:algebra_simps)
    finally show ?thesis .
  qed
  ultimately show 
      "Re lb  Re ub  Re lb  Re x  Re x  Re ub" 
      "Im lb  Im ub  Im lb  Im x  Im x  Im ub" 
    by auto
qed
  
lemma closed_segment_degen_complex:
  "Re lb = Re ub; Im lb  Im ub  
     x  closed_segment lb ub  Re x = Re lb  Im lb  Im x  Im x  Im ub "
  "Im lb = Im ub; Re lb  Re ub  
     x  closed_segment lb ub  Im x = Im lb  Re lb  Re x  Re x  Re ub "  
proof -
  show "x  closed_segment lb ub  Re x = Re lb  Im lb  Im x  Im x  Im ub"
    when "Re lb = Re ub" "Im lb  Im ub"
  proof 
    show "Re x = Re lb  Im lb  Im x  Im x  Im ub" when "x  closed_segment lb ub"
      using closed_segment_imp_Re_Im[OF that] ‹Re lb = Re ub ‹Im lb  Im ub by fastforce
  next
    assume asm:"Re x = Re lb  Im lb  Im x  Im x  Im ub"
    define u where "u=(Im x - Im lb)/ (Im ub - Im lb)"
    have "x = (1 - u) *R lb + u *R ub"
      unfolding u_def using asm ‹Re lb = Re ub ‹Im lb  Im ub
      apply (intro complex_eqI)
      apply (auto simp add:field_simps)
      apply (cases "Im ub - Im lb =0")
      apply (auto simp add:field_simps)
      done
    moreover have "0u" "u1" unfolding u_def 
      using ‹Im lb  Im ub asm
      by (cases "Im ub - Im lb =0",auto simp add:field_simps)+
    ultimately show "x  closed_segment lb ub" unfolding closed_segment_def by auto
  qed
  show "x  closed_segment lb ub  Im x = Im lb  Re lb  Re x  Re x  Re ub"
    when "Im lb = Im ub" "Re lb  Re ub"
  proof 
    show "Im x = Im lb  Re lb  Re x  Re x  Re ub" when "x  closed_segment lb ub"
      using closed_segment_imp_Re_Im[OF that] ‹Im lb = Im ub ‹Re lb  Re ub by fastforce
  next
    assume asm:"Im x = Im lb  Re lb  Re x  Re x  Re ub"
    define u where "u=(Re x - Re lb)/ (Re ub - Re lb)"
    have "x = (1 - u) *R lb + u *R ub"
      unfolding u_def using asm ‹Im lb = Im ub ‹Re lb  Re ub
      apply (intro complex_eqI)
       apply (auto simp add:field_simps)
      apply (cases "Re ub - Re lb =0")
       apply (auto simp add:field_simps)
      done
    moreover have "0u" "u1" unfolding u_def 
      using ‹Re lb  Re ub asm
      by (cases "Re ub - Re lb =0",auto simp add:field_simps)+
    ultimately show "x  closed_segment lb ub" unfolding closed_segment_def by auto
  qed   
qed  
   
lemma complex_box_ne_empty: 
  fixes a b::complex
  shows 
    "cbox a b  {}  (Re a  Re b  Im a  Im b)"
    "box a b  {}  (Re a < Re b  Im a < Im b)"
  by (auto simp add:box_ne_empty Basis_complex_def)  
      
lemma proots_rectangle_code1:
  "proots_rectangle p lb ub = (if Re lb < Re ub  Im lb < Im ub then 
            if p0 then 
            if no_proots_line p lb (Complex (Re ub) (Im lb))
             no_proots_line p (Complex (Re ub) (Im lb)) ub
             no_proots_line p ub (Complex (Re lb) (Im ub))
             no_proots_line p (Complex (Re lb) (Im ub)) lb then  
            (
            let p1 = pcompose p [:lb,  Complex (Re ub - Re lb) 0:];
                pR1 = map_poly Re p1; pI1 = map_poly Im p1; gc1 = gcd pR1 pI1;
                p2 = pcompose p [:Complex (Re ub) (Im lb), Complex 0 (Im ub - Im lb):];
                pR2 = map_poly Re p2; pI2 = map_poly Im p2; gc2 = gcd pR2 pI2;
                p3 = pcompose p [:ub, Complex (Re lb - Re ub) 0:];
                pR3 = map_poly Re p3; pI3 = map_poly Im p3; gc3 = gcd pR3 pI3;
                p4 = pcompose p [:Complex (Re lb) (Im ub), Complex 0 (Im lb - Im ub):];
                pR4 = map_poly Re p4; pI4 = map_poly Im p4; gc4 = gcd pR4 pI4
            in 
              nat (- (changes_alt_itv_smods 0 1 (pR1 div gc1) (pI1 div gc1)
                + changes_alt_itv_smods 0 1 (pR2 div gc2) (pI2 div gc2)
                + changes_alt_itv_smods 0 1 (pR3 div gc3) (pI3 div gc3)
                + changes_alt_itv_smods 0 1 (pR4 div gc4) (pI4 div gc4))  div 4)
            )
            else Code.abort (STR ''proots_rectangle fails when there is a root on the border.'') 
            (λ_. proots_rectangle p lb ub)
            else Code.abort (STR ''proots_rectangle fails when p=0.'') 
            (λ_. proots_rectangle p lb ub)
            else 0)"
proof -
  have ?thesis when "¬ (Re lb < Re ub  Im lb < Im ub)" 
  proof -
    have "box lb ub = {}" using complex_box_ne_empty[of lb ub] that by auto
    then have "proots_rectangle p lb ub = 0" unfolding proots_rectangle_def by auto
    then show ?thesis by (simp add:that) 
  qed
  moreover have ?thesis when "Re lb < Re ub  Im lb < Im ub" "p=0"
    using that by simp
  moreover have ?thesis when     
            "Re lb < Re ub" "Im lb < Im ub" "p0" 
            and no_proots:
            "no_proots_line p lb (Complex (Re ub) (Im lb))"
            "no_proots_line p (Complex (Re ub) (Im lb)) ub"
            "no_proots_line p ub (Complex (Re lb) (Im ub))"
            "no_proots_line p (Complex (Re lb) (Im ub)) lb"
  proof -
    define l1 where "l1 = linepath lb (Complex (Re ub) (Im lb))"
    define l2 where "l2 = linepath (Complex (Re ub) (Im lb)) ub"
    define l3 where "l3 = linepath ub (Complex (Re lb) (Im ub))"
    define l4 where "l4 = linepath (Complex (Re lb) (Im ub)) lb"
    define rec where "rec = l1 +++ l2 +++ l3 +++ l4"
    have valid[simp]:"valid_path rec" and loop[simp]:"pathfinish rec = pathstart rec"
      unfolding rec_def l1_def l2_def l3_def l4_def by auto
    have path_no_proots:"path_image rec  proots p = {}" 
      unfolding rec_def l1_def l2_def l3_def l4_def
      apply (subst path_image_join,simp_all del:Complex_eq)+  
      using no_proots[unfolded no_proots_line_def] by (auto simp del:Complex_eq)
    
    define g1 where "g1 = poly p o l1" 
    define g2 where "g2 = poly p o l2" 
    define g3 where "g3 = poly p o l3" 
    define g4 where "g4 = poly p o l4"    
    have [simp]: "path g1" "path g2" "path g3" "path g4"
      "pathfinish g1 = pathstart g2" "pathfinish g2 = pathstart g3" "pathfinish g3 = pathstart g4"
      "pathfinish g4 = pathstart g1"
      unfolding g1_def g2_def g3_def g4_def l1_def l2_def l3_def l4_def
      by (auto intro!: path_continuous_image continuous_intros 
          simp add:pathfinish_compose pathstart_compose)
    have [simp]: "finite_ReZ_segments g1 0" "finite_ReZ_segments g2 0" 
         "finite_ReZ_segments g3 0" "finite_ReZ_segments g4 0"
      unfolding g1_def l1_def g2_def l2_def g3_def l3_def g4_def l4_def poly_linepath_comp 
      by (rule finite_ReZ_segments_poly_of_real)+
    define p1 pR1 pI1 gc1
           p2 pR2 pI2 gc2
           p3 pR3 pI3 gc3
           p4 pR4 pI4 gc4
      where "p1 = pcompose p [:lb,  Complex (Re ub - Re lb) 0:]"
             and "pR1 = map_poly Re p1" and "pI1 = map_poly Im p1" and "gc1 = gcd pR1 pI1"
             and "p2 = pcompose p [:Complex (Re ub) (Im lb), Complex 0 (Im ub - Im lb):]"
             and "pR2 = map_poly Re p2" and "pI2 = map_poly Im p2" and "gc2 = gcd pR2 pI2"
             and "p3 = pcompose p [:ub, Complex (Re lb - Re ub) 0:]"
             and "pR3 = map_poly Re p3" and "pI3 = map_poly Im p3" and "gc3 = gcd pR3 pI3"
             and "p4 = pcompose p [:Complex (Re lb) (Im ub), Complex 0 (Im lb - Im ub):]"
             and "pR4 = map_poly Re p4" and "pI4 = map_poly Im p4" and "gc4 = gcd pR4 pI4"
    have "gc10" "gc20" "gc30" "gc40"
    proof -
      show "gc10" 
      proof (rule ccontr)
        assume "¬ gc1  0"
        then have "pI1 = 0" "pR1 = 0" unfolding gc1_def by auto
        then have "p1 = 0" unfolding pI1_def pR1_def 
          by (metis cpoly_of_decompose map_poly_0)
        then have "p=0" unfolding p1_def using ‹Re lb < Re ub
          by (auto elim!:pcompose_eq_0 simp add:Complex_eq_0)
        then show False using p0 by simp
      qed
      show "gc20" 
      proof (rule ccontr)
        assume "¬ gc2  0"
        then have "pI2 = 0" "pR2 = 0" unfolding gc2_def by auto
        then have "p2 = 0" unfolding pI2_def pR2_def 
          by (metis cpoly_of_decompose map_poly_0)
        then have "p=0" unfolding p2_def using ‹Im lb < Im ub
          by (auto elim!:pcompose_eq_0 simp add:Complex_eq_0)
        then show False using p0 by simp
      qed 
      show "gc30" 
      proof (rule ccontr)
        assume "¬ gc3  0"
        then have "pI3 = 0" "pR3 = 0" unfolding gc3_def by auto
        then have "p3 = 0" unfolding pI3_def pR3_def 
          by (metis cpoly_of_decompose map_poly_0)
        then have "p=0" unfolding p3_def using ‹Re lb < Re ub
          by (auto elim!:pcompose_eq_0 simp add:Complex_eq_0)
        then show False using p0 by simp
      qed
      show "gc40" 
      proof (rule ccontr)
        assume "¬ gc4  0"
        then have "pI4 = 0" "pR4 = 0" unfolding gc4_def by auto
        then have "p4 = 0" unfolding pI4_def pR4_def 
          by (metis cpoly_of_decompose map_poly_0)
        then have "p=0" unfolding p4_def using ‹Im lb < Im ub
          by (auto elim!:pcompose_eq_0 simp add:Complex_eq_0)
        then show False using p0 by simp
      qed 
    qed
    define sms where 
      "sms = (changes_alt_itv_smods 0 1 (pR1 div gc1) (pI1 div gc1)
                + changes_alt_itv_smods 0 1 (pR2 div gc2) (pI2 div gc2)
                + changes_alt_itv_smods 0 1 (pR3 div gc3) (pI3 div gc3)
                + changes_alt_itv_smods 0 1 (pR4 div gc4) (pI4 div gc4))"
    have "proots_rectangle p lb ub = (rproots p. winding_number rec r * (order r p))" 
    proof -
      have "winding_number rec x * of_nat (order x p) = 0" 
        when "xproots p - proots_within p (box lb ub)" for x 
      proof -
        have *:"cbox lb ub = box lb ub  path_image rec"
        proof -
          have "xcbox lb ub" when "xbox lb ub  path_image rec" for x
            using that ‹Re lb<Re ub ‹Im lb<Im ub
            unfolding box_def cbox_def Basis_complex_def rec_def l1_def l2_def l3_def l4_def
            apply (auto simp add:path_image_join closed_segment_degen_complex)          
                   apply (subst (asm) closed_segment_commute,simp add: closed_segment_degen_complex)+
            done  
          moreover have "xbox lb ub  path_image rec" when "xcbox lb ub" for x
            using that
            unfolding box_def cbox_def Basis_complex_def rec_def l1_def l2_def l3_def l4_def
            apply (auto simp add:path_image_join closed_segment_degen_complex)
             apply (subst (asm) (1 2) closed_segment_commute,simp add:closed_segment_degen_complex)+
            done
          ultimately show ?thesis by auto
        qed
        moreover have "xpath_image rec" 
          using path_no_proots that by auto
        ultimately have "xcbox lb ub" using that by simp
        from winding_number_zero_outside[OF valid_path_imp_path[OF valid] _ loop this,simplified] * 
        have "winding_number rec x = 0" by auto
        then show ?thesis by auto
      qed
      moreover have "of_nat (order x p) = winding_number rec x * of_nat (order x p)" when 
        "x  proots_within p (box lb ub)" for x 
      proof -
        have "xbox lb ub" using that unfolding proots_within_def by auto 
        then have order_asms:"Re lb < Re x" "Re x < Re ub" "Im lb < Im x" "Im x < Im ub"  
          by (auto simp add:box_def Basis_complex_def)  
        have "winding_number rec x = 1"
          unfolding rec_def l1_def l2_def l3_def l4_def 
        proof eval_winding 
          let ?l1 = "linepath lb (Complex (Re ub) (Im lb))"
          and ?l2 = "linepath (Complex (Re ub) (Im lb)) ub"
          and ?l3 = "linepath ub (Complex (Re lb) (Im ub))"
          and ?l4 = "linepath (Complex (Re lb) (Im ub)) lb"
          show l1: "x  path_image ?l1" and l2: "x  path_image ?l2" and 
               l3: "x  path_image ?l3" and l4:"x  path_image ?l4"
            using no_proots that unfolding no_proots_line_def by auto
          show "- of_real (cindex_pathE ?l1 x + (cindex_pathE ?l2 x +
            (cindex_pathE ?l3 x + cindex_pathE ?l4 x))) = 2 * 1"
          proof -
            have "(Im x - Im ub) * (Re ub - Re lb) < 0" 
              using mult_less_0_iff order_asms(1) order_asms(2) order_asms(4) by fastforce
            then have "cindex_pathE ?l3 x = -1" 
              apply (subst cindex_pathE_linepath)
              using l3 by (auto simp add:algebra_simps order_asms)
            moreover have "(Im lb - Im x) * (Re ub - Re lb) <0" 
              using mult_less_0_iff order_asms(1) order_asms(2) order_asms(3) by fastforce
            then have "cindex_pathE ?l1 x = -1"    
              apply (subst cindex_pathE_linepath)
              using l1 by (auto simp add:algebra_simps order_asms)  
            moreover have "cindex_pathE ?l2 x = 0"    
              apply (subst cindex_pathE_linepath)
              using l2 order_asms by auto
            moreover have "cindex_pathE ?l4 x = 0"
              apply (subst cindex_pathE_linepath)
              using l4 order_asms by auto
            ultimately show ?thesis by auto
          qed
        qed
        then show ?thesis by auto
      qed 
      ultimately show ?thesis using p0
      unfolding proots_rectangle_def proots_count_def 
        by (auto intro!: sum.mono_neutral_cong_left[of "proots p" "proots_within p (box lb ub)"])
    qed    
    also have "... = 1/(2 * of_real pi * 𝗂) *contour_integral rec (λx. deriv (poly p) x / poly p x)"
    proof -
      have "contour_integral rec (λx. deriv (poly p) x / poly p x) = 2 * of_real pi * 𝗂 
              * (x | poly p x = 0. winding_number rec x * of_int (zorder (poly p) x))"
      proof (rule argument_principle[of UNIV "poly p" "{}" "λ_. 1" rec,simplified])
        show "connected (UNIV::complex set)" using connected_UNIV[where 'a=complex] .
        show "path_image rec  UNIV - {x. poly p x = 0}" 
          using path_no_proots unfolding proots_within_def by auto
        show "finite {x. poly p x = 0}" by (simp add: poly_roots_finite that(3))
      qed
      also have " ... = 2 * of_real pi * 𝗂 * (xproots p. winding_number rec x * (order x p))" 
        unfolding proots_within_def 
        apply (auto intro!:sum.cong simp add:order_zorder[OF p0] )
        by (metis nat_eq_iff2 of_nat_nat order_root order_zorder that(3))
      finally show ?thesis by auto
    qed  
    also have "... = winding_number (poly p  rec) 0"
    proof -
      have "0  path_image (poly p  rec)" 
        using path_no_proots unfolding path_image_compose proots_within_def by fastforce
      from winding_number_comp[OF _ poly_holomorphic_on _ _ this,of UNIV,simplified]
      show ?thesis by auto
    qed
    also have winding_eq:"... = - cindex_pathE (poly p  rec) 0 / 2"
    proof (rule winding_number_cindex_pathE)
      show "finite_ReZ_segments (poly p  rec) 0"
        unfolding rec_def path_compose_join 
        apply (fold g1_def g2_def g3_def g4_def)
        by (auto intro!: finite_ReZ_segments_joinpaths path_join_imp)
      show "valid_path (poly p  rec)"
        by (rule valid_path_compose_holomorphic[where S=UNIV]) auto
      show "0  path_image (poly p  rec)"
        using path_no_proots unfolding path_image_compose proots_def by fastforce
      show "pathfinish (poly p  rec) = pathstart (poly p  rec)"
        unfolding rec_def pathstart_compose pathfinish_compose  by (auto simp add:l1_def l4_def)   
    qed
    also have cindex_pathE_eq:"... = of_int (- sms) / of_int 4"
    proof -
      have "cindex_pathE (poly p  rec) 0 = cindex_pathE (g1+++g2+++g3+++g4) 0"
        unfolding rec_def path_compose_join g1_def g2_def g3_def g4_def by simp
      also have "... = cindex_pathE g1 0 + cindex_pathE g2 0 + cindex_pathE g3 0