src/HOL/Library/Fundamental_Theorem_Algebra.thy
author wenzelm
Tue Apr 29 21:54:26 2014 +0200 (2014-04-29)
changeset 56795 e8cce2bd23e5
parent 56778 cb0929421ca6
child 56889 48a745e1bde7
permissions -rw-r--r--
tuned proofs;
     1 (* Author: Amine Chaieb, TU Muenchen *)
     2 
     3 header{*Fundamental Theorem of Algebra*}
     4 
     5 theory Fundamental_Theorem_Algebra
     6 imports Polynomial Complex_Main
     7 begin
     8 
     9 subsection {* Square root of complex numbers *}
    10 
    11 definition csqrt :: "complex \<Rightarrow> complex"
    12 where
    13   "csqrt z =
    14     (if Im z = 0 then
    15        if 0 \<le> Re z then Complex (sqrt(Re z)) 0
    16        else Complex 0 (sqrt(- Re z))
    17      else Complex (sqrt((cmod z + Re z) /2)) ((Im z / abs(Im z)) * sqrt((cmod z - Re z) /2)))"
    18 
    19 lemma csqrt[algebra]: "(csqrt z)\<^sup>2 = z"
    20 proof -
    21   obtain x y where xy: "z = Complex x y" by (cases z)
    22   {
    23     assume y0: "y = 0"
    24     {
    25       assume x0: "x \<ge> 0"
    26       then have ?thesis
    27         using y0 xy real_sqrt_pow2[OF x0]
    28         by (simp add: csqrt_def power2_eq_square)
    29     }
    30     moreover
    31     {
    32       assume "\<not> x \<ge> 0"
    33       then have x0: "- x \<ge> 0" by arith
    34       then have ?thesis
    35         using y0 xy real_sqrt_pow2[OF x0]
    36         by (simp add: csqrt_def power2_eq_square)
    37     }
    38     ultimately have ?thesis by blast
    39   }
    40   moreover
    41   {
    42     assume y0: "y \<noteq> 0"
    43     {
    44       fix x y
    45       let ?z = "Complex x y"
    46       from abs_Re_le_cmod[of ?z] have tha: "abs x \<le> cmod ?z"
    47         by auto
    48       then have "cmod ?z - x \<ge> 0" "cmod ?z + x \<ge> 0"
    49         by arith+
    50       then have "(sqrt (x * x + y * y) + x) / 2 \<ge> 0" "(sqrt (x * x + y * y) - x) / 2 \<ge> 0"
    51         by (simp_all add: power2_eq_square)
    52     }
    53     note th = this
    54     have sq4: "\<And>x::real. x\<^sup>2 / 4 = (x / 2)\<^sup>2"
    55       by (simp add: power2_eq_square)
    56     from th[of x y]
    57     have sq4': "sqrt (((sqrt (x * x + y * y) + x)\<^sup>2 / 4)) = (sqrt (x * x + y * y) + x) / 2"
    58       "sqrt (((sqrt (x * x + y * y) - x)\<^sup>2 / 4)) = (sqrt (x * x + y * y) - x) / 2"
    59       unfolding sq4 by simp_all
    60     then have th1: "sqrt ((sqrt (x * x + y * y) + x) * (sqrt (x * x + y * y) + x) / 4) -
    61         sqrt ((sqrt (x * x + y * y) - x) * (sqrt (x * x + y * y) - x) / 4) = x"
    62       unfolding power2_eq_square by simp
    63     have "sqrt 4 = sqrt (2\<^sup>2)"
    64       by simp
    65     then have sqrt4: "sqrt 4 = 2"
    66       by (simp only: real_sqrt_abs)
    67     have th2: "2 * (y * sqrt ((sqrt (x * x + y * y) - x) * (sqrt (x * x + y * y) + x) / 4)) / \<bar>y\<bar> = y"
    68       using iffD2[OF real_sqrt_pow2_iff sum_power2_ge_zero[of x y]] y0
    69       unfolding power2_eq_square
    70       by (simp add: algebra_simps real_sqrt_divide sqrt4)
    71     from y0 xy have ?thesis
    72       apply (simp add: csqrt_def power2_eq_square)
    73       apply (simp add: real_sqrt_sum_squares_mult_ge_zero[of x y]
    74         real_sqrt_pow2[OF th(1)[of x y], unfolded power2_eq_square]
    75         real_sqrt_pow2[OF th(2)[of x y], unfolded power2_eq_square]
    76         real_sqrt_mult[symmetric])
    77       using th1 th2  ..
    78   }
    79   ultimately show ?thesis by blast
    80 qed
    81 
    82 lemma csqrt_Complex: "x \<ge> 0 \<Longrightarrow> csqrt (Complex x 0) = Complex (sqrt x) 0"
    83   by (simp add: csqrt_def)
    84 
    85 lemma csqrt_0 [simp]: "csqrt 0 = 0"
    86   by (simp add: csqrt_def)
    87 
    88 lemma csqrt_1 [simp]: "csqrt 1 = 1"
    89   by (simp add: csqrt_def)
    90 
    91 lemma csqrt_principal: "0 < Re(csqrt(z)) | Re(csqrt(z)) = 0 & 0 \<le> Im(csqrt(z))"
    92 proof (cases z)
    93   case (Complex x y)
    94   then show ?thesis
    95     using real_sqrt_sum_squares_ge1 [of "x" y]
    96           real_sqrt_sum_squares_ge1 [of "-x" y]
    97           real_sqrt_sum_squares_eq_cancel [of x y]
    98     apply (auto simp: csqrt_def intro!: Rings.ordered_ring_class.split_mult_pos_le)
    99     apply (metis add_commute diff_add_cancel le_add_same_cancel1 real_sqrt_sum_squares_ge1)
   100     apply (metis add_commute less_eq_real_def power_minus_Bit0
   101             real_0_less_add_iff real_sqrt_sum_squares_eq_cancel)
   102     done
   103 qed
   104 
   105 lemma Re_csqrt: "0 \<le> Re(csqrt z)"
   106   by (metis csqrt_principal le_less)
   107 
   108 lemma csqrt_square: "0 < Re z \<or> Re z = 0 \<and> 0 \<le> Im z \<Longrightarrow> csqrt (z\<^sup>2) = z"
   109   using csqrt [of "z\<^sup>2"] csqrt_principal [of "z\<^sup>2"]
   110   by (cases z) (auto simp: power2_eq_iff)
   111 
   112 lemma csqrt_eq_0 [simp]: "csqrt z = 0 \<longleftrightarrow> z = 0"
   113   by auto (metis csqrt power_eq_0_iff)
   114 
   115 lemma csqrt_eq_1 [simp]: "csqrt z = 1 \<longleftrightarrow> z = 1"
   116   by auto (metis csqrt power2_eq_1_iff)
   117 
   118 
   119 subsection {* More lemmas about module of complex numbers *}
   120 
   121 lemma complex_of_real_power: "complex_of_real x ^ n = complex_of_real (x^n)"
   122   by (rule of_real_power [symmetric])
   123 
   124 text{* The triangle inequality for cmod *}
   125 lemma complex_mod_triangle_sub: "cmod w \<le> cmod (w + z) + norm z"
   126   using complex_mod_triangle_ineq2[of "w + z" "-z"] by auto
   127 
   128 
   129 subsection {* Basic lemmas about polynomials *}
   130 
   131 lemma poly_bound_exists:
   132   fixes p :: "'a::{comm_semiring_0,real_normed_div_algebra} poly"
   133   shows "\<exists>m. m > 0 \<and> (\<forall>z. norm z \<le> r \<longrightarrow> norm (poly p z) \<le> m)"
   134 proof (induct p)
   135   case 0
   136   then show ?case by (rule exI[where x=1]) simp
   137 next
   138   case (pCons c cs)
   139   from pCons.hyps obtain m where m: "\<forall>z. norm z \<le> r \<longrightarrow> norm (poly cs z) \<le> m"
   140     by blast
   141   let ?k = " 1 + norm c + \<bar>r * m\<bar>"
   142   have kp: "?k > 0"
   143     using abs_ge_zero[of "r*m"] norm_ge_zero[of c] by arith
   144   {
   145     fix z :: 'a
   146     assume H: "norm z \<le> r"
   147     from m H have th: "norm (poly cs z) \<le> m"
   148       by blast
   149     from H have rp: "r \<ge> 0"
   150       using norm_ge_zero[of z] by arith
   151     have "norm (poly (pCons c cs) z) \<le> norm c + norm (z * poly cs z)"
   152       using norm_triangle_ineq[of c "z* poly cs z"] by simp
   153     also have "\<dots> \<le> norm c + r * m"
   154       using mult_mono[OF H th rp norm_ge_zero[of "poly cs z"]]
   155       by (simp add: norm_mult)
   156     also have "\<dots> \<le> ?k"
   157       by simp
   158     finally have "norm (poly (pCons c cs) z) \<le> ?k" .
   159   }
   160   with kp show ?case by blast
   161 qed
   162 
   163 
   164 text{* Offsetting the variable in a polynomial gives another of same degree *}
   165 
   166 definition offset_poly :: "'a::comm_semiring_0 poly \<Rightarrow> 'a \<Rightarrow> 'a poly"
   167   where "offset_poly p h = fold_coeffs (\<lambda>a q. smult h q + pCons a q) p 0"
   168 
   169 lemma offset_poly_0: "offset_poly 0 h = 0"
   170   by (simp add: offset_poly_def)
   171 
   172 lemma offset_poly_pCons:
   173   "offset_poly (pCons a p) h =
   174     smult h (offset_poly p h) + pCons a (offset_poly p h)"
   175   by (cases "p = 0 \<and> a = 0") (auto simp add: offset_poly_def)
   176 
   177 lemma offset_poly_single: "offset_poly [:a:] h = [:a:]"
   178   by (simp add: offset_poly_pCons offset_poly_0)
   179 
   180 lemma poly_offset_poly: "poly (offset_poly p h) x = poly p (h + x)"
   181   apply (induct p)
   182   apply (simp add: offset_poly_0)
   183   apply (simp add: offset_poly_pCons algebra_simps)
   184   done
   185 
   186 lemma offset_poly_eq_0_lemma: "smult c p + pCons a p = 0 \<Longrightarrow> p = 0"
   187   by (induct p arbitrary: a) (simp, force)
   188 
   189 lemma offset_poly_eq_0_iff: "offset_poly p h = 0 \<longleftrightarrow> p = 0"
   190   apply (safe intro!: offset_poly_0)
   191   apply (induct p)
   192   apply simp
   193   apply (simp add: offset_poly_pCons)
   194   apply (frule offset_poly_eq_0_lemma, simp)
   195   done
   196 
   197 lemma degree_offset_poly: "degree (offset_poly p h) = degree p"
   198   apply (induct p)
   199   apply (simp add: offset_poly_0)
   200   apply (case_tac "p = 0")
   201   apply (simp add: offset_poly_0 offset_poly_pCons)
   202   apply (simp add: offset_poly_pCons)
   203   apply (subst degree_add_eq_right)
   204   apply (rule le_less_trans [OF degree_smult_le])
   205   apply (simp add: offset_poly_eq_0_iff)
   206   apply (simp add: offset_poly_eq_0_iff)
   207   done
   208 
   209 definition "psize p = (if p = 0 then 0 else Suc (degree p))"
   210 
   211 lemma psize_eq_0_iff [simp]: "psize p = 0 \<longleftrightarrow> p = 0"
   212   unfolding psize_def by simp
   213 
   214 lemma poly_offset:
   215   fixes p :: "'a::comm_ring_1 poly"
   216   shows "\<exists>q. psize q = psize p \<and> (\<forall>x. poly q x = poly p (a + x))"
   217 proof (intro exI conjI)
   218   show "psize (offset_poly p a) = psize p"
   219     unfolding psize_def
   220     by (simp add: offset_poly_eq_0_iff degree_offset_poly)
   221   show "\<forall>x. poly (offset_poly p a) x = poly p (a + x)"
   222     by (simp add: poly_offset_poly)
   223 qed
   224 
   225 text{* An alternative useful formulation of completeness of the reals *}
   226 lemma real_sup_exists:
   227   assumes ex: "\<exists>x. P x"
   228     and bz: "\<exists>z. \<forall>x. P x \<longrightarrow> x < z"
   229   shows "\<exists>s::real. \<forall>y. (\<exists>x. P x \<and> y < x) \<longleftrightarrow> y < s"
   230 proof
   231   from bz have "bdd_above (Collect P)"
   232     by (force intro: less_imp_le)
   233   then show "\<forall>y. (\<exists>x. P x \<and> y < x) \<longleftrightarrow> y < Sup (Collect P)"
   234     using ex bz by (subst less_cSup_iff) auto
   235 qed
   236 
   237 subsection {* Fundamental theorem of algebra *}
   238 lemma  unimodular_reduce_norm:
   239   assumes md: "cmod z = 1"
   240   shows "cmod (z + 1) < 1 \<or> cmod (z - 1) < 1 \<or> cmod (z + ii) < 1 \<or> cmod (z - ii) < 1"
   241 proof -
   242   obtain x y where z: "z = Complex x y "
   243     by (cases z) auto
   244   from md z have xy: "x\<^sup>2 + y\<^sup>2 = 1"
   245     by (simp add: cmod_def)
   246   {
   247     assume C: "cmod (z + 1) \<ge> 1" "cmod (z - 1) \<ge> 1" "cmod (z + ii) \<ge> 1" "cmod (z - ii) \<ge> 1"
   248     from C z xy have "2 * x \<le> 1" "2 * x \<ge> -1" "2 * y \<le> 1" "2 * y \<ge> -1"
   249       by (simp_all add: cmod_def power2_eq_square algebra_simps)
   250     then have "abs (2 * x) \<le> 1" "abs (2 * y) \<le> 1"
   251       by simp_all
   252     then have "(abs (2 * x))\<^sup>2 \<le> 1\<^sup>2" "(abs (2 * y))\<^sup>2 \<le> 1\<^sup>2"
   253       by - (rule power_mono, simp, simp)+
   254     then have th0: "4 * x\<^sup>2 \<le> 1" "4 * y\<^sup>2 \<le> 1"
   255       by (simp_all add: power_mult_distrib)
   256     from add_mono[OF th0] xy have False by simp
   257   }
   258   then show ?thesis
   259     unfolding linorder_not_le[symmetric] by blast
   260 qed
   261 
   262 text{* Hence we can always reduce modulus of @{text "1 + b z^n"} if nonzero *}
   263 lemma reduce_poly_simple:
   264   assumes b: "b \<noteq> 0"
   265     and n: "n \<noteq> 0"
   266   shows "\<exists>z. cmod (1 + b * z^n) < 1"
   267   using n
   268 proof (induct n rule: nat_less_induct)
   269   fix n
   270   assume IH: "\<forall>m<n. m \<noteq> 0 \<longrightarrow> (\<exists>z. cmod (1 + b * z ^ m) < 1)"
   271   assume n: "n \<noteq> 0"
   272   let ?P = "\<lambda>z n. cmod (1 + b * z ^ n) < 1"
   273   {
   274     assume e: "even n"
   275     then have "\<exists>m. n = 2 * m"
   276       by presburger
   277     then obtain m where m: "n = 2 * m"
   278       by blast
   279     from n m have "m \<noteq> 0" "m < n"
   280       by presburger+
   281     with IH[rule_format, of m] obtain z where z: "?P z m"
   282       by blast
   283     from z have "?P (csqrt z) n"
   284       by (simp add: m power_mult csqrt)
   285     then have "\<exists>z. ?P z n" ..
   286   }
   287   moreover
   288   {
   289     assume o: "odd n"
   290     have th0: "cmod (complex_of_real (cmod b) / b) = 1"
   291       using b by (simp add: norm_divide)
   292     from o have "\<exists>m. n = Suc (2 * m)"
   293       by presburger+
   294     then obtain m where m: "n = Suc (2 * m)"
   295       by blast
   296     from unimodular_reduce_norm[OF th0] o
   297     have "\<exists>v. cmod (complex_of_real (cmod b) / b + v^n) < 1"
   298       apply (cases "cmod (complex_of_real (cmod b) / b + 1) < 1")
   299       apply (rule_tac x="1" in exI)
   300       apply simp
   301       apply (cases "cmod (complex_of_real (cmod b) / b - 1) < 1")
   302       apply (rule_tac x="-1" in exI)
   303       apply simp
   304       apply (cases "cmod (complex_of_real (cmod b) / b + ii) < 1")
   305       apply (cases "even m")
   306       apply (rule_tac x="ii" in exI)
   307       apply (simp add: m power_mult)
   308       apply (rule_tac x="- ii" in exI)
   309       apply (simp add: m power_mult)
   310       apply (cases "even m")
   311       apply (rule_tac x="- ii" in exI)
   312       apply (simp add: m power_mult)
   313       apply (auto simp add: m power_mult)
   314       apply (rule_tac x="ii" in exI)
   315       apply (auto simp add: m power_mult)
   316       done
   317     then obtain v where v: "cmod (complex_of_real (cmod b) / b + v^n) < 1"
   318       by blast
   319     let ?w = "v / complex_of_real (root n (cmod b))"
   320     from odd_real_root_pow[OF o, of "cmod b"]
   321     have th1: "?w ^ n = v^n / complex_of_real (cmod b)"
   322       by (simp add: power_divide complex_of_real_power)
   323     have th2:"cmod (complex_of_real (cmod b) / b) = 1"
   324       using b by (simp add: norm_divide)
   325     then have th3: "cmod (complex_of_real (cmod b) / b) \<ge> 0"
   326       by simp
   327     have th4: "cmod (complex_of_real (cmod b) / b) *
   328         cmod (1 + b * (v ^ n / complex_of_real (cmod b))) <
   329         cmod (complex_of_real (cmod b) / b) * 1"
   330       apply (simp only: norm_mult[symmetric] distrib_left)
   331       using b v
   332       apply (simp add: th2)
   333       done
   334     from mult_less_imp_less_left[OF th4 th3]
   335     have "?P ?w n" unfolding th1 .
   336     then have "\<exists>z. ?P z n" ..
   337   }
   338   ultimately show "\<exists>z. ?P z n" by blast
   339 qed
   340 
   341 text{* Bolzano-Weierstrass type property for closed disc in complex plane. *}
   342 
   343 lemma metric_bound_lemma: "cmod (x - y) \<le> \<bar>Re x - Re y\<bar> + \<bar>Im x - Im y\<bar>"
   344   using real_sqrt_sum_squares_triangle_ineq[of "Re x - Re y" 0 0 "Im x - Im y"]
   345   unfolding cmod_def by simp
   346 
   347 lemma bolzano_weierstrass_complex_disc:
   348   assumes r: "\<forall>n. cmod (s n) \<le> r"
   349   shows "\<exists>f z. subseq f \<and> (\<forall>e >0. \<exists>N. \<forall>n \<ge> N. cmod (s (f n) - z) < e)"
   350 proof-
   351   from seq_monosub[of "Re \<circ> s"]
   352   obtain f where f: "subseq f" "monoseq (\<lambda>n. Re (s (f n)))"
   353     unfolding o_def by blast
   354   from seq_monosub[of "Im \<circ> s \<circ> f"]
   355   obtain g where g: "subseq g" "monoseq (\<lambda>n. Im (s (f (g n))))"
   356     unfolding o_def by blast
   357   let ?h = "f \<circ> g"
   358   from r[rule_format, of 0] have rp: "r \<ge> 0"
   359     using norm_ge_zero[of "s 0"] by arith
   360   have th: "\<forall>n. r + 1 \<ge> \<bar>Re (s n)\<bar>"
   361   proof
   362     fix n
   363     from abs_Re_le_cmod[of "s n"] r[rule_format, of n]
   364     show "\<bar>Re (s n)\<bar> \<le> r + 1" by arith
   365   qed
   366   have conv1: "convergent (\<lambda>n. Re (s (f n)))"
   367     apply (rule Bseq_monoseq_convergent)
   368     apply (simp add: Bseq_def)
   369     apply (metis gt_ex le_less_linear less_trans order.trans th)
   370     apply (rule f(2))
   371     done
   372   have th: "\<forall>n. r + 1 \<ge> \<bar>Im (s n)\<bar>"
   373   proof
   374     fix n
   375     from abs_Im_le_cmod[of "s n"] r[rule_format, of n]
   376     show "\<bar>Im (s n)\<bar> \<le> r + 1"
   377       by arith
   378   qed
   379 
   380   have conv2: "convergent (\<lambda>n. Im (s (f (g n))))"
   381     apply (rule Bseq_monoseq_convergent)
   382     apply (simp add: Bseq_def)
   383     apply (metis gt_ex le_less_linear less_trans order.trans th)
   384     apply (rule g(2))
   385     done
   386 
   387   from conv1[unfolded convergent_def] obtain x where "LIMSEQ (\<lambda>n. Re (s (f n))) x"
   388     by blast
   389   then have x: "\<forall>r>0. \<exists>n0. \<forall>n\<ge>n0. \<bar>Re (s (f n)) - x\<bar> < r"
   390     unfolding LIMSEQ_iff real_norm_def .
   391 
   392   from conv2[unfolded convergent_def] obtain y where "LIMSEQ (\<lambda>n. Im (s (f (g n)))) y"
   393     by blast
   394   then have y: "\<forall>r>0. \<exists>n0. \<forall>n\<ge>n0. \<bar>Im (s (f (g n))) - y\<bar> < r"
   395     unfolding LIMSEQ_iff real_norm_def .
   396   let ?w = "Complex x y"
   397   from f(1) g(1) have hs: "subseq ?h"
   398     unfolding subseq_def by auto
   399   {
   400     fix e :: real
   401     assume ep: "e > 0"
   402     then have e2: "e/2 > 0"
   403       by simp
   404     from x[rule_format, OF e2] y[rule_format, OF e2]
   405     obtain N1 N2 where N1: "\<forall>n\<ge>N1. \<bar>Re (s (f n)) - x\<bar> < e / 2"
   406       and N2: "\<forall>n\<ge>N2. \<bar>Im (s (f (g n))) - y\<bar> < e / 2"
   407       by blast
   408     {
   409       fix n
   410       assume nN12: "n \<ge> N1 + N2"
   411       then have nN1: "g n \<ge> N1" and nN2: "n \<ge> N2"
   412         using seq_suble[OF g(1), of n] by arith+
   413       from add_strict_mono[OF N1[rule_format, OF nN1] N2[rule_format, OF nN2]]
   414       have "cmod (s (?h n) - ?w) < e"
   415         using metric_bound_lemma[of "s (f (g n))" ?w] by simp
   416     }
   417     then have "\<exists>N. \<forall>n\<ge>N. cmod (s (?h n) - ?w) < e"
   418       by blast
   419   }
   420   with hs show ?thesis by blast
   421 qed
   422 
   423 text{* Polynomial is continuous. *}
   424 
   425 lemma poly_cont:
   426   fixes p :: "'a::{comm_semiring_0,real_normed_div_algebra} poly"
   427   assumes ep: "e > 0"
   428   shows "\<exists>d >0. \<forall>w. 0 < norm (w - z) \<and> norm (w - z) < d \<longrightarrow> norm (poly p w - poly p z) < e"
   429 proof -
   430   obtain q where q: "degree q = degree p" "\<And>x. poly q x = poly p (z + x)"
   431   proof
   432     show "degree (offset_poly p z) = degree p"
   433       by (rule degree_offset_poly)
   434     show "\<And>x. poly (offset_poly p z) x = poly p (z + x)"
   435       by (rule poly_offset_poly)
   436   qed
   437   have th: "\<And>w. poly q (w - z) = poly p w"
   438     using q(2)[of "w - z" for w] by simp
   439   show ?thesis unfolding th[symmetric]
   440   proof (induct q)
   441     case 0
   442     then show ?case
   443       using ep by auto
   444   next
   445     case (pCons c cs)
   446     from poly_bound_exists[of 1 "cs"]
   447     obtain m where m: "m > 0" "\<And>z. norm z \<le> 1 \<Longrightarrow> norm (poly cs z) \<le> m"
   448       by blast
   449     from ep m(1) have em0: "e/m > 0"
   450       by (simp add: field_simps)
   451     have one0: "1 > (0::real)"
   452       by arith
   453     from real_lbound_gt_zero[OF one0 em0]
   454     obtain d where d: "d > 0" "d < 1" "d < e / m"
   455       by blast
   456     from d(1,3) m(1) have dm: "d * m > 0" "d * m < e"
   457       by (simp_all add: field_simps)
   458     show ?case
   459     proof (rule ex_forward[OF real_lbound_gt_zero[OF one0 em0]], clarsimp simp add: norm_mult)
   460       fix d w
   461       assume H: "d > 0" "d < 1" "d < e/m" "w \<noteq> z" "norm (w - z) < d"
   462       then have d1: "norm (w-z) \<le> 1" "d \<ge> 0"
   463         by simp_all
   464       from H(3) m(1) have dme: "d*m < e"
   465         by (simp add: field_simps)
   466       from H have th: "norm (w - z) \<le> d"
   467         by simp
   468       from mult_mono[OF th m(2)[OF d1(1)] d1(2) norm_ge_zero] dme
   469       show "norm (w - z) * norm (poly cs (w - z)) < e"
   470         by simp
   471     qed
   472   qed
   473 qed
   474 
   475 text{* Hence a polynomial attains minimum on a closed disc
   476   in the complex plane. *}
   477 lemma poly_minimum_modulus_disc: "\<exists>z. \<forall>w. cmod w \<le> r \<longrightarrow> cmod (poly p z) \<le> cmod (poly p w)"
   478 proof -
   479   {
   480     assume "\<not> r \<ge> 0"
   481     then have ?thesis
   482       by (metis norm_ge_zero order.trans)
   483   }
   484   moreover
   485   {
   486     assume rp: "r \<ge> 0"
   487     from rp have "cmod 0 \<le> r \<and> cmod (poly p 0) = - (- cmod (poly p 0))"
   488       by simp
   489     then have mth1: "\<exists>x z. cmod z \<le> r \<and> cmod (poly p z) = - x"
   490       by blast
   491     {
   492       fix x z
   493       assume H: "cmod z \<le> r" "cmod (poly p z) = - x" "\<not> x < 1"
   494       then have "- x < 0 "
   495         by arith
   496       with H(2) norm_ge_zero[of "poly p z"] have False
   497         by simp
   498     }
   499     then have mth2: "\<exists>z. \<forall>x. (\<exists>z. cmod z \<le> r \<and> cmod (poly p z) = - x) \<longrightarrow> x < z"
   500       by blast
   501     from real_sup_exists[OF mth1 mth2] obtain s where
   502       s: "\<forall>y. (\<exists>x. (\<exists>z. cmod z \<le> r \<and> cmod (poly p z) = - x) \<and> y < x) \<longleftrightarrow> y < s" by blast
   503     let ?m = "- s"
   504     {
   505       fix y
   506       from s[rule_format, of "-y"]
   507       have "(\<exists>z x. cmod z \<le> r \<and> - (- cmod (poly p z)) < y) \<longleftrightarrow> ?m < y"
   508         unfolding minus_less_iff[of y ] equation_minus_iff by blast
   509     }
   510     note s1 = this[unfolded minus_minus]
   511     from s1[of ?m] have s1m: "\<And>z x. cmod z \<le> r \<Longrightarrow> cmod (poly p z) \<ge> ?m"
   512       by auto
   513     {
   514       fix n :: nat
   515       from s1[rule_format, of "?m + 1/real (Suc n)"]
   516       have "\<exists>z. cmod z \<le> r \<and> cmod (poly p z) < - s + 1 / real (Suc n)"
   517         by simp
   518     }
   519     then have th: "\<forall>n. \<exists>z. cmod z \<le> r \<and> cmod (poly p z) < - s + 1 / real (Suc n)" ..
   520     from choice[OF th] obtain g where
   521         g: "\<forall>n. cmod (g n) \<le> r" "\<forall>n. cmod (poly p (g n)) <?m + 1 /real(Suc n)"
   522       by blast
   523     from bolzano_weierstrass_complex_disc[OF g(1)]
   524     obtain f z where fz: "subseq f" "\<forall>e>0. \<exists>N. \<forall>n\<ge>N. cmod (g (f n) - z) < e"
   525       by blast
   526     {
   527       fix w
   528       assume wr: "cmod w \<le> r"
   529       let ?e = "\<bar>cmod (poly p z) - ?m\<bar>"
   530       {
   531         assume e: "?e > 0"
   532         then have e2: "?e/2 > 0"
   533           by simp
   534         from poly_cont[OF e2, of z p] obtain d where
   535             d: "d > 0" "\<forall>w. 0<cmod (w - z)\<and> cmod(w - z) < d \<longrightarrow> cmod(poly p w - poly p z) < ?e/2"
   536           by blast
   537         {
   538           fix w
   539           assume w: "cmod (w - z) < d"
   540           have "cmod(poly p w - poly p z) < ?e / 2"
   541             using d(2)[rule_format, of w] w e by (cases "w = z") simp_all
   542         }
   543         note th1 = this
   544 
   545         from fz(2) d(1) obtain N1 where N1: "\<forall>n\<ge>N1. cmod (g (f n) - z) < d"
   546           by blast
   547         from reals_Archimedean2[of "2/?e"] obtain N2 :: nat where N2: "2/?e < real N2"
   548           by blast
   549         have th2: "cmod (poly p (g (f (N1 + N2))) - poly p z) < ?e/2"
   550           using N1[rule_format, of "N1 + N2"] th1 by simp
   551         {
   552           fix a b e2 m :: real
   553           have "a < e2 \<Longrightarrow> \<bar>b - m\<bar> < e2 \<Longrightarrow> 2 * e2 \<le> \<bar>b - m\<bar> + a \<Longrightarrow> False"
   554             by arith
   555         }
   556         note th0 = this
   557         have ath: "\<And>m x e::real. m \<le> x \<Longrightarrow> x < m + e \<Longrightarrow> \<bar>x - m\<bar> < e"
   558           by arith
   559         from s1m[OF g(1)[rule_format]] have th31: "?m \<le> cmod(poly p (g (f (N1 + N2))))" .
   560         from seq_suble[OF fz(1), of "N1 + N2"]
   561         have th00: "real (Suc (N1 + N2)) \<le> real (Suc (f (N1 + N2)))"
   562           by simp
   563         have th000: "0 \<le> (1::real)" "(1::real) \<le> 1" "real (Suc (N1 + N2)) > 0"
   564           using N2 by auto
   565         from frac_le[OF th000 th00]
   566         have th00: "?m + 1 / real (Suc (f (N1 + N2))) \<le> ?m + 1 / real (Suc (N1 + N2))"
   567           by simp
   568         from g(2)[rule_format, of "f (N1 + N2)"]
   569         have th01:"cmod (poly p (g (f (N1 + N2)))) < - s + 1 / real (Suc (f (N1 + N2)))" .
   570         from order_less_le_trans[OF th01 th00]
   571         have th32: "cmod (poly p (g (f (N1 + N2)))) < ?m + (1/ real(Suc (N1 + N2)))" .
   572         from N2 have "2/?e < real (Suc (N1 + N2))"
   573           by arith
   574         with e2 less_imp_inverse_less[of "2/?e" "real (Suc (N1 + N2))"]
   575         have "?e/2 > 1/ real (Suc (N1 + N2))"
   576           by (simp add: inverse_eq_divide)
   577         with ath[OF th31 th32]
   578         have thc1: "\<bar>cmod (poly p (g (f (N1 + N2)))) - ?m\<bar> < ?e/2"
   579           by arith
   580         have ath2: "\<And>a b c m::real. \<bar>a - b\<bar> \<le> c \<Longrightarrow> \<bar>b - m\<bar> \<le> \<bar>a - m\<bar> + c"
   581           by arith
   582         have th22: "\<bar>cmod (poly p (g (f (N1 + N2)))) - cmod (poly p z)\<bar> \<le>
   583             cmod (poly p (g (f (N1 + N2))) - poly p z)"
   584           by (simp add: norm_triangle_ineq3)
   585         from ath2[OF th22, of ?m]
   586         have thc2: "2 * (?e/2) \<le>
   587             \<bar>cmod(poly p (g (f (N1 + N2)))) - ?m\<bar> + cmod (poly p (g (f (N1 + N2))) - poly p z)"
   588           by simp
   589         from th0[OF th2 thc1 thc2] have False .
   590       }
   591       then have "?e = 0"
   592         by auto
   593       then have "cmod (poly p z) = ?m"
   594         by simp
   595       with s1m[OF wr] have "cmod (poly p z) \<le> cmod (poly p w)"
   596         by simp
   597     }
   598     then have ?thesis by blast
   599   }
   600   ultimately show ?thesis by blast
   601 qed
   602 
   603 lemma "(rcis (sqrt (abs r)) (a/2))\<^sup>2 = rcis (abs r) a"
   604   unfolding power2_eq_square
   605   apply (simp add: rcis_mult)
   606   apply (simp add: power2_eq_square[symmetric])
   607   done
   608 
   609 lemma cispi: "cis pi = -1"
   610   by (simp add: cis_def)
   611 
   612 lemma "(rcis (sqrt (abs r)) ((pi + a) / 2))\<^sup>2 = rcis (- abs r) a"
   613   unfolding power2_eq_square
   614   apply (simp add: rcis_mult add_divide_distrib)
   615   apply (simp add: power2_eq_square[symmetric] rcis_def cispi cis_mult[symmetric])
   616   done
   617 
   618 text {* Nonzero polynomial in z goes to infinity as z does. *}
   619 
   620 lemma poly_infinity:
   621   fixes p:: "'a::{comm_semiring_0,real_normed_div_algebra} poly"
   622   assumes ex: "p \<noteq> 0"
   623   shows "\<exists>r. \<forall>z. r \<le> norm z \<longrightarrow> d \<le> norm (poly (pCons a p) z)"
   624   using ex
   625 proof (induct p arbitrary: a d)
   626   case 0
   627   then show ?case by simp
   628 next
   629   case (pCons c cs a d)
   630   show ?case
   631   proof (cases "cs = 0")
   632     case False
   633     with pCons.hyps obtain r where r: "\<forall>z. r \<le> norm z \<longrightarrow> d + norm a \<le> norm (poly (pCons c cs) z)"
   634       by blast
   635     let ?r = "1 + \<bar>r\<bar>"
   636     {
   637       fix z :: 'a
   638       assume h: "1 + \<bar>r\<bar> \<le> norm z"
   639       have r0: "r \<le> norm z"
   640         using h by arith
   641       from r[rule_format, OF r0] have th0: "d + norm a \<le> 1 * norm(poly (pCons c cs) z)"
   642         by arith
   643       from h have z1: "norm z \<ge> 1"
   644         by arith
   645       from order_trans[OF th0 mult_right_mono[OF z1 norm_ge_zero[of "poly (pCons c cs) z"]]]
   646       have th1: "d \<le> norm(z * poly (pCons c cs) z) - norm a"
   647         unfolding norm_mult by (simp add: algebra_simps)
   648       from norm_diff_ineq[of "z * poly (pCons c cs) z" a]
   649       have th2: "norm (z * poly (pCons c cs) z) - norm a \<le> norm (poly (pCons a (pCons c cs)) z)"
   650         by (simp add: algebra_simps)
   651       from th1 th2 have "d \<le> norm (poly (pCons a (pCons c cs)) z)"
   652         by arith
   653     }
   654     then show ?thesis by blast
   655   next
   656     case True
   657     with pCons.prems have c0: "c \<noteq> 0"
   658       by simp
   659     {
   660       fix z :: 'a
   661       assume h: "(\<bar>d\<bar> + norm a) / norm c \<le> norm z"
   662       from c0 have "norm c > 0"
   663         by simp
   664       from h c0 have th0: "\<bar>d\<bar> + norm a \<le> norm (z * c)"
   665         by (simp add: field_simps norm_mult)
   666       have ath: "\<And>mzh mazh ma. mzh \<le> mazh + ma \<Longrightarrow> \<bar>d\<bar> + ma \<le> mzh \<Longrightarrow> d \<le> mazh"
   667         by arith
   668       from norm_diff_ineq[of "z * c" a] have th1: "norm (z * c) \<le> norm (a + z * c) + norm a"
   669         by (simp add: algebra_simps)
   670       from ath[OF th1 th0] have "d \<le> norm (poly (pCons a (pCons c cs)) z)"
   671         using True by simp
   672     }
   673     then show ?thesis by blast
   674   qed
   675 qed
   676 
   677 text {* Hence polynomial's modulus attains its minimum somewhere. *}
   678 lemma poly_minimum_modulus: "\<exists>z.\<forall>w. cmod (poly p z) \<le> cmod (poly p w)"
   679 proof (induct p)
   680   case 0
   681   then show ?case by simp
   682 next
   683   case (pCons c cs)
   684   show ?case
   685   proof (cases "cs = 0")
   686     case False
   687     from poly_infinity[OF False, of "cmod (poly (pCons c cs) 0)" c]
   688     obtain r where r: "\<And>z. r \<le> cmod z \<Longrightarrow> cmod (poly (pCons c cs) 0) \<le> cmod (poly (pCons c cs) z)"
   689       by blast
   690     have ath: "\<And>z r. r \<le> cmod z \<or> cmod z \<le> \<bar>r\<bar>"
   691       by arith
   692     from poly_minimum_modulus_disc[of "\<bar>r\<bar>" "pCons c cs"]
   693     obtain v where v: "\<And>w. cmod w \<le> \<bar>r\<bar> \<Longrightarrow> cmod (poly (pCons c cs) v) \<le> cmod (poly (pCons c cs) w)"
   694       by blast
   695     {
   696       fix z
   697       assume z: "r \<le> cmod z"
   698       from v[of 0] r[OF z] have "cmod (poly (pCons c cs) v) \<le> cmod (poly (pCons c cs) z)"
   699         by simp
   700     }
   701     note v0 = this
   702     from v0 v ath[of r] show ?thesis
   703       by blast
   704   next
   705     case True
   706     with pCons.hyps show ?thesis by simp
   707   qed
   708 qed
   709 
   710 text{* Constant function (non-syntactic characterization). *}
   711 definition "constant f \<longleftrightarrow> (\<forall>x y. f x = f y)"
   712 
   713 lemma nonconstant_length: "\<not> constant (poly p) \<Longrightarrow> psize p \<ge> 2"
   714   by (induct p) (auto simp: constant_def psize_def)
   715 
   716 lemma poly_replicate_append: "poly (monom 1 n * p) (x::'a::comm_ring_1) = x^n * poly p x"
   717   by (simp add: poly_monom)
   718 
   719 text {* Decomposition of polynomial, skipping zero coefficients
   720   after the first.  *}
   721 
   722 lemma poly_decompose_lemma:
   723   assumes nz: "\<not> (\<forall>z. z \<noteq> 0 \<longrightarrow> poly p z = (0::'a::idom))"
   724   shows "\<exists>k a q. a \<noteq> 0 \<and> Suc (psize q + k) = psize p \<and> (\<forall>z. poly p z = z^k * poly (pCons a q) z)"
   725   unfolding psize_def
   726   using nz
   727 proof (induct p)
   728   case 0
   729   then show ?case by simp
   730 next
   731   case (pCons c cs)
   732   show ?case
   733   proof (cases "c = 0")
   734     case True
   735     from pCons.hyps pCons.prems True show ?thesis
   736       apply (auto)
   737       apply (rule_tac x="k+1" in exI)
   738       apply (rule_tac x="a" in exI, clarsimp)
   739       apply (rule_tac x="q" in exI)
   740       apply auto
   741       done
   742   next
   743     case False
   744     show ?thesis
   745       apply (rule exI[where x=0])
   746       apply (rule exI[where x=c], auto simp add: False)
   747       done
   748   qed
   749 qed
   750 
   751 lemma poly_decompose:
   752   assumes nc: "\<not> constant (poly p)"
   753   shows "\<exists>k a q. a \<noteq> (0::'a::idom) \<and> k \<noteq> 0 \<and>
   754                psize q + k + 1 = psize p \<and>
   755               (\<forall>z. poly p z = poly p 0 + z^k * poly (pCons a q) z)"
   756   using nc
   757 proof (induct p)
   758   case 0
   759   then show ?case
   760     by (simp add: constant_def)
   761 next
   762   case (pCons c cs)
   763   {
   764     assume C: "\<forall>z. z \<noteq> 0 \<longrightarrow> poly cs z = 0"
   765     {
   766       fix x y
   767       from C have "poly (pCons c cs) x = poly (pCons c cs) y"
   768         by (cases "x = 0") auto
   769     }
   770     with pCons.prems have False
   771       by (auto simp add: constant_def)
   772   }
   773   then have th: "\<not> (\<forall>z. z \<noteq> 0 \<longrightarrow> poly cs z = 0)" ..
   774   from poly_decompose_lemma[OF th]
   775   show ?case
   776     apply clarsimp
   777     apply (rule_tac x="k+1" in exI)
   778     apply (rule_tac x="a" in exI)
   779     apply simp
   780     apply (rule_tac x="q" in exI)
   781     apply (auto simp add: psize_def split: if_splits)
   782     done
   783 qed
   784 
   785 text{* Fundamental theorem of algebra *}
   786 
   787 lemma fundamental_theorem_of_algebra:
   788   assumes nc: "\<not> constant (poly p)"
   789   shows "\<exists>z::complex. poly p z = 0"
   790   using nc
   791 proof (induct "psize p" arbitrary: p rule: less_induct)
   792   case less
   793   let ?p = "poly p"
   794   let ?ths = "\<exists>z. ?p z = 0"
   795 
   796   from nonconstant_length[OF less(2)] have n2: "psize p \<ge> 2" .
   797   from poly_minimum_modulus obtain c where c: "\<forall>w. cmod (?p c) \<le> cmod (?p w)"
   798     by blast
   799 
   800   show ?ths
   801   proof (cases "?p c = 0")
   802     case True
   803     then show ?thesis by blast
   804   next
   805     case False
   806     note pc0 = this
   807     from poly_offset[of p c] obtain q where q: "psize q = psize p" "\<forall>x. poly q x = ?p (c + x)"
   808       by blast
   809     {
   810       assume h: "constant (poly q)"
   811       from q(2) have th: "\<forall>x. poly q (x - c) = ?p x"
   812         by auto
   813       {
   814         fix x y
   815         from th have "?p x = poly q (x - c)"
   816           by auto
   817         also have "\<dots> = poly q (y - c)"
   818           using h unfolding constant_def by blast
   819         also have "\<dots> = ?p y"
   820           using th by auto
   821         finally have "?p x = ?p y" .
   822       }
   823       with less(2) have False
   824         unfolding constant_def by blast
   825     }
   826     then have qnc: "\<not> constant (poly q)"
   827       by blast
   828     from q(2) have pqc0: "?p c = poly q 0"
   829       by simp
   830     from c pqc0 have cq0: "\<forall>w. cmod (poly q 0) \<le> cmod (?p w)"
   831       by simp
   832     let ?a0 = "poly q 0"
   833     from pc0 pqc0 have a00: "?a0 \<noteq> 0"
   834       by simp
   835     from a00 have qr: "\<forall>z. poly q z = poly (smult (inverse ?a0) q) z * ?a0"
   836       by simp
   837     let ?r = "smult (inverse ?a0) q"
   838     have lgqr: "psize q = psize ?r"
   839       using a00
   840       unfolding psize_def degree_def
   841       by (simp add: poly_eq_iff)
   842     {
   843       assume h: "\<And>x y. poly ?r x = poly ?r y"
   844       {
   845         fix x y
   846         from qr[rule_format, of x] have "poly q x = poly ?r x * ?a0"
   847           by auto
   848         also have "\<dots> = poly ?r y * ?a0"
   849           using h by simp
   850         also have "\<dots> = poly q y"
   851           using qr[rule_format, of y] by simp
   852         finally have "poly q x = poly q y" .
   853       }
   854       with qnc have False
   855         unfolding constant_def by blast
   856     }
   857     then have rnc: "\<not> constant (poly ?r)"
   858       unfolding constant_def by blast
   859     from qr[rule_format, of 0] a00 have r01: "poly ?r 0 = 1"
   860       by auto
   861     {
   862       fix w
   863       have "cmod (poly ?r w) < 1 \<longleftrightarrow> cmod (poly q w / ?a0) < 1"
   864         using qr[rule_format, of w] a00 by (simp add: divide_inverse mult_ac)
   865       also have "\<dots> \<longleftrightarrow> cmod (poly q w) < cmod ?a0"
   866         using a00 unfolding norm_divide by (simp add: field_simps)
   867       finally have "cmod (poly ?r w) < 1 \<longleftrightarrow> cmod (poly q w) < cmod ?a0" .
   868     }
   869     note mrmq_eq = this
   870     from poly_decompose[OF rnc] obtain k a s where
   871       kas: "a \<noteq> 0" "k \<noteq> 0" "psize s + k + 1 = psize ?r"
   872         "\<forall>z. poly ?r z = poly ?r 0 + z^k* poly (pCons a s) z" by blast
   873     {
   874       assume "psize p = k + 1"
   875       with kas(3) lgqr[symmetric] q(1) have s0: "s = 0"
   876         by auto
   877       {
   878         fix w
   879         have "cmod (poly ?r w) = cmod (1 + a * w ^ k)"
   880           using kas(4)[rule_format, of w] s0 r01 by (simp add: algebra_simps)
   881       }
   882       note hth = this [symmetric]
   883       from reduce_poly_simple[OF kas(1,2)] have "\<exists>w. cmod (poly ?r w) < 1"
   884         unfolding hth by blast
   885     }
   886     moreover
   887     {
   888       assume kn: "psize p \<noteq> k + 1"
   889       from kn kas(3) q(1) lgqr have k1n: "k + 1 < psize p"
   890         by simp
   891       have th01: "\<not> constant (poly (pCons 1 (monom a (k - 1))))"
   892         unfolding constant_def poly_pCons poly_monom
   893         using kas(1)
   894         apply simp
   895         apply (rule exI[where x=0])
   896         apply (rule exI[where x=1])
   897         apply simp
   898         done
   899       from kas(1) kas(2) have th02: "k + 1 = psize (pCons 1 (monom a (k - 1)))"
   900         by (simp add: psize_def degree_monom_eq)
   901       from less(1) [OF k1n [simplified th02] th01]
   902       obtain w where w: "1 + w^k * a = 0"
   903         unfolding poly_pCons poly_monom
   904         using kas(2) by (cases k) (auto simp add: algebra_simps)
   905       from poly_bound_exists[of "cmod w" s] obtain m where
   906         m: "m > 0" "\<forall>z. cmod z \<le> cmod w \<longrightarrow> cmod (poly s z) \<le> m" by blast
   907       have w0: "w \<noteq> 0"
   908         using kas(2) w by (auto simp add: power_0_left)
   909       from w have "(1 + w ^ k * a) - 1 = 0 - 1"
   910         by simp
   911       then have wm1: "w^k * a = - 1"
   912         by simp
   913       have inv0: "0 < inverse (cmod w ^ (k + 1) * m)"
   914         using norm_ge_zero[of w] w0 m(1)
   915         by (simp add: inverse_eq_divide zero_less_mult_iff)
   916       with real_lbound_gt_zero[OF zero_less_one] obtain t where
   917         t: "t > 0" "t < 1" "t < inverse (cmod w ^ (k + 1) * m)" by blast
   918       let ?ct = "complex_of_real t"
   919       let ?w = "?ct * w"
   920       have "1 + ?w^k * (a + ?w * poly s ?w) = 1 + ?ct^k * (w^k * a) + ?w^k * ?w * poly s ?w"
   921         using kas(1) by (simp add: algebra_simps power_mult_distrib)
   922       also have "\<dots> = complex_of_real (1 - t^k) + ?w^k * ?w * poly s ?w"
   923         unfolding wm1 by simp
   924       finally have "cmod (1 + ?w^k * (a + ?w * poly s ?w)) =
   925         cmod (complex_of_real (1 - t^k) + ?w^k * ?w * poly s ?w)"
   926         by metis
   927       with norm_triangle_ineq[of "complex_of_real (1 - t^k)" "?w^k * ?w * poly s ?w"]
   928       have th11: "cmod (1 + ?w^k * (a + ?w * poly s ?w)) \<le> \<bar>1 - t^k\<bar> + cmod (?w^k * ?w * poly s ?w)"
   929         unfolding norm_of_real by simp
   930       have ath: "\<And>x t::real. 0 \<le> x \<Longrightarrow> x < t \<Longrightarrow> t \<le> 1 \<Longrightarrow> \<bar>1 - t\<bar> + x < 1"
   931         by arith
   932       have "t * cmod w \<le> 1 * cmod w"
   933         apply (rule mult_mono)
   934         using t(1,2)
   935         apply auto
   936         done
   937       then have tw: "cmod ?w \<le> cmod w"
   938         using t(1) by (simp add: norm_mult)
   939       from t inv0 have "t * (cmod w ^ (k + 1) * m) < 1"
   940         by (simp add: inverse_eq_divide field_simps)
   941       with zero_less_power[OF t(1), of k] have th30: "t^k * (t* (cmod w ^ (k + 1) * m)) < t^k * 1"
   942         by (metis comm_mult_strict_left_mono)
   943       have "cmod (?w^k * ?w * poly s ?w) = t^k * (t* (cmod w ^ (k + 1) * cmod (poly s ?w)))"
   944         using w0 t(1)
   945         by (simp add: algebra_simps power_mult_distrib norm_power norm_mult)
   946       then have "cmod (?w^k * ?w * poly s ?w) \<le> t^k * (t* (cmod w ^ (k + 1) * m))"
   947         using t(1,2) m(2)[rule_format, OF tw] w0
   948         by auto
   949       with th30 have th120: "cmod (?w^k * ?w * poly s ?w) < t^k"
   950         by simp
   951       from power_strict_mono[OF t(2), of k] t(1) kas(2) have th121: "t^k \<le> 1"
   952         by auto
   953       from ath[OF norm_ge_zero[of "?w^k * ?w * poly s ?w"] th120 th121]
   954       have th12: "\<bar>1 - t^k\<bar> + cmod (?w^k * ?w * poly s ?w) < 1" .
   955       from th11 th12 have "cmod (1 + ?w^k * (a + ?w * poly s ?w)) < 1"
   956         by arith
   957       then have "cmod (poly ?r ?w) < 1"
   958         unfolding kas(4)[rule_format, of ?w] r01 by simp
   959       then have "\<exists>w. cmod (poly ?r w) < 1"
   960         by blast
   961     }
   962     ultimately have cr0_contr: "\<exists>w. cmod (poly ?r w) < 1"
   963       by blast
   964     from cr0_contr cq0 q(2) show ?thesis
   965       unfolding mrmq_eq not_less[symmetric] by auto
   966   qed
   967 qed
   968 
   969 text {* Alternative version with a syntactic notion of constant polynomial. *}
   970 
   971 lemma fundamental_theorem_of_algebra_alt:
   972   assumes nc: "\<not> (\<exists>a l. a \<noteq> 0 \<and> l = 0 \<and> p = pCons a l)"
   973   shows "\<exists>z. poly p z = (0::complex)"
   974   using nc
   975 proof (induct p)
   976   case 0
   977   then show ?case by simp
   978 next
   979   case (pCons c cs)
   980   show ?case
   981   proof (cases "c = 0")
   982     case True
   983     then show ?thesis by auto
   984   next
   985     case False
   986     {
   987       assume nc: "constant (poly (pCons c cs))"
   988       from nc[unfolded constant_def, rule_format, of 0]
   989       have "\<forall>w. w \<noteq> 0 \<longrightarrow> poly cs w = 0" by auto
   990       then have "cs = 0"
   991       proof (induct cs)
   992         case 0
   993         then show ?case by simp
   994       next
   995         case (pCons d ds)
   996         show ?case
   997         proof (cases "d = 0")
   998           case True
   999           then show ?thesis using pCons.prems pCons.hyps by simp
  1000         next
  1001           case False
  1002           from poly_bound_exists[of 1 ds] obtain m where
  1003             m: "m > 0" "\<forall>z. \<forall>z. cmod z \<le> 1 \<longrightarrow> cmod (poly ds z) \<le> m" by blast
  1004           have dm: "cmod d / m > 0"
  1005             using False m(1) by (simp add: field_simps)
  1006           from real_lbound_gt_zero[OF dm zero_less_one] obtain x where
  1007             x: "x > 0" "x < cmod d / m" "x < 1" by blast
  1008           let ?x = "complex_of_real x"
  1009           from x have cx: "?x \<noteq> 0"  "cmod ?x \<le> 1"
  1010             by simp_all
  1011           from pCons.prems[rule_format, OF cx(1)]
  1012           have cth: "cmod (?x*poly ds ?x) = cmod d"
  1013             by (simp add: eq_diff_eq[symmetric])
  1014           from m(2)[rule_format, OF cx(2)] x(1)
  1015           have th0: "cmod (?x*poly ds ?x) \<le> x*m"
  1016             by (simp add: norm_mult)
  1017           from x(2) m(1) have "x * m < cmod d"
  1018             by (simp add: field_simps)
  1019           with th0 have "cmod (?x*poly ds ?x) \<noteq> cmod d"
  1020             by auto
  1021           with cth show ?thesis
  1022             by blast
  1023         qed
  1024       qed
  1025     }
  1026     then have nc: "\<not> constant (poly (pCons c cs))"
  1027       using pCons.prems False by blast
  1028     from fundamental_theorem_of_algebra[OF nc] show ?thesis .
  1029   qed
  1030 qed
  1031 
  1032 
  1033 subsection{* Nullstellensatz, degrees and divisibility of polynomials *}
  1034 
  1035 lemma nullstellensatz_lemma:
  1036   fixes p :: "complex poly"
  1037   assumes "\<forall>x. poly p x = 0 \<longrightarrow> poly q x = 0"
  1038     and "degree p = n"
  1039     and "n \<noteq> 0"
  1040   shows "p dvd (q ^ n)"
  1041   using assms
  1042 proof (induct n arbitrary: p q rule: nat_less_induct)
  1043   fix n :: nat
  1044   fix p q :: "complex poly"
  1045   assume IH: "\<forall>m<n. \<forall>p q.
  1046                  (\<forall>x. poly p x = (0::complex) \<longrightarrow> poly q x = 0) \<longrightarrow>
  1047                  degree p = m \<longrightarrow> m \<noteq> 0 \<longrightarrow> p dvd (q ^ m)"
  1048     and pq0: "\<forall>x. poly p x = 0 \<longrightarrow> poly q x = 0"
  1049     and dpn: "degree p = n"
  1050     and n0: "n \<noteq> 0"
  1051   from dpn n0 have pne: "p \<noteq> 0" by auto
  1052   let ?ths = "p dvd (q ^ n)"
  1053   {
  1054     fix a
  1055     assume a: "poly p a = 0"
  1056     {
  1057       assume oa: "order a p \<noteq> 0"
  1058       let ?op = "order a p"
  1059       from pne have ap: "([:- a, 1:] ^ ?op) dvd p" "\<not> [:- a, 1:] ^ (Suc ?op) dvd p"
  1060         using order by blast+
  1061       note oop = order_degree[OF pne, unfolded dpn]
  1062       {
  1063         assume q0: "q = 0"
  1064         then have ?ths using n0
  1065           by (simp add: power_0_left)
  1066       }
  1067       moreover
  1068       {
  1069         assume q0: "q \<noteq> 0"
  1070         from pq0[rule_format, OF a, unfolded poly_eq_0_iff_dvd]
  1071         obtain r where r: "q = [:- a, 1:] * r" by (rule dvdE)
  1072         from ap(1) obtain s where s: "p = [:- a, 1:] ^ ?op * s"
  1073           by (rule dvdE)
  1074         have sne: "s \<noteq> 0" using s pne by auto
  1075         {
  1076           assume ds0: "degree s = 0"
  1077           from ds0 obtain k where kpn: "s = [:k:]"
  1078             by (cases s) (auto split: if_splits)
  1079           from sne kpn have k: "k \<noteq> 0" by simp
  1080           let ?w = "([:1/k:] * ([:-a,1:] ^ (n - ?op))) * (r ^ n)"
  1081           have "q ^ n = p * ?w"
  1082             apply (subst r)
  1083             apply (subst s)
  1084             apply (subst kpn)
  1085             using k oop [of a]
  1086             apply (subst power_mult_distrib)
  1087             apply simp
  1088             apply (subst power_add [symmetric])
  1089             apply simp
  1090             done
  1091           then have ?ths
  1092             unfolding dvd_def by blast
  1093         }
  1094         moreover
  1095         {
  1096           assume ds0: "degree s \<noteq> 0"
  1097           from ds0 sne dpn s oa
  1098             have dsn: "degree s < n"
  1099               apply auto
  1100               apply (erule ssubst)
  1101               apply (simp add: degree_mult_eq degree_linear_power)
  1102               done
  1103             {
  1104               fix x assume h: "poly s x = 0"
  1105               {
  1106                 assume xa: "x = a"
  1107                 from h[unfolded xa poly_eq_0_iff_dvd] obtain u where u: "s = [:- a, 1:] * u"
  1108                   by (rule dvdE)
  1109                 have "p = [:- a, 1:] ^ (Suc ?op) * u"
  1110                   apply (subst s)
  1111                   apply (subst u)
  1112                   apply (simp only: power_Suc mult_ac)
  1113                   done
  1114                 with ap(2)[unfolded dvd_def] have False
  1115                   by blast
  1116               }
  1117               note xa = this
  1118               from h have "poly p x = 0"
  1119                 by (subst s) simp
  1120               with pq0 have "poly q x = 0"
  1121                 by blast
  1122               with r xa have "poly r x = 0"
  1123                 by auto
  1124             }
  1125             note impth = this
  1126             from IH[rule_format, OF dsn, of s r] impth ds0
  1127             have "s dvd (r ^ (degree s))"
  1128               by blast
  1129             then obtain u where u: "r ^ (degree s) = s * u" ..
  1130             then have u': "\<And>x. poly s x * poly u x = poly r x ^ degree s"
  1131               by (simp only: poly_mult[symmetric] poly_power[symmetric])
  1132             let ?w = "(u * ([:-a,1:] ^ (n - ?op))) * (r ^ (n - degree s))"
  1133             from oop[of a] dsn have "q ^ n = p * ?w"
  1134               apply -
  1135               apply (subst s)
  1136               apply (subst r)
  1137               apply (simp only: power_mult_distrib)
  1138               apply (subst mult_assoc [where b=s])
  1139               apply (subst mult_assoc [where a=u])
  1140               apply (subst mult_assoc [where b=u, symmetric])
  1141               apply (subst u [symmetric])
  1142               apply (simp add: mult_ac power_add [symmetric])
  1143               done
  1144             then have ?ths
  1145               unfolding dvd_def by blast
  1146         }
  1147         ultimately have ?ths by blast
  1148       }
  1149       ultimately have ?ths by blast
  1150     }
  1151     then have ?ths using a order_root pne by blast
  1152   }
  1153   moreover
  1154   {
  1155     assume exa: "\<not> (\<exists>a. poly p a = 0)"
  1156     from fundamental_theorem_of_algebra_alt[of p] exa
  1157     obtain c where ccs: "c \<noteq> 0" "p = pCons c 0"
  1158       by blast
  1159     then have pp: "\<And>x. poly p x = c"
  1160       by simp
  1161     let ?w = "[:1/c:] * (q ^ n)"
  1162     from ccs have "(q ^ n) = (p * ?w)"
  1163       by simp
  1164     then have ?ths
  1165       unfolding dvd_def by blast
  1166   }
  1167   ultimately show ?ths by blast
  1168 qed
  1169 
  1170 lemma nullstellensatz_univariate:
  1171   "(\<forall>x. poly p x = (0::complex) \<longrightarrow> poly q x = 0) \<longleftrightarrow>
  1172     p dvd (q ^ (degree p)) \<or> (p = 0 \<and> q = 0)"
  1173 proof -
  1174   {
  1175     assume pe: "p = 0"
  1176     then have eq: "(\<forall>x. poly p x = (0::complex) \<longrightarrow> poly q x = 0) \<longleftrightarrow> q = 0"
  1177       by (auto simp add: poly_all_0_iff_0)
  1178     {
  1179       assume "p dvd (q ^ (degree p))"
  1180       then obtain r where r: "q ^ (degree p) = p * r" ..
  1181       from r pe have False by simp
  1182     }
  1183     with eq pe have ?thesis by blast
  1184   }
  1185   moreover
  1186   {
  1187     assume pe: "p \<noteq> 0"
  1188     {
  1189       assume dp: "degree p = 0"
  1190       then obtain k where k: "p = [:k:]" "k \<noteq> 0" using pe
  1191         by (cases p) (simp split: if_splits)
  1192       then have th1: "\<forall>x. poly p x \<noteq> 0"
  1193         by simp
  1194       from k dp have "q ^ (degree p) = p * [:1/k:]"
  1195         by (simp add: one_poly_def)
  1196       then have th2: "p dvd (q ^ (degree p))" ..
  1197       from th1 th2 pe have ?thesis
  1198         by blast
  1199     }
  1200     moreover
  1201     {
  1202       assume dp: "degree p \<noteq> 0"
  1203       then obtain n where n: "degree p = Suc n "
  1204         by (cases "degree p") auto
  1205       {
  1206         assume "p dvd (q ^ (Suc n))"
  1207         then obtain u where u: "q ^ (Suc n) = p * u" ..
  1208         {
  1209           fix x
  1210           assume h: "poly p x = 0" "poly q x \<noteq> 0"
  1211           then have "poly (q ^ (Suc n)) x \<noteq> 0"
  1212             by simp
  1213           then have False using u h(1)
  1214             by (simp only: poly_mult) simp
  1215         }
  1216       }
  1217       with n nullstellensatz_lemma[of p q "degree p"] dp
  1218       have ?thesis by auto
  1219     }
  1220     ultimately have ?thesis by blast
  1221   }
  1222   ultimately show ?thesis by blast
  1223 qed
  1224 
  1225 text {* Useful lemma *}
  1226 
  1227 lemma constant_degree:
  1228   fixes p :: "'a::{idom,ring_char_0} poly"
  1229   shows "constant (poly p) \<longleftrightarrow> degree p = 0" (is "?lhs = ?rhs")
  1230 proof
  1231   assume l: ?lhs
  1232   from l[unfolded constant_def, rule_format, of _ "0"]
  1233   have th: "poly p = poly [:poly p 0:]"
  1234     by auto
  1235   then have "p = [:poly p 0:]"
  1236     by (simp add: poly_eq_poly_eq_iff)
  1237   then have "degree p = degree [:poly p 0:]"
  1238     by simp
  1239   then show ?rhs
  1240     by simp
  1241 next
  1242   assume r: ?rhs
  1243   then obtain k where "p = [:k:]"
  1244     by (cases p) (simp split: if_splits)
  1245   then show ?lhs
  1246     unfolding constant_def by auto
  1247 qed
  1248 
  1249 lemma divides_degree:
  1250   assumes pq: "p dvd (q:: complex poly)"
  1251   shows "degree p \<le> degree q \<or> q = 0"
  1252   by (metis dvd_imp_degree_le pq)
  1253 
  1254 text {* Arithmetic operations on multivariate polynomials. *}
  1255 
  1256 lemma mpoly_base_conv:
  1257   fixes x :: "'a::comm_ring_1"
  1258   shows "0 = poly 0 x" "c = poly [:c:] x" "x = poly [:0,1:] x"
  1259   by simp_all
  1260 
  1261 lemma mpoly_norm_conv:
  1262   fixes x :: "'a::comm_ring_1"
  1263   shows "poly [:0:] x = poly 0 x" "poly [:poly 0 y:] x = poly 0 x"
  1264   by simp_all
  1265 
  1266 lemma mpoly_sub_conv:
  1267   fixes x :: "'a::comm_ring_1"
  1268   shows "poly p x - poly q x = poly p x + -1 * poly q x"
  1269   by simp
  1270 
  1271 lemma poly_pad_rule: "poly p x = 0 \<Longrightarrow> poly (pCons 0 p) x = 0"
  1272   by simp
  1273 
  1274 lemma poly_cancel_eq_conv:
  1275   fixes x :: "'a::field"
  1276   shows "x = 0 \<Longrightarrow> a \<noteq> 0 \<Longrightarrow> y = 0 \<longleftrightarrow> a * y - b * x = 0"
  1277   by auto
  1278 
  1279 lemma poly_divides_pad_rule:
  1280   fixes p:: "('a::comm_ring_1) poly"
  1281   assumes pq: "p dvd q"
  1282   shows "p dvd (pCons 0 q)"
  1283 proof -
  1284   have "pCons 0 q = q * [:0,1:]" by simp
  1285   then have "q dvd (pCons 0 q)" ..
  1286   with pq show ?thesis by (rule dvd_trans)
  1287 qed
  1288 
  1289 lemma poly_divides_conv0:
  1290   fixes p:: "'a::field poly"
  1291   assumes lgpq: "degree q < degree p"
  1292     and lq: "p \<noteq> 0"
  1293   shows "p dvd q \<longleftrightarrow> q = 0" (is "?lhs \<longleftrightarrow> ?rhs")
  1294 proof
  1295   assume r: ?rhs
  1296   then have "q = p * 0" by simp
  1297   then show ?lhs ..
  1298 next
  1299   assume l: ?lhs
  1300   show ?rhs
  1301   proof (cases "q = 0")
  1302     case True
  1303     then show ?thesis by simp
  1304   next
  1305     assume q0: "q \<noteq> 0"
  1306     from l q0 have "degree p \<le> degree q"
  1307       by (rule dvd_imp_degree_le)
  1308     with lgpq show ?thesis by simp
  1309   qed
  1310 qed
  1311 
  1312 lemma poly_divides_conv1:
  1313   fixes p :: "'a::field poly"
  1314   assumes a0: "a \<noteq> 0"
  1315     and pp': "p dvd p'"
  1316     and qrp': "smult a q - p' = r"
  1317   shows "p dvd q \<longleftrightarrow> p dvd r" (is "?lhs \<longleftrightarrow> ?rhs")
  1318 proof
  1319   from pp' obtain t where t: "p' = p * t" ..
  1320   {
  1321     assume l: ?lhs
  1322     then obtain u where u: "q = p * u" ..
  1323     have "r = p * (smult a u - t)"
  1324       using u qrp' [symmetric] t by (simp add: algebra_simps)
  1325     then show ?rhs ..
  1326   next
  1327     assume r: ?rhs
  1328     then obtain u where u: "r = p * u" ..
  1329     from u [symmetric] t qrp' [symmetric] a0
  1330     have "q = p * smult (1/a) (u + t)" by (simp add: algebra_simps)
  1331     then show ?lhs ..
  1332   }
  1333 qed
  1334 
  1335 lemma basic_cqe_conv1:
  1336   "(\<exists>x. poly p x = 0 \<and> poly 0 x \<noteq> 0) \<longleftrightarrow> False"
  1337   "(\<exists>x. poly 0 x \<noteq> 0) \<longleftrightarrow> False"
  1338   "(\<exists>x. poly [:c:] x \<noteq> 0) \<longleftrightarrow> c \<noteq> 0"
  1339   "(\<exists>x. poly 0 x = 0) \<longleftrightarrow> True"
  1340   "(\<exists>x. poly [:c:] x = 0) \<longleftrightarrow> c = 0"
  1341   by simp_all
  1342 
  1343 lemma basic_cqe_conv2:
  1344   assumes l: "p \<noteq> 0"
  1345   shows "\<exists>x. poly (pCons a (pCons b p)) x = (0::complex)"
  1346 proof -
  1347   {
  1348     fix h t
  1349     assume h: "h \<noteq> 0" "t = 0" and "pCons a (pCons b p) = pCons h t"
  1350     with l have False by simp
  1351   }
  1352   then have th: "\<not> (\<exists> h t. h \<noteq> 0 \<and> t = 0 \<and> pCons a (pCons b p) = pCons h t)"
  1353     by blast
  1354   from fundamental_theorem_of_algebra_alt[OF th] show ?thesis
  1355     by auto
  1356 qed
  1357 
  1358 lemma  basic_cqe_conv_2b: "(\<exists>x. poly p x \<noteq> (0::complex)) \<longleftrightarrow> p \<noteq> 0"
  1359   by (metis poly_all_0_iff_0)
  1360 
  1361 lemma basic_cqe_conv3:
  1362   fixes p q :: "complex poly"
  1363   assumes l: "p \<noteq> 0"
  1364   shows "(\<exists>x. poly (pCons a p) x = 0 \<and> poly q x \<noteq> 0) \<longleftrightarrow> \<not> (pCons a p) dvd (q ^ psize p)"
  1365 proof -
  1366   from l have dp: "degree (pCons a p) = psize p"
  1367     by (simp add: psize_def)
  1368   from nullstellensatz_univariate[of "pCons a p" q] l
  1369   show ?thesis
  1370     by (metis dp pCons_eq_0_iff)
  1371 qed
  1372 
  1373 lemma basic_cqe_conv4:
  1374   fixes p q :: "complex poly"
  1375   assumes h: "\<And>x. poly (q ^ n) x = poly r x"
  1376   shows "p dvd (q ^ n) \<longleftrightarrow> p dvd r"
  1377 proof -
  1378   from h have "poly (q ^ n) = poly r"
  1379     by auto
  1380   then have "(q ^ n) = r"
  1381     by (simp add: poly_eq_poly_eq_iff)
  1382   then show "p dvd (q ^ n) \<longleftrightarrow> p dvd r"
  1383     by simp
  1384 qed
  1385 
  1386 lemma poly_const_conv:
  1387   fixes x :: "'a::comm_ring_1"
  1388   shows "poly [:c:] x = y \<longleftrightarrow> c = y"
  1389   by simp
  1390 
  1391 end