Theory Subresultant
section ‹Subresultants and the subresultant PRS›
text ‹This theory contains most of the soundness proofs of the subresultant PRS algorithm,
where we closely follow the papers of Brown \<^cite>‹"Brown"› and Brown and Traub \<^cite>‹"BrownTraub"›.
This is in contrast to a similar Coq formalization of Mahboubi \<^cite>‹"Mahboubi06"› which
is based on polynomial determinants.
Whereas the current file only contains an algorithm to compute the resultant of two polynomials
efficiently, there is another theory ``Subresultant-Gcd'' which also contains the algorithm
to compute the GCD of two polynomials via the subresultant algorithm. In both algorithms we
integrate Lazard's optimization in the dichotomous version,
but not the second optmization described by Ducos \<^cite>‹"Ducos"›.›
theory Subresultant
imports
Resultant_Prelim
Dichotomous_Lazard
Binary_Exponentiation
More_Homomorphisms
Coeff_Int
begin
subsection ‹Algorithm›
locale div_exp_param =
fixes div_exp :: "'a :: idom_divide ⇒ 'a ⇒ nat ⇒ 'a"
begin
partial_function(tailrec) subresultant_prs_main where
"subresultant_prs_main f g c = (let
m = degree f;
n = degree g;
lf = lead_coeff f;
lg = lead_coeff g;
δ = m - n;
d = div_exp lg c δ;
h = pseudo_mod f g
in if h = 0 then (g,d)
else subresultant_prs_main g (sdiv_poly h ((-1) ^ (δ + 1) * lf * (c ^ δ))) d)"
definition subresultant_prs where
"subresultant_prs f g = (let
h = pseudo_mod f g;
δ = (degree f - degree g);
d = lead_coeff g ^ δ
in if h = 0 then (g,d)
else subresultant_prs_main g ((- 1) ^ (δ + 1) * h) d)"
definition resultant_impl_main where
"resultant_impl_main G1 G2 = (if G2 = 0 then (if degree G1 = 0 then 1 else 0) else
case subresultant_prs G1 G2 of
(Gk,hk) ⇒ (if degree Gk = 0 then hk else 0))"
definition resultant_impl where
"resultant_impl f g =
(if length (coeffs f) ≥ length (coeffs g) then resultant_impl_main f g
else let res = resultant_impl_main g f in
if even (degree f) ∨ even (degree g) then res else - res)"
end
locale div_exp_sound = div_exp_param +
assumes div_exp: "⋀ x y n.
(to_fract x)^n / (to_fract y)^(n-1) ∈ range to_fract
⟹ to_fract (div_exp x y n) = (to_fract x)^n / (to_fract y)^(n-1)"
definition basic_div_exp :: "'a :: idom_divide ⇒ 'a ⇒ nat ⇒ 'a" where
"basic_div_exp x y n = x^n div y^(n-1)"
text ‹We have an instance for arbitrary integral domains.›
lemma basic_div_exp: "div_exp_sound basic_div_exp"
by (unfold_locales, unfold basic_div_exp_def, rule sym, rule div_divide_to_fract, auto simp: hom_distribs)
text ‹Lazard's optimization is only proven for factorial rings.›
lemma dichotomous_Lazard: "div_exp_sound (dichotomous_Lazard :: 'a :: factorial_ring_gcd ⇒ _)"
by (unfold_locales, rule dichotomous_Lazard)
subsection ‹Soundness Proof for @{term "div_exp_param.resultant_impl div_exp = resultant"}›
abbreviation pdivmod :: "'a::field poly ⇒ 'a poly ⇒ 'a poly × 'a poly"
where
"pdivmod p q ≡ (p div q, p mod q)"
lemma even_sum_list: assumes "⋀ x. x ∈ set xs ⟹ even (f x) = even (g x)"
shows "even (sum_list (map f xs)) = even (sum_list (map g xs))"
using assms by (induct xs, auto)
lemma for_all_Suc: "P i ⟹ (∀ j ≥ Suc i. P j) = (∀ j ≥ i. P j)" for P
by (metis (full_types) Suc_le_eq less_le)
lemma pseudo_mod_left_0[simp]: "pseudo_mod 0 x = 0"
unfolding pseudo_mod_def pseudo_divmod_def
by (cases "x = 0"; cases "length (coeffs x)", auto)
lemma pseudo_mod_right_0[simp]: "pseudo_mod x 0 = x"
unfolding pseudo_mod_def pseudo_divmod_def by simp
lemma snd_pseudo_divmod_main_cong:
assumes "a1 = b1" "a3 = b3" "a4 = b4" "a5 = b5" "a6 = b6"
shows "snd (pseudo_divmod_main a1 a2 a3 a4 a5 a6) = snd (pseudo_divmod_main b1 b2 b3 b4 b5 b6)"
using assms snd_pseudo_divmod_main by metis
lemma snd_pseudo_mod_smult_invar_right:
shows "(snd (pseudo_divmod_main (x * lc) q r (smult x d) dr n))
= snd (pseudo_divmod_main lc q' (smult (x^n) r) d dr n)"
proof(induct n arbitrary: q q' r dr)
case (Suc n)
let ?q = "smult (x * lc) q + monom (coeff r dr) n"
let ?r = "smult (x * lc) r - (smult x (monom (coeff r dr) n * d))"
let ?dr = "dr - 1"
let ?rec_lhs = "pseudo_divmod_main (x * lc) ?q ?r (smult x d) ?dr n"
let ?rec_rhs = "pseudo_divmod_main lc q' (smult (x^n) ?r) d ?dr n"
have [simp]: "⋀ n. x ^ n * (x * lc) = lc * (x * x ^ n)"
"⋀ n c. x ^ n * (x * c) = x * x ^ n * c"
"⋀ n. x * x ^ n * lc = lc * (x * x ^ n)"
by (auto simp: ac_simps)
have "snd (pseudo_divmod_main (x * lc) q r (smult x d) dr (Suc n)) = snd ?rec_lhs"
by (auto simp:Let_def)
also have "… = snd ?rec_rhs" using Suc by auto
also have "… = snd (pseudo_divmod_main lc q' (smult (x^Suc n) r) d dr (Suc n))"
unfolding pseudo_divmod_main.simps Let_def
proof(rule snd_pseudo_divmod_main_cong,goal_cases)
case 2
show ?case by (auto simp:smult_add_right smult_diff_right smult_monom smult_monom_mult)
qed auto
finally show ?case by auto
qed auto
lemma snd_pseudo_mod_smult_invar_left:
shows "snd (pseudo_divmod_main lc q (smult x r) d dr n)
= smult x (snd (pseudo_divmod_main lc q' r d dr n))"
proof(induct n arbitrary:x lc q q' r d dr)
case (Suc n)
have sm:"smult lc (smult x r) - monom (coeff (smult x r) dr) n * d
=smult x (smult lc r - monom (coeff r dr) n * d) "
by (auto simp:smult_diff_right smult_monom smult_monom_mult mult.commute[of lc x])
let ?q' = "smult lc q' + monom (coeff r dr) n"
show ?case unfolding pseudo_divmod_main.simps Let_def Suc(1)[of lc _ _ _ _ _ ?q'] sm by auto
qed auto
lemma snd_pseudo_mod_smult_left[simp]:
shows "snd (pseudo_divmod (smult (x::'a::idom) p) q) = (smult x (snd (pseudo_divmod p q)))"
unfolding pseudo_divmod_def
by (auto simp:snd_pseudo_mod_smult_invar_left[of _ _ _ _ _ _ _ 0] Polynomial.coeffs_smult)
lemma pseudo_mod_smult_right:
assumes "(x::'a::idom)≠0" "q≠0"
shows "(pseudo_mod p (smult (x::'a::idom) q)) = (smult (x^(Suc (length (coeffs p)) - length (coeffs q))) (pseudo_mod p q))"
unfolding pseudo_divmod_def pseudo_mod_def
by (auto simp:snd_pseudo_mod_smult_invar_right[of _ _ _ _ _ _ _ 0]
snd_pseudo_mod_smult_invar_left[of _ _ _ _ _ _ _ 0] Polynomial.coeffs_smult assms)
lemma pseudo_mod_zero[simp]:
"pseudo_mod 0 f = (0::'a :: {idom} poly)"
"pseudo_mod f 0 = f"
unfolding pseudo_mod_def snd_pseudo_mod_smult_left[of 0 _ f,simplified]
unfolding pseudo_divmod_def by auto
lemma prod_combine:
assumes "j ≤ i"
shows "f i * (∏l←[j..<i]. (f l :: 'a::comm_monoid_mult)) = prod_list (map f [j..<Suc i])"
proof(subst prod_list_map_remove1[of i "[j..<Suc i]" f],goal_cases)
case 2
have "remove1 i ([j..<i] @ [i]) = [j..<i]" by (simp add: remove1_append)
thus ?case by auto
qed (insert assms, auto)
lemma prod_list_minus_1_exp: "prod_list (map (λ i. (-1)^(f i)) xs)
= (-1)^(sum_list (map f xs))"
by (induct xs, auto simp: power_add)
lemma minus_1_power_even: "(- (1 :: 'b :: comm_ring_1))^ k = (if even k then 1 else (-1))"
by auto
lemma minus_1_even_eqI: assumes "even k = even l" shows
"(- (1 :: 'b :: comm_ring_1))^k = (- 1)^l"
unfolding minus_1_power_even assms by auto
lemma (in comm_monoid_mult) prod_list_multf:
"(∏x←xs. f x * g x) = prod_list (map f xs) * prod_list (map g xs)"
by (induct xs) (simp_all add: algebra_simps)
lemma inverse_prod_list: "inverse (prod_list xs) = prod_list (map inverse (xs :: 'a :: field list))"
by (induct xs, auto)
definition pow_int :: "'a :: field ⇒ int ⇒ 'a" where
"pow_int x e = (if e < 0 then 1 / (x ^ (nat (-e))) else x ^ (nat e))"
lemma pow_int_0[simp]: "pow_int x 0 = 1" unfolding pow_int_def by auto
lemma pow_int_1[simp]: "pow_int x 1 = x" unfolding pow_int_def by auto
lemma exp_pow_int: "x ^ n = pow_int x n"
unfolding pow_int_def by auto
lemma pow_int_add: assumes x: "x ≠ 0" shows "pow_int x (a + b) = pow_int x a * pow_int x b"
proof -
have *:
"¬ a + b < 0 ⟹ a < 0 ⟹ nat b = nat (a + b) + nat (-a)"
"¬ a + b < 0 ⟹ b < 0 ⟹ nat a = nat (a + b) + nat (-b)"
"a + b < 0 ⟹ ¬ a < 0 ⟹ nat (-b) = nat a + nat (-a -b) "
"a + b < 0 ⟹ ¬ b < 0 ⟹ nat (-a) = nat b + nat (-a -b) "
by auto
have pow_eq: "l = m ⟹ (x ^ l = x ^ m)" for l m by auto
from x show ?thesis unfolding pow_int_def
by (auto split: if_splits simp: power_add[symmetric] simp: * intro!: pow_eq, auto simp: power_add)
qed
lemma pow_int_mult: "pow_int (x * y) a = pow_int x a * pow_int y a"
unfolding pow_int_def by (cases "a < 0", auto simp: power_mult_distrib)
lemma pow_int_base_1[simp]: "pow_int 1 a = 1"
unfolding pow_int_def by (cases "a < 0", auto)
lemma pow_int_divide: "a / pow_int x b = a * pow_int x (-b)"
unfolding pow_int_def by (cases b rule: linorder_cases[of _ 0], auto)
lemma divide_prod_assoc: "x / (y * z :: 'a :: field) = x / y / z" by (simp add: field_simps)
lemma minus_1_inverse_pow[simp]: "x / (-1)^n = (x :: 'a :: field) * (-1)^n"
by (simp add: minus_1_power_even)
definition subresultant_mat :: "nat ⇒ 'a :: comm_ring_1 poly ⇒ 'a poly ⇒ 'a poly mat" where
"subresultant_mat J F G = (let
dg = degree G; df = degree F; f = coeff_int F; g = coeff_int G; n = (df - J) + (dg - J)
in mat n n (λ (i,j). if j < dg - J then
if i = n - 1 then monom 1 (dg - J - 1 - j) * F else [: f (df - int i + int j) :]
else let jj = j - (dg - J) in
if i = n - 1 then monom 1 (df - J - 1 - jj) * G else [: g (dg - int i + int jj) :]))"
lemma subresultant_mat_dim[simp]:
fixes j p q
defines "S ≡ subresultant_mat j p q"
shows "dim_row S = (degree p - j) + (degree q - j)" and "dim_col S = (degree p - j) + (degree q - j)"
unfolding S_def subresultant_mat_def Let_def by auto
definition subresultant'_mat :: "nat ⇒ nat ⇒ 'a :: comm_ring_1 poly ⇒ 'a poly ⇒ 'a mat" where
"subresultant'_mat J l F G = (let
γ = degree G; φ = degree F; f = coeff_int F; g = coeff_int G; n = (φ - J) + (γ - J)
in mat n n (λ (i,j). if j < γ - J then
if i = n - 1 then (f (l - int (γ - J - 1) + int j)) else (f (φ - int i + int j))
else let jj = j - (γ - J) in
if i = n - 1 then (g (l - int (φ - J - 1) + int jj)) else (g (γ - int i + int jj))))"
lemma subresultant_index_mat:
fixes F G
assumes i: "i < (degree F - J) + (degree G - J)" and j: "j < (degree F - J) + (degree G - J)"
shows "subresultant_mat J F G $$ (i,j) =
(if j < degree G - J then
if i = (degree F - J) + (degree G - J) - 1 then monom 1 (degree G - J - 1 - j) * F else ([: coeff_int F ( degree F - int i + int j) :])
else let jj = j - (degree G - J) in
if i = (degree F - J) + (degree G - J) - 1 then monom 1 ( degree F - J - 1 - jj) * G else ([: coeff_int G (degree G - int i + int jj) :]))"
unfolding subresultant_mat_def Let_def
unfolding index_mat(1)[OF i j] split by auto
definition subresultant :: "nat ⇒ 'a :: comm_ring_1 poly ⇒ 'a poly ⇒ 'a poly" where
"subresultant J F G = det (subresultant_mat J F G)"
lemma subresultant_smult_left: assumes "(c :: 'a :: {comm_ring_1, semiring_no_zero_divisors}) ≠ 0"
shows "subresultant J (smult c f) g = smult (c ^ (degree g - J)) (subresultant J f g)"
proof -
let ?df = "degree f"
let ?dg = "degree g"
let ?n = "(?df - J) + (?dg - J)"
let ?m = "?dg - J"
let ?M = "mat ?n ?n (λ (i,j). if i = j then if i < ?m then [:c:] else 1 else 0)"
from ‹c ≠ 0› have deg: "degree (smult c f) = ?df" by simp
let ?S = "subresultant_mat J f g"
let ?cS = "subresultant_mat J (smult c f) g"
have dim: "dim_row ?S = ?n" "dim_col ?S = ?n" "dim_row ?cS = ?n" "dim_col ?cS = ?n" using deg by auto
hence C: "?S ∈ carrier_mat ?n ?n" "?cS ∈ carrier_mat ?n ?n" "?M ∈ carrier_mat ?n ?n" by auto
have dim': "dim_row (?S * ?M) = ?n" "dim_col (?S * ?M) = ?n" using dim (1,2) by simp_all
define S where "S = ?S"
have "?cS = ?S * ?M"
proof (rule eq_matI, unfold dim' dim)
fix i j
assume ij: "i < ?n" "j < ?n"
have "(?S * ?M) $$ (i,j) = row ?S i ∙ col ?M j"
by (rule index_mult_mat, insert ij dim, auto)
also have "… = (∑k = 0..<?n. row S i $ k * col ?M j $ k)" unfolding scalar_prod_def S_def[symmetric]
by simp
also have "… = (∑k = 0..<?n. S $$ (i,k) * ?M $$ (k,j))"
by (rule sum.cong, insert ij, auto simp: S_def)
also have "… = S $$ (i,j) * ?M $$ (j,j) + sum (λ k. S $$ (i,k) * ?M $$ (k,j)) ({0..<?n} - {j})"
by (rule sum.remove, insert ij, auto)
also have "… = S $$ (i,j) * ?M $$ (j,j)"
by (subst sum.neutral, insert ij, auto)
also have "… = ?cS $$ (i,j)" unfolding subresultant_index_mat[OF ij] S_def
by (subst subresultant_index_mat, unfold deg, insert ij, auto)
finally show "?cS $$ (i,j) = (?S * ?M) $$ (i,j)" by simp
qed auto
from arg_cong[OF this, of det] det_mult[OF C(1) C(3)]
have "subresultant J (smult c f) g = subresultant J f g * det ?M"
unfolding subresultant_def by auto
also have "det ?M = [:c ^ ?m :]"
proof (subst det_upper_triangular[OF _ C(3)])
show "upper_triangular ?M"
by (rule upper_triangularI, auto)
have "prod_list (diag_mat ?M) = (∏k = 0..<?n. (?M $$ (k,k)))"
unfolding prod_list_diag_prod by simp
also have "… = (∏k = 0..<?m. ?M $$ (k,k)) * (∏k = ?m..<?n. ?M $$ (k,k))"
by (subst prod.union_disjoint[symmetric], (auto)[3], rule prod.cong, auto)
also have "(∏k = 0..<?m. ?M $$ (k,k)) = (∏k = 0..<?m. [: c :])"
by (rule prod.cong, auto)
also have "(∏k = 0..<?m. [: c :]) = [: c :] ^ ?m" by simp
also have "(∏k = ?m..<?n. ?M $$ (k,k)) = (∏k = ?m..<?n. 1)"
by (rule prod.cong, auto)
finally show "prod_list (diag_mat ?M) = [: c^?m :]" unfolding poly_const_pow by simp
qed
finally show ?thesis by simp
qed
lemma subresultant_swap:
shows "subresultant J f g = smult ((- 1) ^ ((degree f - J) * (degree g - J))) (subresultant J g f)"
proof -
let ?A = "subresultant_mat J f g"
let ?k = "degree f - J"
let ?n = "degree g - J"
have nk: "?n + ?k = ?k + ?n" by simp
have change: "j < ?k + ?n ⟹ ((if j < ?k then j + ?n else j - ?k) < ?n)
= (¬ (j < ?k))" for j by auto
have "subresultant J f g = det ?A" unfolding subresultant_def by simp
also have "… = (-1)^(?k * ?n) * det (mat (?k + ?n) (?k + ?n) (λ (i,j).
?A $$ (i,(if j < ?k then j + ?n else j - ?k))))" (is "_ = _ * det ?B")
by (rule det_swap_cols, auto simp: subresultant_mat_def Let_def)
also have "?B = subresultant_mat J g f"
unfolding subresultant_mat_def Let_def
by (rule eq_matI, unfold dim_row_mat dim_col_mat nk index_mat split,
subst index_mat, (auto)[2], unfold split, subst change, force,
unfold if_conn, rule if_cong[OF refl if_cong if_cong], auto)
also have "det … = subresultant J g f" unfolding subresultant_def ..
also have "(-1)^(?k * ?n) * … = [: (-1)^(?k * ?n) :] * … " by (unfold hom_distribs, simp)
also have "… = smult ((-1)^(?k * ?n)) (subresultant J g f)" by simp
finally show ?thesis .
qed
lemma subresultant_smult_right:assumes "(c :: 'a :: {comm_ring_1, semiring_no_zero_divisors}) ≠ 0"
shows "subresultant J f (smult c g) = smult (c ^ (degree f - J)) (subresultant J f g)"
unfolding subresultant_swap[of _ f] subresultant_smult_left[OF assms]
degree_smult_eq using assms by (simp add: ac_simps)
lemma coeff_subresultant: "coeff (subresultant J F G) l =
(if degree F - J + (degree G - J) = 0 ∧ l ≠ 0 then 0 else det (subresultant'_mat J l F G))"
proof (cases "degree F - J + (degree G - J) = 0")
case True
show ?thesis unfolding subresultant_def subresultant_mat_def subresultant'_mat_def Let_def True
by simp
next
case False
let ?n = "degree F - J + (degree G - J)"
define n where "n = ?n"
from False have n: "n ≠ 0" unfolding n_def by auto
hence id: "{0..<n} = insert (n - 1) {0..< (n - 1)}" by (cases n, auto)
have idn: "(x = x) = True" for x :: nat by simp
let ?M = "subresultant_mat J F G"
define M where "M = ?M"
let ?L = "subresultant'_mat J l F G"
define L where "L = ?L"
{
fix p
assume p: "p permutes {0..<n}"
from n p have n1: "n - 1 < n" "p (n - 1) < n" by auto
have "coeff_int (∏i = 0..<n. M $$ (i, p i)) l =
(∏i = 0 ..< (n - 1). coeff_int (M $$ (i, p i)) 0) * coeff_int (M $$ (n - 1, p (n - 1))) l"
unfolding id
proof (rule coeff_int_prod_const, (auto)[2])
fix i
assume "i ∈ {0 ..< n - 1}"
with p have i: "i ≠ n - 1" and "i < n" "p i < n" by (auto simp: n_def)
note id = subresultant_index_mat[OF this(2-3)[unfolded n_def], folded M_def n_def]
show "degree (M $$ (i, p i)) = 0" unfolding id Let_def using i
by (simp split: if_splits)
qed
also have "(∏i = 0 ..< (n - 1). coeff_int (M $$ (i, p i)) 0)
= (∏i = 0 ..< (n - 1). L $$ (i, p i))"
proof (rule prod.cong[OF refl])
fix i
assume "i ∈ {0 ..< n - 1}"
with p have i: "i ≠ n - 1" and ii: "i < n" "p i < n" by (auto simp: n_def)
note id = subresultant_index_mat[OF this(2-3)[unfolded n_def], folded M_def n_def]
note id' = L_def[unfolded subresultant'_mat_def Let_def, folded n_def] index_mat[OF ii]
show "coeff_int (M $$ (i, p i)) 0 = L $$ (i, p i)"
unfolding id id' split using i proof (simp add: if_splits Let_def)
qed
qed
also have "coeff_int (M $$ (n - 1, p (n - 1))) l =
(if p (n - 1) < degree G - J then
coeff_int (monom 1 (degree G - J - 1 - p (n - 1)) * F) l
else coeff_int (monom 1 (degree F - J - 1 - (p (n - 1) - (degree G - J))) * G) l)"
using subresultant_index_mat[OF n1[unfolded n_def], folded M_def n_def, unfolded idn if_True Let_def]
by simp
also have "… = (if p (n - 1) < degree G - J
then coeff_int F (int l - int (degree G - J - 1 - p (n - 1)))
else coeff_int G (int l - int (degree F - J - 1 - (p (n - 1) - (degree G - J)))))"
unfolding coeff_int_monom_mult by simp
also have "… = (if p (n - 1) < degree G - J
then coeff_int F (int l - int (degree G - J - 1) + p (n - 1))
else coeff_int G (int l - int (degree F - J - 1) + (p (n - 1) - (degree G - J))))"
proof (cases "p (n - 1) < degree G - J")
case True
hence "int (degree G - J - 1 - p (n - 1)) = int (degree G - J - 1) - p (n - 1)" by simp
hence id: "int l - int (degree G - J - 1 - p (n - 1)) = int l - int (degree G - J - 1) + p (n - 1)" by simp
show ?thesis using True unfolding id by simp
next
case False
from n1 False have "degree F - J - 1 ≥ p (n - 1) - (degree G - J)"
unfolding n_def by linarith
hence "int (degree F - J - 1 - (p (n - 1) - (degree G - J))) = int (degree F - J - 1) - (p (n - 1) - (degree G - J))"
by linarith
hence id: "int l - int (degree F - J - 1 - (p (n - 1) - (degree G - J))) =
int l - int (degree F - J - 1) + (p (n - 1) - (degree G - J))" by simp
show ?thesis unfolding id using False by simp
qed
also have "… = L $$ (n - 1, p (n - 1))"
unfolding L_def subresultant'_mat_def Let_def n_def[symmetric] using n1 by simp
also have "(∏i = 0..<n - 1. L $$ (i, p i)) * … = (∏i = 0..<n. L $$ (i, p i))"
unfolding id by simp
finally have "coeff_int (∏i = 0..<n. M $$ (i, p i)) (int l) = (∏i = 0..<n. L $$ (i, p i))" .
} note * = this
have "coeff_int (subresultant J F G) l =
(∑p∈{p. p permutes {0..<n}}. signof p * coeff_int (∏i = 0..<n. M $$ (i, p i)) l)"
unfolding subresultant_def det_def subresultant_mat_dim idn if_True n_def[symmetric] M_def
coeff_int_sum coeff_int_signof_mult by simp
also have "… = (∑p∈{p. p permutes {0..<n}}. signof p * (∏i = 0..<n. L $$ (i, p i)))"
by (rule sum.cong[OF refl], insert *, simp)
also have "… = det L"
proof -
have id: "dim_row (subresultant'_mat J l F G) = n"
"dim_col (subresultant'_mat J l F G) = n" unfolding subresultant'_mat_def Let_def n_def
by auto
show ?thesis unfolding det_def L_def id by simp
qed
finally show ?thesis unfolding L_def coeff_int_def using False by auto
qed
lemma subresultant'_zero_ge: assumes "(degree f - J) + (degree g - J) ≠ 0" and "k ≥ degree f + (degree g - J)"
shows "det (subresultant'_mat J k f g) = 0"
proof -
obtain dg where dg: "degree g - J = dg" by simp
obtain df where df: "degree f - J = df" by simp
obtain ddf where ddf: "degree f = ddf" by simp
note * = assms(2)[unfolded ddf dg] assms(1)
define M where "M = (λ i j. if j < dg
then coeff_int f (degree f - int i + int j)
else coeff_int g (degree g - int i + int (j - dg)))"
let ?M = "subresultant'_mat J k f g"
have M: "det ?M = det (mat (df + dg) (df + dg)
(λ(i, j).
if i = df + dg - 1 then
if j < dg
then coeff_int f (int k - int (dg - 1) + int j)
else coeff_int g (int k - int (df - 1) + int (j - dg))
else M i j))" (is "_ = det ?N")
unfolding subresultant'_mat_def Let_def M_def
by (rule arg_cong[of _ _ det], rule eq_matI, auto simp: df dg)
also have "?N = mat (df + dg) (df + dg)
(λ(i, j).
if i = df + dg - 1 then 0
else M i j)"
by (rule cong_mat[OF refl refl], unfold split, rule if_cong[OF refl _ refl],
auto simp add: coeff_int_def df dg ddf intro!: coeff_eq_0, insert *(1),
unfold ddf[symmetric] dg[symmetric] df[symmetric], linarith+)
also have "… = mat⇩r (df + dg) (df + dg) (λi. if i = df + dg - 1 then 0⇩v (df + dg) else
vec (df + dg) (λ j. M i j))"
by (rule eq_matI, auto)
also have "det … = 0"
by (rule det_row_0, insert *, auto simp: df[symmetric] dg[symmetric] ddf[symmetric])
finally show ?thesis .
qed
lemma subresultant'_zero_lt: assumes
J: "J ≤ degree f" "J ≤ degree g" "J < k"
and k: "k < degree f + (degree g - J)"
shows "det (subresultant'_mat J k f g) = 0"
proof -
obtain dg where dg: "dg = degree g - J" by simp
obtain df where df: "df = degree f - J" by simp
note * = assms[folded df dg]
define M where "M = (λ i j. if j < dg
then coeff_int f (degree f - int i + int j)
else coeff_int g (degree g - int i + int (j - dg)))"
define N where "N = (λ j. if j < dg
then coeff_int f (int k - int (dg - 1) + int j)
else coeff_int g (int k - int (df - 1) + int (j - dg)))"
let ?M = "subresultant'_mat J k f g"
have M: "?M = mat (df + dg) (df + dg)
(λ(i, j).
if i = df + dg - 1 then N j
else M i j)"
unfolding subresultant'_mat_def Let_def
by (rule eq_matI, auto simp: df dg M_def N_def)
also have "… = mat (df + dg) (df + dg)
(λ(i, j).
if i = df + dg - 1 then N j
else if i = degree f + dg - 1 - k then N j else M i j)" (is "_ = ?N")
unfolding N_def
by (rule cong_mat[OF refl refl], unfold split, rule if_cong[OF refl refl], unfold M_def N_def,
insert J k, auto simp: df dg intro!: arg_cong[of _ _ "coeff_int _"])
finally have id: "?M = ?N" .
have deg: "degree f + dg - 1 - k < df + dg" "df + dg - 1 < df + dg"
using k J unfolding df dg by auto
have id: "row ?M (degree f + dg - 1 - k) = row ?M (df + dg - 1)"
unfolding arg_cong[OF id, of row]
by (rule eq_vecI, insert deg, auto)
show ?thesis
by (rule det_identical_rows[OF _ _ _ _ id, of "df + dg"], insert deg assms,
auto simp: subresultant'_mat_def Let_def df dg)
qed
lemma subresultant'_mat_sylvester_mat: "transpose_mat (subresultant'_mat 0 0 f g) = sylvester_mat f g"
proof -
obtain dg where dg: "degree g = dg" by simp
obtain df where df: "degree f = df" by simp
let ?M = "transpose_mat (subresultant'_mat 0 0 f g)"
let ?n = "degree f + degree g"
have dim: "dim_row ?M = ?n" "dim_col ?M = ?n" by (auto simp: subresultant'_mat_def Let_def)
show ?thesis
proof (rule eq_matI, unfold sylvester_mat_dim dim df dg, goal_cases)
case ij: (1 i j)
have "?M $$ (i,j) = (if i < dg
then if j = df + dg - 1
then coeff_int f (- int (dg - 1) + int i)
else coeff_int f (int df - int j + int i)
else if j = df + dg - 1
then coeff_int g (- int (df - 1) + int (i - dg))
else coeff_int g (int dg - int j + int (i - dg)))"
using ij unfolding subresultant'_mat_def Let_def by (simp add: if_splits df dg)
also have "… = (if i < dg
then coeff_int f (int df - int j + int i)
else coeff_int g (int dg - int j + int (i - dg)))"
proof -
have cong: "(b ⟹ x = z) ⟹ (¬ b ⟹ y = z) ⟹ (if b then coeff_int f x else coeff_int f y) = coeff_int f z"
for b x y z and f :: "'a poly" by auto
show ?thesis
by (rule if_cong[OF refl cong cong], insert ij, auto)
qed
also have "… = sylvester_mat f g $$ (i,j)"
proof -
have *: "i ≤ j ⟹ j - i ≤ df ⟹ nat (int df - int j + int i) = df - (j - i)" for j i df
by simp
show ?thesis unfolding sylvester_index_mat[OF ij[folded df dg]] df dg
proof (rule if_cong[OF refl])
assume i: "i < dg"
have "int df - int j + int i < 0 ⟶ ¬ j - i ≤ df" by auto
thus "coeff_int f (int df - int j + int i) = (if i ≤ j ∧ j - i ≤ df then coeff f (df + i - j) else 0)"
using i ij by (simp add: coeff_int_def *, intro impI coeff_eq_0[of f, unfolded df], linarith)
next
assume i: "¬ i < dg"
hence **: "i - dg ≤ j ⟹ dg - (j + dg - i) = i - j" using ij by linarith
have "int dg - int j + int (i - dg) < 0 ⟶ ¬ j ≤ i" by auto
thus "coeff_int g (int dg - int j + int (i - dg)) = (if i - dg ≤ j ∧ j ≤ i then coeff g (i - j) else 0)"
using ij i by (simp add: coeff_int_def * **, intro impI coeff_eq_0[of g, unfolded dg], linarith)
qed
qed
finally show ?case .
qed auto
qed
lemma coeff_subresultant_0_0_resultant: "coeff (subresultant 0 f g) 0 = resultant f g"
proof -
let ?M = "transpose_mat (subresultant'_mat 0 0 f g)"
have "det (subresultant'_mat 0 0 f g) = det ?M"
by (subst det_transpose, auto simp: subresultant'_mat_def Let_def)
also have "?M = sylvester_mat f g"
by (rule subresultant'_mat_sylvester_mat)
finally show ?thesis by (simp add: coeff_subresultant resultant_def)
qed
lemma subresultant_zero_ge: assumes "k ≥ degree f + (degree g - J)"
and "(degree f - J) + (degree g - J) ≠ 0"
shows "coeff (subresultant J f g) k = 0"
unfolding coeff_subresultant
by (subst subresultant'_zero_ge[OF assms(2,1)], simp)
lemma subresultant_zero_lt: assumes "k < degree f + (degree g - J)"
and "J ≤ degree f" "J ≤ degree g" "J < k"
shows "coeff (subresultant J f g) k = 0"
unfolding coeff_subresultant
by (subst subresultant'_zero_lt[OF assms(2,3,4,1)], simp)
lemma subresultant_resultant: "subresultant 0 f g = [: resultant f g :]"
proof (cases "degree f + degree g = 0")
case True
thus ?thesis unfolding subresultant_def subresultant_mat_def resultant_def Let_def
sylvester_mat_def sylvester_mat_sub_def
by simp
next
case 0: False
show ?thesis
proof (rule poly_eqI)
fix k
show "coeff (subresultant 0 f g) k = coeff [:resultant f g:] k"
proof (cases "k = 0")
case True
thus ?thesis using coeff_subresultant_0_0_resultant[of f g] by auto
next
case False
hence "0 < k ∧ k < degree f + degree g ∨ k ≥ degree f + degree g" by auto
thus ?thesis using subresultant_zero_ge[of f g 0 k] 0
subresultant_zero_lt[of k f g 0] 0 False by (cases k, auto)
qed
qed
qed
lemma (in inj_comm_ring_hom) subresultant_hom:
"map_poly hom (subresultant J f g) = subresultant J (map_poly hom f) (map_poly hom g)"
proof -
note d = subresultant_mat_def Let_def
interpret p: map_poly_inj_comm_ring_hom hom..
show ?thesis unfolding subresultant_def unfolding p.hom_det[symmetric]
proof (rule arg_cong[of _ _ det])
show "p.mat_hom (subresultant_mat J f g) =
subresultant_mat J (map_poly hom f) (map_poly hom g)"
proof (rule eq_matI, goal_cases)
case (1 i j)
hence ij: "i < degree f - J + (degree g - J)" "j < degree f - J + (degree g - J)"
unfolding d degree_map_poly by auto
show ?case
by (auto simp add: coeff_int_def d map_mat_def index_mat(1)[OF ij] hom_distribs)
qed (auto simp: d)
qed
qed
text ‹We now derive properties of the resultant via the connection to subresultants.›
lemma resultant_smult_left: assumes "(c :: 'a :: idom) ≠ 0"
shows "resultant (smult c f) g = c ^ degree g * resultant f g"
unfolding coeff_subresultant_0_0_resultant[symmetric] subresultant_smult_left[OF assms] coeff_smult
by simp
lemma resultant_smult_right: assumes "(c :: 'a :: idom) ≠ 0"
shows "resultant f (smult c g) = c ^ degree f * resultant f g"
unfolding coeff_subresultant_0_0_resultant[symmetric] subresultant_smult_right[OF assms] coeff_smult
by simp
lemma resultant_swap: "resultant f g = (-1)^(degree f * degree g) * (resultant g f)"
unfolding coeff_subresultant_0_0_resultant[symmetric]
unfolding arg_cong[OF subresultant_swap[of 0 f g], of "λ x. coeff x 0"] coeff_smult by simp
text ‹The following equations are taken from Brown-Traub ``On Euclid's Algorithm and
the Theory of Subresultant'' (BT)›
lemma fixes F B G H :: "'a :: idom poly" and J :: nat
defines df: "df ≡ degree F"
and dg: "dg ≡ degree G"
and dh: "dh ≡ degree H"
and db: "db ≡ degree B"
defines
n: "n ≡ (df - J) + (dg - J)"
and f: "f ≡ coeff_int F"
and b: "b ≡ coeff_int B"
and g: "g ≡ coeff_int G"
and h: "h ≡ coeff_int H"
assumes FGH: "F + B * G = H"
and dfg: "df ≥ dg"
and choice: "dg > dh ∨ H = 0 ∧ F ≠ 0 ∧ G ≠ 0"
shows BT_eq_18: "subresultant J F G = smult ((-1)^((df - J) * (dg - J))) (det (mat n n
(λ (i,j).
if j < df - J
then if i = n - 1 then monom 1 ((df - J) - 1 - j) * G
else [:g (int dg - int i + int j):]
else if i = n - 1 then monom 1 ((dg - J) - 1 - (j - (df - J))) * H
else [:h (int df - int i + int (j - (df - J))):])))"
(is "_ = smult ?m1 ?right")
and BT_eq_19: "dh ≤ J ⟹ J < dg ⟹ subresultant J F G = smult (
(-1)^((df - J) * (dg - J)) * lead_coeff G ^ (df - J) * coeff H J ^ (dg - J - 1)) H"
(is "_ ⟹ _ ⟹ _ = smult (_ * ?G * ?H) H")
and BT_lemma_1_12: "J < dh ⟹ subresultant J F G = smult (
(-1)^((df - J) * (dg - J)) * lead_coeff G ^ (df - dh)) (subresultant J G H)"
and BT_lemma_1_13': "J = dh ⟹ dg > dh ∨ H ≠ 0 ⟹ subresultant dh F G = smult (
(-1)^((df - dh) * (dg - dh)) * lead_coeff G ^ (df - dh) * lead_coeff H ^ (dg - dh - 1)) H"
and BT_lemma_1_14: "dh < J ⟹ J < dg - 1 ⟹ subresultant J F G = 0"
and BT_lemma_1_15': "J = dg - 1 ⟹ dg > dh ∨ H ≠ 0 ⟹ subresultant (dg - 1) F G = smult (
(-1)^(df - dg + 1) * lead_coeff G ^ (df - dg + 1)) H"
proof -
define dfj where "dfj = df - J"
define dgj where "dgj = dg - J"
note d = df dg dh db
have F0: "F ≠ 0" using dfg choice df by auto
have G0: "G ≠ 0" using choice dg by auto
have dgh: "dg ≥ dh" using choice unfolding dh by auto
have B0: "B ≠ 0" using FGH dfg dgh choice F0 G0 unfolding d by auto
have dfh: "df ≥ dh" using dfg dgh by auto
have "df = degree (B * G)"
proof (cases "H = 0")
case False
with choice dfg have dfh: "df > dh" by auto
show ?thesis using dfh[folded arg_cong[OF FGH, of degree, folded dh]] choice
unfolding df by (metis ‹degree (F + B * G) < df› degree_add_eq_left degree_add_eq_right df nat_neq_iff)
next
case True
have "F = - B * G" using arg_cong[OF FGH[unfolded True], of "λ x. x - B * G"] by auto
thus ?thesis using F0 G0 B0 unfolding df by simp
qed
hence dfbg: "df = db + dg" using degree_mult_eq[OF B0 G0] by (simp add: d)
hence dbfg: "db = df - dg" by simp
let ?dfj = "df - J"
let ?dgj = "dg - J"
have norm: "?dgj + ?dfj = ?dfj + ?dgj" by simp
let ?bij = "λ i j. b (db - int i + int (j - dfj))"
let ?M = "mat n n (λ (i,j). if i = j then 1 else if j < dfj then 0 else if i < j
then [: ?bij i j :] else 0)"
let ?GF = "λ i j.
if j < dfj
then if i = n - 1 then monom 1 (dfj - 1 - j) * G
else [:g (int dg - int i + int j):]
else if i = n - 1 then monom 1 (dgj - 1 - (j - dfj)) * F
else [:f (int df - int i + int (j - dfj)):]"
let ?G_F = "mat n n (λ (i,j). ?GF i j)"
let ?GH = "λ i j.
if j < dfj
then if i = n - 1 then monom 1 (dfj - 1 - j) * G
else [:g (int dg - int i + int j):]
else if i = n - 1 then monom 1 (dgj - 1 - (j - dfj)) * H
else [:h (int df - int i + int (j - dfj)):]"
let ?G_H = "mat n n (λ (i,j). ?GH i j)"
have hfg: "h i = f i + coeff_int (B * G) i" for i
unfolding FGH[symmetric] f g h unfolding coeff_int_def by simp
have dM1: "det ?M = 1"
by (subst det_upper_triangular, (auto)[2], subst prod_list_diag_prod, auto)
have "subresultant J F G = smult ?m1 (subresultant J G F)"
unfolding subresultant_swap[of _ F] d by simp
also have "subresultant J G F = det ?G_F"
unfolding subresultant_def n norm subresultant_mat_def g f Let_def d[symmetric] dfj_def dgj_def by simp
also have "… = det (?G_F * ?M)"
by (subst det_mult[of _ n], unfold dM1, auto)
also have "?G_F * ?M = ?G_H"
proof (rule eq_matI, unfold dim_col_mat dim_row_mat)
fix i j
assume i: "i < n" and j: "j < n"
have "(?G_F * ?M) $$ (i,j) = row (?G_F) i ∙ col ?M j"
using i j by simp
also have "… = ?GH i j"
proof (cases "j < dfj")
case True
have id: "col ?M j = unit_vec n j"
by (rule eq_vecI, insert True i j, auto)
show ?thesis unfolding id using True i j by simp
next
case False
define d where "d = j - dfj"
from False have jd: "j = d + dfj" unfolding d_def by auto
hence idj: "{0 ..< j} = {0 ..< dfj} ∪ {dfj ..< dfj + d}" by auto
let ?H = "(if i = n - 1 then monom 1 (dgj - Suc d) * H else [:h (int df - int i + int d):])"
have idr: "?GH i j = ?H" unfolding d_def using jd by auto
let ?bi = "λ i. b (db - int i + int d)"
let ?m = "λ i. if i = j then 1 else if i < j then [:?bij i j:] else 0"
let ?P = "λ k. (?GF i k * ?m k)"
let ?Q = "λ k. ?GF i k * [: ?bi k :]"
let ?G = "λ k. if i = n - 1 then monom 1 (dfj - 1 - k) * G else [:g (int dg - int i + int k):]"
let ?Gb = "λ k. ?G k * [:?bi k:]"
let ?off = "- (int db - int dfj + 1 + int d)"
have off0: "?off ≥ 0" using False dfg j unfolding dfj_def d_def dbfg n by simp
from nat_0_le[OF this]
obtain off where off: "int off = ?off" by blast
have "int off ≤ int dfj" unfolding off by auto
hence "off ≤ dfj" by simp
hence split1: "{0 ..< dfj} = {0 ..< off} ∪ {off ..< dfj}" by auto
have "int off + Suc db ≤ dfj" unfolding off by auto
hence split2: "{off ..< dfj} = {off .. off + db} ∪ {off + Suc db ..< dfj} " by auto
let ?g_b = "λk. (if i = n - 1 then monom 1 k * G else [:g (int dg - int i + int (dfj - Suc k)):]) *
[:b (k - int off):]"
let ?gb = "λk. (if i = n - 1 then monom 1 (k + off) * G else [:g (int dg - int i + int (dfj - Suc k - off)):]) *
[:coeff B k:]"
let ?F = "λ k. if i = n - 1 then monom 1 (dgj - 1 - (k - dfj)) * F
else [:f (int df - int i + int (k - dfj)):]"
let ?Fb = "λ k. ?F k * [:?bi k:]"
let ?Pj = "if i = n - 1 then monom 1 (dgj - Suc d) * F else [:f (int df - int i + int d):]"
from False have id: "col ?M j = vec n ?m"
using j i by (intro eq_vecI, auto)
have "row ?G_F i ∙ col ?M j = sum ?P {0 ..< n}"
using i j unfolding id by (simp add: scalar_prod_def)
also have "{0 ..< n} = {0 ..< j} ∪ {j} ∪ {Suc j ..< n}" using j by auto
also have "sum ?P … = sum ?P {0 ..< j} + ?P j + sum ?P {Suc j ..< n}"
by (simp add: sum.union_disjoint)
also have "sum ?P {Suc j ..< n} = 0" by (rule sum.neutral, auto)
also have "?P j = ?Pj"
unfolding d_def using jd by simp
also have "sum ?P {0 ..< j} = sum ?Q {0 ..< j}"
by (rule sum.cong[OF refl], unfold d_def, insert jd, auto)
also have "sum ?Q {0 ..< j} = sum ?Q {0 ..< dfj} + sum ?Q {dfj ..< dfj+d}" unfolding idj
by (simp add: sum.union_disjoint)
also have "sum ?Q {0 ..< dfj} = sum ?Gb {0 ..< dfj}"
by (rule sum.cong, auto)
also have "sum ?Q {dfj ..< dfj+d} = sum ?Fb {dfj ..< dfj+d}"
by (rule sum.cong, auto)
also have "… = 0"
proof (rule sum.neutral, intro ballI)
fix k
assume k: "k ∈ {dfj ..< dfj+d}"
hence k: "db + d < k" using k j False unfolding n db[symmetric] dfbg dfj_def d_def by auto
let ?k = "(int db - int k + int d)"
have "?k < 0" using k by auto
hence "b ?k = 0" unfolding b by (intro coeff_int_eq_0, auto)
thus "?Fb k = 0" by simp
qed
also have "sum ?Gb {0 ..< dfj} = sum ?g_b {0 ..< dfj}"
proof (rule sum.reindex_cong[of "λ k. dfj - Suc k"], (auto simp: inj_on_def off)[2], goal_cases)
case (1 k)
hence "k = dfj - (Suc (dfj - Suc k))" and "(dfj - Suc k) ∈ {0..<dfj}" by auto
thus ?case by blast
next
case (2 k)
hence [simp]: "dfj - Suc (dfj - Suc k) = k"
"int db - int (dfj - Suc k) + int d = int k - off" by (auto simp: off)
show ?case by auto
qed
also have "… = sum ?g_b {0 ..< off} + sum ?g_b {off ..< dfj}" unfolding split1
by (simp add: sum.union_disjoint)
also have "sum ?g_b {0 ..< off} = 0"
by (rule sum.neutral, intro ballI, auto simp: b coeff_int_def)
also have "sum ?g_b {off ..< dfj} = sum ?g_b {off .. off + db} + sum ?g_b {off + Suc db ..< dfj}"
unfolding split2 by (rule sum.union_disjoint, auto)
also have "sum ?g_b {off + Suc db ..< dfj} = 0"
proof (rule sum.neutral, intro ballI, goal_cases)
case (1 k)
hence "b (int k - int off) = 0" unfolding b db
by (intro coeff_int_eq_0, auto)
thus ?case by simp
qed
also have "sum ?g_b {off .. off + db} = sum ?gb {0 .. db}"
using sum.atLeastAtMost_shift_bounds [of ?g_b 0 off db]
by (auto intro: sum.cong simp add: b ac_simps)
finally have id: "row ?G_F i ∙ col ?M j - ?H = ?Pj + sum ?gb {0 .. db} - ?H"
(is "_ = ?E")
by (simp add: ac_simps)
define E where "E = ?E"
let ?b = "coeff B"
have Bsum: "(∑k = 0..db. monom (?b k) k) = B" unfolding db
using atMost_atLeast0 poly_as_sum_of_monoms by auto
have "E = 0"
proof (cases "i = n - 1")
case i_n: False
hence id: "(i = n - 1) = False" by simp
with i have i: "i < n - 1" by auto
let ?ii = "int df - int i + int d"
have "?thesis = ([:f ?ii:] +
(∑k = 0..db.
[:g (int dg - int i + int (dfj - Suc k - off)):] * [:?b k:]) -
[:h ?ii:] = 0)" (is "_ = (?e = 0)") unfolding E_def id if_False by simp
also have "?e = [: f ?ii +
(∑k = 0..db.
g (int dg - int i + int (dfj - Suc k - off)) * ?b k) -
h ?ii:]" (is "_ = [: ?e :]")
proof (rule poly_eqI, goal_cases)
case (1 n)
show ?case unfolding coeff_diff coeff_add coeff_sum coeff_const
by (cases n, auto simp: ac_simps)
qed
also have "[: ?e :] = 0 ⟷ ?e = 0" by simp
also have "?e = (∑k = 0..db. g (int dg - int i + int (dfj - Suc k - off)) * ?b k)
- coeff_int (B * G) ?ii"
unfolding hfg by simp
also have "(B * G) = (∑k = 0..db. monom (?b k) k) * G" unfolding Bsum by simp
also have "… = (∑k = 0..db. monom (?b k) k * G)" by (rule sum_distrib_right)
also have "coeff_int … ?ii = (∑k = 0..db. g (?ii - k) * ?b k)"
unfolding coeff_int_sum coeff_int_monom_mult g by (simp add: ac_simps)
also have "… = (∑k = 0..db. g (int dg - int i + int (dfj - Suc k - off)) * ?b k)"
proof (rule sum.cong[OF refl], goal_cases)
case (1 k)
hence "k ≤ db" by simp
hence id: "int dg - int i + int (dfj - Suc k - off) = ?ii - k"
using False i j off dfg
unfolding dbfg d_def dfj_def n by linarith
show ?case unfolding id ..
qed
finally show ?thesis by simp
next
case True
let ?jj = "dgj - Suc d"
have zero: "int off - (dgj - Suc d) = 0" using dfg False j unfolding off dbfg dfj_def d_def dgj_def n
by linarith
from True have "E = monom 1 ?jj * F + (∑k = 0.. db.
monom 1 (k + off) * G * [: ?b k :]) - monom 1 ?jj * H"
(is "_ = ?A + ?sum - ?mon") unfolding id E_def by simp
also have "?mon = monom 1 ?jj * F + monom 1 ?jj * (B * G)"
unfolding FGH[symmetric] by (simp add: ring_distribs)
also have "?A + ?sum - … = ?sum - (monom 1 ?jj * G) * B" (is "_ = _ - ?GB * B") by simp
also have "?sum = (∑k = 0..db.
(monom 1 ?jj * G) * (monom 1 (k + off - ?jj) * [: ?b k :]))"
proof (rule sum.cong[OF refl], goal_cases)
case (1 k)
let ?one = "1 :: 'a"
have "int off ≥ int ?jj" using j False i True
unfolding off d_def dfj_def dgj_def dfbg n by linarith
hence "k + off = ?jj + (k + off - ?jj)" by linarith
hence id: "monom ?one (k + off) = monom (1 * 1) (?jj + (k + off - ?jj))" by simp
show ?case unfolding id[folded mult_monom] by (simp add: ac_simps)
qed
also have "… = (monom 1 ?jj * G) * (∑k = 0..db. monom 1 (k + off - ?jj) * [:?b k:])"
(is "_ = _ * ?sum")
unfolding sum_distrib_left ..
also have "… - (monom 1 ?jj * G) * B = (monom 1 ?jj * G) * (?sum - B)" by (simp add: ring_distribs)
also have "?sum = (∑k = 0..db. monom 1 k * [:?b k:])"
by (rule sum.cong[OF refl], insert zero, auto)
also have "… = (∑k = 0..db. monom (?b k) k)"
by (rule sum.cong[OF refl], rule poly_eqI, auto)
also have "… = B" unfolding Bsum ..
finally show ?thesis by simp
qed
from id[folded E_def, unfolded this]
show ?thesis using False unfolding d_def by simp
qed
also have "… = ?G_H $$ (i,j)" using i j by simp
finally show "(?G_F * ?M) $$ (i,j) = ?G_H $$ (i,j)" .
qed auto
finally show eq_18: "subresultant J F G = smult ?m1 (det ?G_H)" unfolding dfj_def dgj_def .
{
fix i j
assume ij: "i < j" and j: "j < n"
with dgh have "int dg - int i + int j > int dg" by auto
hence "g (int dg - int i + int j) = 0" unfolding g dg by (intro coeff_int_eq_0, auto)
} note g0 = this
{
assume *: "dh ≤ J" "J < dg"
have n_dfj: "n > dfj" using * unfolding n dfj_def by auto
note eq_18
also have "det ?G_H = prod_list (diag_mat ?G_H)"
proof (rule det_lower_triangular[of n])
fix i j
assume ij: "i < j" and j: "j < n"
from ij j have if_e: "i = n - 1 ⟷ False" by auto
have "?G_H $$ (i,j) = ?GH i j" using ij j by auto
also have "… = 0"
proof (cases "j < dfj")
case True
with True g0[OF ij j] show ?thesis unfolding if_e by simp
next
case False
have "h (int df - int i + int (j - dfj)) = 0" unfolding h
by (rule coeff_int_eq_0, insert False * ij j dfg, unfold dfj_def dh[symmetric], auto)
with False show ?thesis unfolding if_e by auto
qed
finally show "?G_H $$ (i,j) = 0" .
qed auto
also have "… = (∏i = 0..<n. ?GH i i)"
by (subst prod_list_diag_prod, simp)
also have "{0 ..< n} = {0 ..< dfj} ∪ {dfj ..< n}" unfolding n dfj_def by auto
also have "prod (λ i. ?GH i i) … = prod (λ i. ?GH i i) {0 ..< dfj} * prod (λ i. ?GH i i) {dfj ..< n}"
by (simp add: prod.union_disjoint)
also have "prod (λ i. ?GH i i) {0 ..< dfj} = prod (λ i. [: lead_coeff G :]) {0 ..< dfj}"
proof -
show ?thesis
by (rule prod.cong[OF refl], insert n_dfj, auto simp: g coeff_int_def dg)
qed
also have "… = [: (lead_coeff G)^dfj :]" by (simp add: poly_const_pow)
also have "{dfj ..< n} = {dfj ..< n-1} ∪ {n - 1}" using n_dfj by auto
also have "prod (λ i. ?GH i i) … = prod (λ i. ?GH i i) {dfj ..< n-1} * ?GH (n - 1) (n - 1)"
by (simp add: prod.union_disjoint)
also have "?GH (n - 1) (n - 1) = H"
proof -
have "dgj - 1 - (n - 1 - dfj) = 0" using n_dfj unfolding dgj_def dfj_def n by auto
with n_dfj show ?thesis by auto
qed
also have "prod (λ i. ?GH i i) {dfj ..< n-1} = prod (λ i. [:h (int df - dfj):]) {dfj ..< n-1}"
by (rule prod.cong[OF refl], auto intro!: arg_cong[of _ _ h])
also have "… = [: h (int df - dfj) ^ (n - 1 - dfj) :]"
unfolding prod_constant by (simp add: poly_const_pow)
also have "n - 1 - dfj = dg - J - 1" unfolding n dfj_def by simp
also have "int df - dfj = J" using * dfg unfolding dfj_def by auto
also have "h J = coeff H J" unfolding h coeff_int_def by simp
finally show "subresultant J F G = smult (?m1 * ?G * ?H) H" by (simp add: dfj_def ac_simps)
} note eq_19 = this
{
assume J: "J < dh"
define dhj where "dhj = dh - J"
have n_add: "n = (df - dh) + (dhj + dgj)" unfolding dhj_def dgj_def n using J dfg dgh by auto
let ?split = "split_block ?G_H (df - dh) (df - dh)"
have dim: "dim_row ?G_H = (df - dh) + (dhj + dgj)"
"dim_col ?G_H = (df - dh) + (dhj + dgj)"
unfolding n_add by auto
obtain UL UR LL LR where spl: "?split = (UL,UR,LL,LR)" by (cases ?split, auto)
note spl' = spl[unfolded split_block_def Let_def, simplified]
let ?LR = "subresultant_mat J G H"
have "LR = mat (dgj + dhj) (dgj + dhj)
(λ (i,j). ?GH (i + (df - dh)) (j + (df - dh)))"
using spl' by (auto simp: n_add)
also have "… = ?LR"
unfolding subresultant_mat_def Let_def dhj_def dgj_def d[symmetric]
proof (rule eq_matI, unfold dim_row_mat dim_col_mat index_mat split dfj_def, goal_cases)
case (1 i j)
hence id1: "(j + (df - dh) < df - J) = (j < dh - J)" using dgh dfg J by auto
have id2: "(i + (df - dh) = n - 1) = (i = dg - J + (dh - J) - 1)"
unfolding n_add dhj_def dgj_def using dgh dfg J by auto
have id3: "(df - J - 1 - (j + (df - dh))) = (dh - J - 1 - j)"
and id4: "(int dg - int (i + (df - dh)) + int (j + (df - dh))) = (int dg - int i + int j)"
and id5: "(dg - J - 1 - (j + (df - dh) - (df - J))) = (dg - J - 1 - (j - (dh - J)))"
and id6: "(int df - int (i + (df - dh)) + int (j + (df - dh) - (df - J))) = (int dh - int i + int (j - (dh - J)))"
using dgh dfg J by auto
show ?case unfolding g[symmetric] h[symmetric] id3 id4 id5 id6
by (rule if_cong[OF id1 if_cong[OF id2 refl refl] if_cong[OF id2 refl refl]])
qed auto
finally have "LR = ?LR" .
note spl = spl[unfolded this]
let ?UR = "0⇩m (df - dh) (dgj + dhj)"
have "UR = mat (df - dh) (dgj + dhj)
(λ (i,j). ?GH i (j + (df - dh)))"
using spl' by (auto simp: n_add)
also have "… = ?UR"
proof (rule eq_matI, unfold dim_row_mat dim_col_mat index_mat split dfj_def index_zero_mat, goal_cases)
case (1 i j)
hence in1: "i ≠ n - 1" using J unfolding dgj_def dhj_def n_add by auto
{
assume "j + (df - dh) < df - J"
hence "dg < int dg - int i + int (j + (df - dh))" using 1 J unfolding dgj_def dhj_def by auto
hence "g … = 0" unfolding dg g by (intro coeff_int_eq_0, auto)
} note g = this
{
assume "¬ (j + (df - dh) < df - J)"
hence "dh < int df - int i + int (j + (df - dh) - (df - J))" using 1 J unfolding dgj_def dhj_def by auto
hence "h … = 0" unfolding dh h by (intro coeff_int_eq_0, auto)
} note h = this
show ?case using in1 g h by auto
qed auto
finally have "UR = ?UR" .
note spl = spl[unfolded this]
let ?G = "λ (i,j). if i = j then [:lead_coeff G:] else if i < j then 0 else ?GH i j"
let ?UL = "mat (df - dh) (df - dh) ?G"
have "UL = mat (df - dh) (df - dh) (λ (i,j). ?GH i j)"
using spl' by (auto simp: n_add)
also have "… = ?UL"
proof (rule eq_matI, unfold dim_row_mat dim_col_mat index_mat split, goal_cases)
case (1 i j)
{
assume "i = j"
hence "int dg - int i + int j = dg" using 1 by auto
hence "g (int dg - int i + int j) = lead_coeff G"
unfolding g dg coeff_int_def by simp
} note eq = this
{
assume "i < j"
hence "dg < int dg - int i + int j" using 1 by auto
hence "g (int dg - int i + int j) = 0"
unfolding g dg by (intro coeff_int_eq_0, auto)
} note lt = this
from 1 have *: "j < dfj" "i ≠ n - 1" using J unfolding n_add dhj_def dgj_def dfj_def by auto
hence "?GH i j = [:g (int dg - int i + int j):]" by simp
also have "… = (if i = j then [: lead_coeff G :] else if i < j then 0 else ?GH i j)"
using eq lt * by auto
finally show ?case by simp
qed auto
finally have "UL = ?UL" .
note spl = spl[unfolded this]
from split_block[OF spl dim]
have GH: "?G_H = four_block_mat ?UL ?UR LL ?LR"
and C: "?UL ∈ carrier_mat (df - dh) (df - dh)"
"?UR ∈ carrier_mat (df - dh) (dhj + dgj)"
"LL ∈ carrier_mat (dhj + dgj) (df - dh)"
"?LR ∈ carrier_mat (dhj + dgj) (dhj + dgj)" by auto
from arg_cong[OF GH, of det]
have "det ?G_H = det (four_block_mat ?UL ?UR LL ?LR)" unfolding GH[symmetric] ..
also have "… = det ?UL * det ?LR"
by (rule det_four_block_mat_upper_right_zero[OF _ refl], insert C, auto simp: ac_simps)
also have "det ?LR = subresultant J G H" unfolding subresultant_def by simp
also have "det ?UL = prod_list (diag_mat ?UL)"
by (rule det_lower_triangular[of "df - dh"], auto)
also have "… = (∏i = 0..< (df - dh). [: lead_coeff G :])" unfolding prod_list_diag_prod by simp
also have "… = [: lead_coeff G ^ (df - dh) :]" by (simp add: poly_const_pow)
finally have det: "det ?G_H = [:lead_coeff G ^ (df - dh):] * subresultant J G H" by auto
show "subresultant J F G = smult (?m1 * lead_coeff G ^ (df - dh)) (subresultant J G H)"
unfolding eq_18 det by simp
}
{
assume J: "dh < J" "J < dg - 1"
hence "dh ≤ J" "J < dg" by auto
from eq_19[OF this]
have "subresultant J F G = smult ((- 1) ^ ((df - J) * (dg - J)) * lead_coeff G ^ (df - J) * coeff H J ^ (dg - J - 1)) H"
by simp
also have "coeff H J = 0" by (rule coeff_eq_0, insert J, auto simp: dh)
also have "… ^ (dg - J - 1) = 0" using J by auto
finally show "subresultant J F G = 0" by simp
}
{
assume J: "J = dh" and "dg > dh ∨ H ≠ 0"
with choice have dgh: "dg > dh" by auto
show "subresultant dh F G = smult (
(-1)^((df - dh) * (dg - dh)) * lead_coeff G ^ (df - dh) * lead_coeff H ^ (dg - dh - 1)) H"
unfolding eq_19[unfolded J, OF le_refl dgh] unfolding dh by simp
}
{
assume J: "J = dg - 1" and "dg > dh ∨ H ≠ 0"
with choice have dgh: "dg > dh" by auto
have *: "dh ≤ dg - 1" "dg - 1 < dg" using dgh by auto
have **: "df - (dg - 1) = df - dg + 1" "dg - (dg - 1) - 1 = 0" "dg - (dg - 1) = 1"
using dfg dgh by linarith+
show "subresultant (dg - 1) F G = smult (
(-1)^(df - dg + 1) * lead_coeff G ^ (df - dg + 1)) H"
unfolding eq_19[unfolded J, OF *] unfolding ** by simp
}
qed
lemmas BT_lemma_1_13 = BT_lemma_1_13'[OF _ _ _ refl]
lemmas BT_lemma_1_15 = BT_lemma_1_15'[OF _ _ _ refl]
lemma subresultant_product: fixes F :: "'a :: idom poly"
assumes "F = B * G"
and FG: "degree F ≥ degree G"
shows "subresultant J F G = (if J < degree G then 0 else
if J < degree F then smult (lead_coeff G ^ (degree F - J - 1)) G else 1)"
proof (cases "J < degree G")
case J: True
from assms have eq: "F + (-B) * G = 0" by auto
from J have lt: "degree 0 < degree G ∨ b" for b by auto
from BT_lemma_1_13[OF eq FG lt lt]
have "subresultant 0 F G = 0" using J by auto
with BT_lemma_1_14[OF eq FG lt, of J] have 00: "J = 0 ∨ J < degree G - 1 ⟹ subresultant J F G = 0" by auto
from BT_lemma_1_15[OF eq FG lt lt] J have 01: "subresultant (degree G - 1) F G = 0" by simp
from J have "(J = 0 ∨ J < degree G - 1) ∨ J = degree G - 1" by linarith
with 00 01 have "subresultant J F G = 0" by auto
thus ?thesis using J by simp
next
case J: False
hence dg: "degree G - J = 0" by simp
let ?n = "degree F - J"
have *: "(j :: nat) < 0 ⟷ False" "j - 0 = j" for j by auto
let ?M = "mat ?n ?n
(λ(i, j).
if i = ?n - 1 then monom 1 (?n - 1 - j) * G
else [:coeff_int G (int (degree G) - int i + int j):])"
have "subresultant J F G = det ?M"
unfolding subresultant_def subresultant_mat_def Let_def dg * by auto
also have "det ?M = prod_list (diag_mat ?M)"
by (rule det_lower_triangular[of ?n], auto intro: coeff_int_eq_0)
also have "… = (∏i = 0..< ?n. ?M $$ (i,i))" unfolding prod_list_diag_prod by simp
also have "… = (∏i = 0..< ?n. if i = ?n - 1 then G else [: lead_coeff G :])"
by (rule prod.cong[OF refl], auto simp: coeff_int_def)
also have "… = (if J < degree F then smult (lead_coeff G ^ (?n - 1)) G else 1)"
proof (cases "J < degree F")
case True
hence id: "{ 0 ..< ?n} = { 0 ..< ?n - 1} ∪ {?n - 1}" by auto
have "(∏i = 0..< ?n. if i = ?n - 1 then G else [: lead_coeff G :])
= (∏i = 0 ..< ?n - 1. if i = ?n - 1 then G else [: lead_coeff G :]) * G" (is "_ = ?P * G")
unfolding id
by (subst prod.union_disjoint, auto)
also have "?P = (∏i = 0 ..< ?n - 1. [: lead_coeff G :])"
by (rule prod.cong, auto)
also have "… = [: lead_coeff G ^ (?n - 1) :]"
by (simp add: poly_const_pow)
finally show ?thesis by auto
qed auto
finally have "subresultant J F G =
(if J < degree F then smult (lead_coeff G ^ (degree F - J - 1)) G else 1)" .
thus ?thesis using J by simp
qed
lemma resultant_pseudo_mod_0: assumes "pseudo_mod f g = (0 :: 'a :: idom_divide poly)"
and dfg: "degree f ≥ degree g"
and f: "f ≠ 0" and g: "g ≠ 0"
shows "resultant f g = (if degree g = 0 then lead_coeff g^degree f else 0)"
proof -
let ?df = "degree f" let ?dg = "degree g"
obtain d r where pd: "pseudo_divmod f g = (d,r)" by force
from pd have r: "r = pseudo_mod f g" unfolding pseudo_mod_def by simp
with assms pd have pd: "pseudo_divmod f g = (d,0)" by auto
from pseudo_divmod[OF g pd] g
obtain a q where prod: "smult a f = g * q" and a: "a ≠ 0" "a = lead_coeff g ^ (Suc ?df - ?dg)"
by auto
from a dfg have dfg: "degree g ≤ degree (smult a f)" by auto
have g0: "degree g = 0 ⟹ coeff g 0 = 0 ⟹ g = 0"
using leading_coeff_0_iff by fastforce
from prod have "smult a f = q * g" by simp
from arg_cong[OF subresultant_product[OF this dfg, of 0, unfolded subresultant_resultant
resultant_smult_left[OF a(1)]], of "λ x. coeff x 0"]
show ?thesis using a g0 by (cases "degree f", auto)
qed
locale primitive_remainder_sequence =
fixes F :: "nat ⇒ 'a :: idom_divide poly"
and n :: "nat ⇒ nat"
and δ :: "nat ⇒ nat"
and f :: "nat ⇒ 'a"
and k :: nat
and β :: "nat ⇒ 'a"
assumes f: "⋀ i. f i = lead_coeff (F i)"
and n: "⋀ i. n i = degree (F i)"
and δ: "⋀ i. δ i = n i - n (Suc i)"
and n12: "n 1 ≥ n 2"
and F12: "F 1 ≠ 0" "F 2 ≠ 0"
and F0: "⋀ i. i ≠ 0 ⟹ F i = 0 ⟷ i > k"
and β0: "⋀ i. β i ≠ 0"
and pmod: "⋀ i. i ≥ 3 ⟹ i ≤ Suc k ⟹ smult (β i) (F i) = pseudo_mod (F (i - 2)) (F (i - 1))"
begin
lemma f10: "f 1 ≠ 0" and f20: "f 2 ≠ 0" unfolding f using F12 by auto
lemma f0: "i ≠ 0 ⟹ f i = 0 ⟷ i > k"
using F0[of i] unfolding f by auto
lemma n_gt: assumes "2 ≤ i" "i < k"
shows "n i > n (Suc i)"
proof -
from assms have "3 ≤ Suc i" "Suc i ≤ Suc k" by auto
note pmod = pmod[OF this]
from assms F0 have "F (Suc i - 1) ≠ 0" "F (Suc i) ≠ 0" by auto
from pseudo_mod(2)[OF this(1), of "F (Suc i - 2)", folded pmod] this(2)
show ?thesis unfolding n using β0 by auto
qed
lemma n_ge: assumes "1 ≤ i" "i < k"
shows "n i ≥ n (Suc i)"
using n12 n_gt[OF _ assms(2)] assms(1) by (cases "i = 1", auto simp: numeral_2_eq_2)
lemma n_ge_trans: assumes "1 ≤ i" "i ≤ j" "j ≤ k"
shows "n i ≥ n j"
proof -
from assms(2) have "j = i + (j - i)" by simp
then obtain jj where j: "j = i + jj" by blast
from assms(3)[unfolded j] show ?thesis unfolding j
proof (induct jj)
case (Suc j)
from Suc(2) have "i + j ≤ k" by simp
from Suc(1)[OF this] have IH: "n (i + j) ≤ n i" .
have "n (Suc (i + j)) ≤ n (i + j)"
by (rule n_ge, insert assms(1) Suc(2), auto)
with IH show ?case by auto
qed auto
qed
lemma delta_gt: assumes "2 ≤ i" "i < k"
shows "δ i > 0" using n_gt[OF assms] unfolding δ by auto
lemma k2:"2 ≤ k"
by (metis le_cases linorder_not_le F0 F12(2) zero_order(2))
lemma k0: "k ≠ 0" using k2 by auto
lemma ni2:"3 ≤ i ⟹ i ≤ k ⟹ n i ≠ n 2"
by (metis Suc_numeral δ delta_gt k2 le_imp_less_Suc le_less n_ge_trans not_le one_le_numeral
semiring_norm(5) zero_less_diff)
end
locale subresultant_prs_locale = primitive_remainder_sequence F n δ f k β for
F :: "nat ⇒ 'a :: idom_divide fract poly"
and n :: "nat ⇒ nat"
and δ :: "nat ⇒ nat"
and f :: "nat ⇒ 'a fract"
and k :: nat
and β :: "nat ⇒ 'a fract" +
fixes G1 G2 :: "'a poly"
assumes F1: "F 1 = map_poly to_fract G1"
and F2: "F 2 = map_poly to_fract G2"
begin
definition "α i = (f (i - 1))^(Suc (δ (i - 2)))"
lemma α0: "i > 1 ⟹ α i = 0 ⟷ (i - 1) > k"
unfolding α_def using f0[of "i - 1"] by auto
lemma α_char:
assumes "3 ≤ i" "i < k + 2"
shows "α i = (f (i - 1)) ^ (Suc (length (coeffs (F (i - 2)))) - length (coeffs (F (i - 1))))"
proof (cases "i = 3")
case True
have triv:"Suc (Suc 0) = 2" by auto
have l:"length (coeffs (F 2)) ≠ 0" "length (coeffs (F 1)) ≠ 0" using F12 by auto
hence "length (coeffs (F 2)) ≤ length (coeffs (F (Suc 0)))" using n12
unfolding n degree_eq_length_coeffs One_nat_def by linarith
hence "Suc (length (coeffs (F 1)) - 1 - (length (coeffs (F 2)) - 1)) =
(Suc (length (coeffs (F 1))) - length (coeffs (F (3 - 1))))" using l by simp
thus ?thesis unfolding True α_def n δ degree_eq_length_coeffs by (simp add:triv)
next
case False
hence assms:"2 ≤ i - 2" "i - 2 < k" using assms by auto
have i:"i - 2 ≠ 0" "i - 1 ≠ 0" using assms by auto
hence [simp]: "Suc (i - 2) = i - 1" by auto
from assms(2) F0[OF i(2)] have "F (i - 1) ≠ 0" by auto
then have "length (coeffs (F (i - 1))) > 0" by (cases "F (i - 1)") auto
with delta_gt[unfolded δ n degree_eq_length_coeffs,OF assms]
have * : "Suc (δ (i - 2)) = Suc (length (coeffs (F (i - 2)))) - (length (coeffs (F (Suc (i - 2)))))"
by (auto simp:δ n degree_eq_length_coeffs)
show ?thesis unfolding α_def * by simp
qed
definition Q :: "nat ⇒ 'a fract poly" where
"Q i ≡ smult (α i) (fst (pdivmod (F (i - 2)) (F (i - 1))))"
lemma beta_F_as_sum:
assumes "3 ≤ i" "i ≤ Suc k"
shows "smult (β i) (F i) = smult (α i) (F (i - 2)) + - Q i * F (i - 1)" (is ?t1)
proof -
have ik:"i < k + 2" using assms by auto
have f0:"F (i - 1) = 0 ⟷ False" "F (i - Suc 0) = 0 ⟷ False"
using F0[of "i - 1"] assms by auto
hence f0_b:"(inverse (coeff (F (i - 1)) (degree (F (i - 1))))) ≠ 0" "F (i - 1) ≠ 0" by auto
have i:"i - 2 ≠ 0" "Suc (i - 2) = i - 1" "(k < i - 2) ⟷ False" using assms by auto
have "F (i - 2) ≠ 0" using F0[of "i - 2"] assms by auto
let ?c = "(inverse (f (i - 1)) ^ (Suc (length (coeffs (F (i - 2)))) - length (coeffs (F (i - 1)))))"
have inv:"inverse (α i) = ?c" unfolding α_char[OF assms(1) ik] power_inverse by auto
have alpha0:"α i ≠ 0" unfolding α_def f using f0 by auto
have α_inv[simp]:"α i * inverse (α i) = 1"
using field_class.field_inverse[OF alpha0] mult.commute by metis
with field_class.field_inverse[OF alpha0,unfolded inv]
have c_times_Q:"smult ?c (Q i) = fst (pdivmod (F (i - 2)) (F (i - 1)))"
unfolding Q_def by auto
have "pdivmod (F (i - 2)) (F (i - 1)) = (smult ?c (Q i), smult ?c (smult (β i) (F i)))"
unfolding c_times_Q
unfolding pdivmod_via_pseudo_divmod pmod[OF assms] f n c_times_Q
pseudo_mod_smult_right[OF f0_b, of "F (i - 2)",symmetric] f0 if_False Let_def
unfolding pseudo_mod_def by (auto split:prod.split)
from this [symmetric]
have pr: ‹F (i - 2) = smult ?c (Q i) * F (i - 1) + smult ?c ( smult (β i) (F i))›
by (simp only: prod_eq_iff fst_conv snd_conv div_mult_mod_eq)
then have "F (i - 2) = smult (inverse (α i)) (Q i) * F (i - 1)
+ smult (inverse (α i)) ( smult (β i) (F i))" (is "?l = ?r" is "_ = ?t + _")
unfolding inv.
hence eq:"smult (α i) (?l - ?t) = smult (α i) (?r - ?t)" by auto
have " smult (α i) (F (i - 2)) - Q i * (F (i - 1)) = smult (α i) (?l - ?t)"
unfolding smult_diff_right by auto
also have "… = smult (α i) (?r - ?t)" unfolding eq..
also have "… = smult (β i) (F i)" by (auto simp:mult.assoc[symmetric])
finally show ?t1 by auto
qed
lemma assumes "3 ≤ i" "i ≤ k" shows
BT_lemma_2_21: "j < n i ⟹ smult (α i ^ (n (i - 1) - j)) (subresultant j (F (i - 2)) (F (i - 1)))
= smult ((- 1) ^ ((n (i - 2) - j) * (n (i - 1) - j)) * (f (i - 1)) ^ (δ (i - 2) + δ (i - 1)) * (β i) ^ (n (i - 1) - j)) (subresultant j (F (i - 1)) (F i))"
(is "_ ⟹ ?eq_21") and
BT_lemma_2_22: "smult (α i ^ (δ (i - 1))) (subresultant (n i) (F (i - 2)) (F (i - 1)))
= smult ((- 1) ^ ((δ (i - 2) + δ (i - 1)) * δ (i - 1)) * f (i - 1) ^ (δ (i - 2) + δ (i - 1)) * f i ^ (δ (i - 1) - 1) * (β i) ^ δ (i - 1)) (F i)"
(is "?eq_22") and
BT_lemma_2_23: "n i < j ⟹ j < n (i - 1) - 1 ⟹ subresultant j (F (i - 2)) (F (i - 1)) = 0"
(is "_ ⟹ _ ⟹ ?eq_23") and
BT_lemma_2_24: "smult (α i) (subresultant (n (i - 1) - 1) (F (i - 2)) (F (i - 1)))
= smult ((- 1) ^ (δ (i - 2) + 1) * f (i - 1) ^ (δ (i - 2) + 1) * β i) (F i)" (is "?eq_24")
proof -
from assms have ik:"i ≤ Suc k" by auto
note beta_F_as_sum = beta_F_as_sum[OF assms(1) ik, symmetric]
have s[simp]:"Suc (i - 2) = i - 1" "Suc (i - 1) = i" using assms by auto
have α0:"α i ≠ 0" using assms f0[of "i - 1"] unfolding α_def f by auto
hence α0pow:"⋀ x. α i ^ x ≠ 0" by auto
have df:"degree (F (i - 1)) ≤ degree (smult (α i) (F (i - 2)))"
"degree (smult (β i) (F i)) < degree (F (i - 1)) ∨ b" for b
using n_ge[of "i-2"] n_gt[of "i-1"] assms α0 β0 unfolding n by auto
have degree_smult_eq:"⋀ c f. (c::_::{idom_divide}) ≠ 0 ⟹ degree (smult c f) = degree f" by auto
have n_lt:"n i < n (i - 1)" using n_gt[of "i-1"] assms unfolding n by auto
from semiring_normalization_rules(30) mult.commute
have *:"⋀ x y q. (x * (y::'a fract)) ^ q = y ^ q * x ^ q" by metis
have "n (i - 1) - n i > 0" using n_lt by auto
hence **:"β i ^ (n (i - 1) - n i - 1) * β i = β i ^ (n (i - 1) - n i)"
by (subst power_minus_mult) auto
have "max (n (i - 2)) (n (i - 1)) = n (i - 2)" using n_ge[of "i-2"] assms
unfolding max_def by auto
with diff_add_assoc[OF n_ge[of "i-1"],symmetric] assms
have ns : "n (i - 2) - n (i - 1) + (n (i - 1) - n i) = n (i - 2) - n i"
by (auto simp:nat_minus_add_max)
{ assume "j < n i"
hence j:"j < degree (smult (β i) (F i))" using β0 unfolding n by auto
from BT_lemma_1_12[OF beta_F_as_sum df j]
show ?eq_21
unfolding subresultant_smult_right[OF β0] subresultant_smult_left[OF α0]
degree_smult_eq[OF α0] degree_smult_eq[OF β0] n[symmetric] f[symmetric] δ s ns
using f n
by auto}
{ from BT_lemma_1_13[OF beta_F_as_sum df df(2)]
show ?eq_22
unfolding subresultant_smult_left[OF α0] lead_coeff_smult smult_smult
degree_smult_eq[OF α0] degree_smult_eq[OF β0] n[symmetric] f[symmetric] δ s ns
by (metis (no_types, lifting) * ** coeff_smult f mult.assoc n)}
{ assume "n i < j" "j < n (i - 1) - 1"
hence j:"degree (smult (β i) (F i)) < j" "j < degree (F (i - 1)) - 1"
using β0 unfolding n by auto
from BT_lemma_1_14[OF beta_F_as_sum df j]
show ?eq_23 unfolding subresultant_smult_left[OF α0] smult_eq_0_iff using α0pow by auto}
{ have ***: "n (i - 1) - (n (i - 1) - 1) = 1" using n_lt by auto
from BT_lemma_1_15[OF beta_F_as_sum df df(2)]
show ?eq_24
unfolding subresultant_smult_left[OF α0] *** degree_smult_eq[OF α0] n[symmetric] f δ
by (auto simp:mult.commute)}
qed
lemma BT_eq_30: "3 ≤ i ⟹ i ≤ k + 1 ⟹ j < n (i - 1) ⟹
smult (∏l←[3..<i]. α l ^ (n (l - 1) - j)) (subresultant j (F 1) (F 2))
= smult (∏l←[3..<i]. β l ^ (n (l - 1) - j) * f (l - 1) ^ (δ (l - 2) + δ (l - 1))
* (- 1) ^ ((n (l - 2) - j) * (n (l - 1) - j))) (subresultant j (F (i - 2)) (F (i - 1)))"
proof (induct "i - 3" arbitrary:i)
case (Suc x)
from Suc.hyps(2) Suc.prems(1-2)
have prems:"x = (i - 1) - 3" "3 ≤ i - 1" "i - 1 ≤ k + 1" "2 ≤ i - 1 - 1" "i - 1 - 1 < k"
"i - 1 ≤ k" by auto
from prems(2) have inset:"i - 1 ∈ set [3..<i]" by auto
have r1:"remove1 (i - 1) [3..<i] = [3..<i-1]" by (induct i,auto simp:remove1_append)
from Suc.prems(1) have "Suc (i - 1 - 1) = i - 1" by auto
from n_gt[OF prems(4,5),unfolded this] Suc.prems(3) have j:"j < n (i - 1 - 1)" by auto
have *:"⋀ c d e x. smult c d = e ⟹ smult (x * c) d = smult x e" by auto
have **:"⋀ c d e x. smult c d = e ⟹ smult c (smult x d) = smult x e" by (auto simp:mult.commute)
show ?case unfolding prod_list_map_remove1[OF inset(1),unfolded r1]
*[OF Suc.hyps(1)[OF prems(1-3) j]]
**[OF BT_lemma_2_21[OF prems(2,6) Suc.prems(3)]]
by (auto simp: numeral_2_eq_2 ac_simps)
qed auto
lemma nonzero_alphaprod: assumes "i ≤ k + 1" shows "(∏l←[3..<i]. α l ^ (p l)) ≠ 0"
unfolding prod_list_zero_iff using assms by (auto simp: α0)
lemma BT_eq_30': assumes i: "3 ≤ i" "i ≤ k + 1" "j < n (i - 1)"
shows "subresultant j (F 1) (F 2)
= smult ((- 1) ^ (∑l←[3..<i]. (n (l - 2) - j) * (n (l - 1) - j))
* (∏l←[3..<i]. (β l / α l) ^ (n (l - 1) - j)) * (∏l←[3..<i]. f (l - 1) ^ (δ (l - 2) + δ (l - 1)))) (subresultant j (F (i - 2)) (F (i - 1)))"
(is "_ = smult (?mm * ?b * ?f) _")
proof -
let ?a = "∏l←[3..<i]. α l ^ (n (l - 1) - j)"
let ?d = "∏l←[3..<i]. β l ^ (n (l - 1) - j) * f (l - 1) ^ (δ (l - 2) + δ (l - 1)) *
(- 1) ^ ((n (l - 2) - j) * (n (l - 1) - j))"
let ?m = "∏l←[3..<i]. (- 1) ^ ((n (l - 2) - j) * (n (l - 1) - j))"
have a0: "?a ≠ 0" by (rule nonzero_alphaprod, rule i)
with arg_cong[OF BT_eq_30[OF i], of "smult (inverse ?a)", unfolded smult_smult]
have "subresultant j (F 1) (F 2) = smult (inverse ?a * ?d)
(subresultant j (F (i - 2)) (F (i - 1)))"
by simp
also have "inverse ?a * ?d = ?b * ?f * ?m" unfolding prod_list_multf inverse_prod_list map_map o_def
power_inverse[symmetric] power_mult_distrib divide_inverse_commute
by simp
also have "?m = ?mm"
unfolding prod_list_minus_1_exp by simp
finally show ?thesis by (simp add: ac_simps)
qed
text ‹For defining the subresultant PRS, we mainly follow Brown's ``The Subresultant PRS Algorithm'' (B).›
definition "R j = (if j = n 2 then sdiv_poly (smult ((lead_coeff G2)^(δ 1)) G2) (lead_coeff G2) else subresultant j G1 G2)"
abbreviation "ff i ≡ to_fract (i :: 'a)"
abbreviation "ffp ≡ map_poly ff"
sublocale map_poly_hom: map_poly_inj_idom_hom to_fract..
definition "σ i = (∑l←[3..<Suc i]. (n (l - 2) + n (i - 1) + 1) * (n (l - 1) + n (i - 1) + 1))"
definition "τ i = (∑l←[3..<Suc i]. (n (l - 2) + n i) * (n (l - 1) + n i))"
definition "γ i = (-1)^(σ i) * pow_int ( f (i - 1)) (1 - int (δ (i - 1))) * (∏l←[3..<Suc i].
(β l / α l)^(n (l - 1) - n (i - 1) + 1) * (f (l - 1))^(δ (l - 2) + δ (l - 1)))"
definition "Θ i = (-1)^(τ i) * pow_int (f i) (int (δ (i - 1)) - 1) * (∏l←[3..<Suc i].
(β l / α l)^(n (l - 1) - n i) * (f (l - 1))^(δ (l - 2) + δ (l - 1)))"
lemma fundamental_theorem_eq_4: assumes i: "3 ≤ i" "i ≤ k"
shows "ffp (R (n (i - 1) - 1)) = smult (γ i) (F i)"
proof -
have "n (i - 1) ≤ n 2" by (rule n_ge_trans, insert i, auto)
with n_gt[of "i - 1"] i have "n (i - 1) - 1 < n 2"
and lt: "n (i - 1) - 1 < n (i - 1)" by linarith+
hence "R (n (i - 1) - 1) = subresultant (n (i - 1) - 1) G1 G2"
unfolding R_def by auto
from arg_cong[OF this, of ffp, unfolded to_fract_hom.subresultant_hom, folded F1 F2]
have id1: "ffp (R (n (i - 1) - 1)) = subresultant (n (i - 1) - 1) (F 1) (F 2)" .
note eq_24 = BT_lemma_2_24[OF i]
let ?o = "(- 1) :: 'a fract"
let ?m1 = "(δ (i - 2) + 1)"
let ?d1 = "f (i - 1) ^ (δ (i - 2) + 1) * β i"
let ?c1 = "?o ^ ?m1 * ?d1"
let ?c0 = "α i"
have "?c0 ≠ 0" using α0[of i] i by auto
with arg_cong[OF eq_24, of "smult (inverse ?c0)"]
have id2: "subresultant (n (i - 1) - 1) (F (i - 2)) (F (i - 1)) =
smult (inverse ?c0 * ?c1) (F i)"
by (auto intro: poly_eqI)
from i have "3 ≤ i" "i ≤ k + 1" by auto
note id3 = BT_eq_30'[OF this lt]
let ?f = "λ l. f (l - 1) ^ (δ (l - 2) + δ (l - 1))"
let ?b = "λ l. (β l / α l) ^ (n (l - 1) - (n (i - 1) - 1))"
let ?b' = "λ l. (β l / α l) ^ (n (l - 1) - n (i - 1) + 1)"
let ?m = "λ l. (n (l - 2) - (n (i - 1) - 1)) * (n (l - 1) - (n (i - 1) - 1))"
let ?m' = "λ l. (n (l - 2) + n (i - 1) + 1) * (n (l - 1) + n (i - 1) + 1)"
let ?m2 = "(∑l←[3..<i]. ?m l)"
let ?b2 = "(∏l←[3..<i]. ?b l)"
let ?f2 = "(∏l←[3..<i]. ?f l)"
let ?f1 = "pow_int ( f (i - 1)) (1 - int (δ (i - 1)))"
have id4: "γ i = ?o ^ (?m1 + ?m2) * (inverse ?c0 * ?d1 * ?b2 * ?f2)"
proof -
have id: "γ i = (-1)^(σ i) * (?f1 * (∏l←[3..<Suc i]. ?b' l) * (∏l←[3..<Suc i]. ?f l))"
unfolding γ_def prod_list_multf by simp
have cong: "even m1 = even m2 ⟹ c1 = c2 ⟹ ?o^m1 * c1 = ?o^m2 * c2" for m1 m2 c1 c2
unfolding minus_1_power_even by auto
show ?thesis unfolding id
proof (rule cong)
from n_gt[of "i - 1"] i have n1: "n (i - 1) ≠ 0" by linarith
{
fix l
assume "2 ≤ l" "l ≤ i"
hence l: "l ≥ 2" "l - 1 ≤ i - 1" "l ≤ k" using i by auto
from n_ge_trans[OF _ l(2)] l i have n2: "n (i - 1) ≤ n (l - 1)" by auto
from n1 n2 have id: "n (l - 1) - (n (i - 1) - 1) = n (l - 1) - n (i - 1) + 1" by auto
have "even (n (l - 1) - (n (i - 1) - 1)) = even (n (l - 1) + n (i - 1) + 1)"
unfolding id using n2 by auto
note id n2 this
} note diff = this
have f0: "f (i - 1) ≠ 0" using f0[of "i - 1"] i by auto
have "(∏l←[3..<Suc i]. ?b' l) = (∏l←[3..<Suc i]. ?b l)"
by (rule arg_cong, rule map_cong, use diff(1) in auto)
also have "… = ?b2 * ?b i" using i by auto
finally have "?f1 * (∏l←[3..<Suc i]. ?b' l) * (∏l←[3..<Suc i]. ?f l) =
(?b2 * ?f2) * (?f1 * ?b i * ?f i)" using i by simp
also have "?f1 * ?b i * ?f i = (?f1 * ?f i) * β i * inverse ?c0" using n1 by (simp add: divide_inverse)
also have "?f1 * ?f i = f (i - 1) ^ (δ (i - 2) + 1)"
unfolding exp_pow_int pow_int_add[OF f0, symmetric] by simp
finally
show "?f1 * (∏l←[3..<Suc i]. ?b' l) * (∏l←[3..<Suc i]. ?f l)
= inverse ?c0 * ?d1 * ?b2 * ?f2" by simp
have "even (σ i) = even ((∑l←[3..<i]. ?m' l) + ?m' i)" unfolding σ_def using i by simp
also have "… = (even (∑l←[3..<i]. ?m' l) = even (?m' i))" by simp
also have "even (∑l←[3..<i]. ?m' l) = even ?m2"
proof (rule even_sum_list, goal_cases)
case (1 l)
hence l: "l ≥ 2" "l ≤ i" and l1: "l - 1 ≥ 2" "l - 1 ≤ i" by auto
have l2: "l - 2 = l - 1 - 1" by simp
show ?case using diff(3) [OF l] diff(3) [OF l1] l2
by auto
qed
also have "even (?m' i) = even ?m1"
proof -
from i have id: "Suc (i - 1 - 1) = i - 1" "i - 2 = i - 1 - 1" by auto
have "even ?m1 = even (n (i - 2) + n (i - 1) + 1)" unfolding δ id
using diff[of "i - 1"] i by auto
also have "… = even (?m' i)" by auto
finally show ?thesis by simp
qed
also have "(even ?m2 = even ?m1) = even (?m2 + ?m1)" unfolding even_add by simp
also have "?m2 + ?m1 = ?m1 + ?m2" by simp
finally show "even (σ i) = even (?m1 + ?m2)" .
qed
qed
show ?thesis unfolding id1 id3 id2 smult_smult id4 by (simp add: ac_simps power_add)
qed
lemma fundamental_theorem_eq_5: assumes i: "3 ≤ i" "i ≤ k" "n i < j" "j < n (i - 1) - 1"
shows "R j = 0"
proof -
from BT_lemma_2_23[OF i] have id1: "subresultant j (F (i - 2)) (F (i - 1)) = 0" .
have "n (i - 1) ≤ n 2" by (rule n_ge_trans, insert i, auto)
with n_gt[of "i - 1"] i have "n (i - 1) - 1 < n 2"
and lt: "j < n (i - 1)" by linarith+
with i have "R j = subresultant j G1 G2" unfolding R_def by auto
from arg_cong[OF this, of ffp, unfolded to_fract_hom.subresultant_hom, folded F1 F2]
have id2: "ffp (R j) = subresultant j (F 1) (F 2)" .
from i have "3 ≤ i" "i ≤ k + 1" by auto
note eq_30 = BT_eq_30[OF this lt]
let ?c3 = "∏l←[3..<i]. α l ^ (n (l - 1) - j)"
let ?c2 = "∏l←[3..<i]. β l ^ (n (l - 1) - j) * f (l - 1) ^ (δ (l - 2) + δ (l - 1)) *
(- 1) ^ ((n (l - 2) - j) * (n (l - 1) - j))"
have "?c3 ≠ 0" by (rule nonzero_alphaprod, insert i, auto)
with arg_cong[OF eq_30, of "smult (inverse ?c3)"]
have id3: "subresultant j (F 1) (F 2) = smult (inverse ?c3 * ?c2)
(subresultant j (F (i - 2)) (F (i - 1)))"
by (auto intro: poly_eqI)
have "ffp (R j) = 0" unfolding id1 id2 id3 by simp
thus ?thesis by simp
qed
lemma fundamental_theorem_eq_6: assumes "3 ≤ i" "i ≤ k" shows "ffp (R (n i)) = smult (Θ i) (F i)"
(is "?lhs=?rhs")
proof -
from assms have i1:"1 ≤ i" by auto
from assms have nlt:"i ≤ k + 1" "n i < n (i - 1)" using n_gt[of "i - 1"] by auto
from assms have αnz:"α i ^ δ (i - 1) ≠ 0" using α0 by auto
have *:"⋀ a f b. a ≠ 0 ⟹ smult a f = b ⟹ f = smult (inverse (a::'a fract)) b" by auto
have **:"⋀ f g xs c. c * prod_list (map f xs) * prod_list (map g xs)
= c * (∏x←xs. f x * (g:: _ ⇒ (_ :: comm_monoid_mult)) x)"
by (auto simp:ac_simps prod_list_multf)
have ***:"⋀ c. β i ^ δ (i - Suc 0) * (inverse (α i ^ δ (i - Suc 0)) * c) = (β i / α i) ^ δ (i - 1) * c"
by (auto simp:inverse_eq_divide power_divide)
have ****:"int (n (i - Suc 0) - n i) - 1 = int (n (i - 1) - Suc (n i))"
using assms nlt by auto
from assms n_ge[of "i-2"] nlt n_ge[of i]
have nge:"n (i - Suc 0) ≤ n (i - 2)" "n i < n (i - Suc 0)" "n i < n (i - 1)" "Suc (i - 2) = i - 1"
by (cases i,auto simp:numeral_2_eq_2 numeral_3_eq_3)
have *****:"(- 1 :: 'a fract) ^ ((n (i - Suc 0) - n i) * (n (i - Suc 0) - n i + (n (i - 2) - n (Suc (i - 2)))))
= (- 1) ^ ((n i + n (i - Suc 0)) * (n i + n (i - 2)))"
"(- 1 :: 'a fract) ^ (∑l←[3..<i]. (n (l - Suc 0) - n i) * (n (l - 2) - n i))
= (- 1) ^ (∑l←[3..<i]. (n i + n (l - Suc 0)) * (n i + n (l - 2))) "
using nge apply (intro minus_1_even_eqI,auto)
apply (intro minus_1_even_eqI)
apply (intro even_sum_list)
proof(goal_cases) case (1 x)
with n_ge_trans assms
have "n i ≤ n (x - Suc 0)" "n (x - 2) ≥ n i" by auto
with 1 show ?case by auto
qed
have "ffp (R (n i)) = subresultant (n i) (F 1) (F 2)" unfolding R_def F1 F2
by (auto simp: to_fract_hom.subresultant_hom ni2[OF assms])
also have "… = smult
((- 1) ^ (∑l←[3..<i]. (n (l - 2) - n i) * (n (l - 1) - n i)) *
(∏x←[3..<i]. (β x / α x) ^ (n (x - 1) - n i) * f (x - 1) ^ (δ (x - 1) + δ (x - 2))) *
(((β i / α i) ^ δ (i - 1)) * f (i - 1) ^ (δ (i - 1) + δ (i - 2))) *
((- 1) ^ ((δ (i - 2) + δ (i - 1)) * δ (i - 1)) * f i ^ (δ (i - 1) - 1)
))
(F i)"
unfolding BT_eq_30'[OF assms(1) nlt] **
*[OF αnz BT_lemma_2_22[OF assms]] smult_smult by (auto simp:ac_simps *** )
also have "… = ?rhs" unfolding Θ_def τ_def
using prod_combine[OF assms(1)] δ assms
by (auto simp:ac_simps exp_pow_int[symmetric] power_add ***** ****)
finally show ?thesis.
qed
lemma fundamental_theorem_eq_7: assumes j: "j < n k" shows "R j = 0"
proof -
let ?P = "pseudo_divmod (F (k - 1)) (F k)"
from F0[of k] k2 have Fk: "F k ≠ 0" by auto
from pmod[of "Suc k"] k2 F0[of "Suc k"]
have "pseudo_mod (F (k - 1)) (F k) = 0" by auto
then obtain Q where "?P = (Q,0)"
unfolding pseudo_mod_def by (cases ?P, auto)
from pseudo_divmod(1)[OF Fk this] Fk obtain c where id: "smult c (F (k - 1)) = F k * Q"
and c: "c ≠ 0" by auto
from id have id: "smult c (F (k - 1)) = Q * F k" by auto
from n_ge[unfolded n, of "k - 1"] k2 c have "degree (F k) ≤ degree (smult c (F (k - 1)))" by auto
from subresultant_product[OF id this, unfolded subresultant_smult_left[OF c], of j] j
have *:"subresultant j (F (k + 1 - 2)) (F (k + 1 - 1)) = 0" using c unfolding n by simp
from assms have **:"j ≠ n 2"
by (meson k2 n_ge_trans not_le one_le_numeral order_refl)
from k2 assms have "3 ≤ k + 1" "k + 1 ≤ k + 1" "j < n (k + 1 - 1)" by auto
from BT_eq_30[OF this,unfolded *] nonzero_alphaprod[OF le_refl] ** F1 F2
show ?thesis by (auto simp:R_def F0 to_fract_hom.subresultant_hom[symmetric])
qed
definition "G i = R (n (i - 1) - 1)"
definition "H i = R (n i)"
lemma gamma_delta_beta_3: "γ 3 = (- 1) ^ (δ 1 + 1) * β 3"
proof -
have "γ 3 = (- 1) ^ σ 3 * pow_int (f 2) (1 - int (δ 2)) *
(β 3 / (f 2 ^ Suc (δ 1)) * f 2 ^ (δ 1 + δ 2))"
unfolding γ_def δ α_def by (simp add: δ)
also have "f 2 ^ (δ 1 + δ 2) = pow_int (f 2) (int (δ 1 + δ 2))"
unfolding pow_int_def nat_int by auto
also have "int (δ 1 + δ 2) = int (Suc (δ 1)) + (int (δ 2) - 1)" by simp
also have "pow_int (f 2) … = pow_int (f 2) (Suc (δ 1)) * pow_int (f 2) (int (δ 2) - 1)"
by (rule pow_int_add, insert f20, auto)
also have "pow_int (f 2) (Suc (δ 1)) = f 2 ^ (Suc (δ 1))" unfolding pow_int_def nat_int by simp
also have "β 3 / (f 2 ^ Suc (δ 1)) *
(f 2 ^ Suc (δ 1) * pow_int (f 2) (int (δ 2) - 1))
= (β 3 / (f 2 ^ Suc (δ 1)) * f 2 ^ Suc (δ 1) * pow_int (f 2) (int (δ 2) - 1))" by simp
also have "β 3 / (f 2 ^ Suc (δ 1)) * f 2 ^ Suc (δ 1) = β 3" using f20 by auto
finally have "γ 3 = ((- 1) ^ σ 3 * β 3) * (pow_int (f 2) (1 - int (δ 2)) * pow_int (f 2) (int (δ 2) - 1))"
by simp
also have "pow_int (f 2) (1 - int (δ 2)) * pow_int (f 2) (int (δ 2) - 1)
= 1"
by (subst pow_int_add[symmetric], insert f20, auto)
finally have "γ 3 = (- 1) ^ σ 3 * β 3" by simp
also have "σ 3 = (n 1 + n 2 + 1) * (n 2 + n 2 + 1)" unfolding σ_def
by simp
also have "(- (1 :: 'a fract)) ^ … = (- 1) ^ (n 1 - n 2 + 1)"
by (rule minus_1_even_eqI, insert n12, auto)
also have "… = (- 1)^(δ 1 + 1)" unfolding δ by (simp add: numeral_2_eq_2)
finally show "γ 3 = (- 1) ^ (δ 1 + 1) * β 3" .
qed
fun h :: "nat ⇒ 'a fract" where
"h i = (if (i ≤ 1) then 1 else if i = 2 then (f 2 ^ δ 1) else (f i ^ δ (i - 1) / (h (i - 1) ^ (δ (i - 1) - 1))))"
lemma smult_inverse_sdiv_poly: assumes ffp: "p ∈ range ffp"
and p: "p = smult (inverse x) q"
and p': "p' = sdiv_poly q' x'"
and xx: "x = ff x'"
and qq: "q = ffp q'"
shows "p = ffp p'"
proof (rule poly_eqI)
fix i
have "coeff p i = coeff q i / x" unfolding p by (simp add: field_simps)
also have "… = ff (coeff q' i) / ff x'" unfolding qq xx by simp
finally have cpi: "coeff p i = ff (coeff q' i) / ff x'" .
from ffp obtain r where pr: "p = ffp r" by auto
from arg_cong[OF this, of "λ p. coeff p i", unfolded cpi]
have "ff (coeff q' i) / ff x' ∈ range ff" by auto
hence id: "ff (coeff q' i) / ff x' = ff (coeff q' i div x')"
by (rule div_divide_to_fract, auto)
show "coeff p i = coeff (ffp p') i" unfolding cpi id p'
by (simp add: sdiv_poly_def coeff_map_poly)
qed
end
locale subresultant_prs_locale2 = subresultant_prs_locale F n δ f k β G1 G2 for
F :: "nat ⇒ 'a :: idom_divide fract poly"
and n :: "nat ⇒ nat"
and δ :: "nat ⇒ nat"
and f :: "nat ⇒ 'a fract"
and k :: nat
and β :: "nat ⇒ 'a fract"
and G1 G2 :: "'a poly" +
assumes β3: "β 3 = (-1)^(δ 1 + 1)"
and βi: "⋀ i. 4 ≤ i ⟹ i ≤ Suc k ⟹ β i = (-1)^(δ (i - 2) + 1) * f (i - 2) * h (i - 2) ^ (δ (i - 2))"
begin
lemma B_eq_17_main: "2 ≤ i ⟹ i ≤ k ⟹
h i = (-1) ^ (n 1 + n i + i + 1) / f i
* (∏l←[3..< Suc (Suc i)]. (α l / β l)) ∧ h i ≠ 0"
proof (induct i rule: less_induct)
case (less i)
from less(2-) have fi0: "f i ≠ 0" using f0[of i] by simp
have 1: "(- 1) ≠ (0 :: 'a fract)" by simp
show ?case (is "h i = ?r i ∧ _")
proof (cases "i = 2")
case True
have f20: "f 2 ≠ 0" using f20 by auto
have hi: "h i = f 2 ^ δ 1" unfolding True h.simps[of 2] by simp
have id: "int (δ 1) = int (n 1) - int (n 2)" using n12 unfolding δ numeral_2_eq_2 by simp
have "?r i = (- 1) ^ (1 + n 1 + n 2)
* ((f 2 ^ Suc (δ 1)) / (β 3)) / pow_int (f 2) 1" unfolding True α_def by simp
also have "β 3 = (- 1) ^ (δ 1 + 1)" by (rule β3)
also have "f 2 ^ Suc (δ 1) / … = … * f 2 ^ Suc (δ 1)" by simp
finally have "?r i = ((- 1) ^ (1 + n 1 + n 2) * ((- 1) ^ (δ 1 + 1))) *
pow_int (f 2) (int (Suc (δ 1)) + (-1))" (is "_ = ?a * _")
unfolding pow_int_divide exp_pow_int power_add pow_int_add[OF f20] by (simp add: ac_simps pow_int_add)
also have "?a = (-1)^(1 + n 1 + n 2 + δ 1 + 1)" unfolding power_add by simp
also have "… = (-1)^0"
by (rule minus_1_even_eqI, insert n12, auto simp: δ numeral_2_eq_2, presburger)
finally have ri: "?r i = pow_int (f 2) (int (δ 1))" by simp
show ?thesis unfolding ri hi exp_pow_int[symmetric] using f20 by simp
next
case False
hence i: "i ≥ 3" and ii: "i - 1 < i" "2 ≤ i - 1" "i - 1 ≤ k" using less(2-) by auto
from i less(2-) have cc: "4 ≤ Suc i" "Suc i ≤ Suc k" by auto
define P where "P = (∏l←[3..< Suc i]. α l / β l)"
define Q where "Q = P * pow_int (h (i - 1)) (- int (δ (i - 1)))"
define R where "R = f i ^ δ (i - 1)"
define S where "S = pow_int (f (i - 1)) (- 1)"
note IH = less(1)[OF ii]
hence hi0: "h (i - 1) ≠ 0" by auto
have hii: "h i = f i ^ δ (i - 1) / h (i - 1) ^ (δ (i - 1) - 1)"
unfolding h.simps[of i] using i by simp
also have "… = f i ^ δ (i - 1) * pow_int (h (i - 1)) (- int (δ (i - 1) - 1))"
unfolding exp_pow_int pow_int_divide by simp
also have "int (δ (i - 1) - 1) = int (δ (i - 1)) - 1"
proof -
have "δ (i - 1) > 0" unfolding δ[of "i - 1"] using n_gt[OF ii(2)] less(2-) by auto
thus ?thesis by simp
qed
also have "- (int (δ (i - 1)) - 1) = 1 + (- int (δ (i - 1)))" by simp
finally have hi: "h i = (- 1) ^ (n 1 + n (i - 1) + i) * (R * Q * S)"
unfolding pow_int_add[OF hi0] P_def Q_def pow_int_divide[symmetric] R_def S_def using IH i by (simp add: ac_simps)
from i have id: "[3..<Suc (Suc i)] = [3 ..< Suc i] @ [Suc i]" by simp
have "?r i = (- 1) ^ (n 1 + n i + i + 1)
* pow_int (f i) (- 1) * P * α (Suc i) / β (Suc i)"
unfolding pow_int_divide[symmetric] P_def id Fract_conv_to_fract by simp
also have "β (Suc i) = (- 1) ^ (δ (i - 1) + 1) * f (i - 1) * h (i - 1) ^ δ (i - 1)"
using βi[OF cc] by simp
also have "α (Suc i) = f i ^ Suc (δ (i - 1))" unfolding α_def by simp
finally have "?r i = (- 1) ^ (n 1 + n i + i + 1) * pow_int (f i) (- 1) * P * (f i ^ Suc (δ (i - 1))) /
(- 1) ^ (δ (i - 1) + 1) * pow_int (f (i - 1)) (- 1) / h (i - 1) ^ δ (i - 1)"
(is "_ = ?a1 * ?fi1 * P * ?fi2 / ?a2 * ?b / ?c")
unfolding exp_pow_int pow_int_divide[symmetric] by simp
also have "… = (?a1 / ?a2) * (?fi1 * ?fi2) * (P / ?c) * ?b" by (simp add: ac_simps)
also have "?a1 / ?a2 = (- 1) ^ (n 1 + n i + i + 1 + δ (i - 1) + 1)"
by (simp add: power_add)
also have "… = (-1) ^ (n 1 + n i + i + δ (i - 1))"
by (rule minus_1_even_eqI, auto)
also have "n 1 + n i + i + δ (i - 1) = n 1 + n (i - 1) + i"
unfolding δ using i less(2-) n_ge[of "i - 1"] by simp
also have "?fi1 * ?fi2 = pow_int (f i) (-1 + int (Suc (δ (i - 1))))"
unfolding exp_pow_int pow_int_add[OF fi0] by simp
also have "… = pow_int (f i) (int (δ (i - 1)))" by simp
also have "P / ?c = Q" unfolding Q_def exp_pow_int pow_int_divide by simp
also have "?b = S" unfolding S_def by simp
finally have ri: "?r i = (-1)^(n 1 + n (i - 1) + i)
* (R * Q * S)" by (simp add: exp_pow_int R_def)
have id: "h i = ?r i" unfolding hi ri ..
show ?thesis
by (rule conjI[OF id], unfold hii, insert IH fi0, auto)
qed
qed
lemma B_eq_17: "2 ≤ i ⟹ i ≤ k ⟹
h i = (-1) ^ (n 1 + n i + i + 1) / f i * (∏l←[3..< Suc (Suc i)]. (α l / β l))"
using B_eq_17_main by blast
lemma B_theorem_2: "3 ≤ i ⟹ i ≤ Suc k ⟹ γ i = 1"
proof (induct i rule: less_induct)
case (less i)
show ?case
proof (cases "i = 3")
case True
show ?thesis unfolding True unfolding gamma_delta_beta_3 β3 by simp
next
case False
with less(2-)
have i: "i ≥ 4" and ii: "i - 1 < i" "3 ≤ i - 1" "i - 1 ≤ Suc k"
and iii: "4 ≤ i" "i ≤ Suc k"
and iv: "2 ≤ i - 2" "i - 2 ≤ k" by auto
from less(1)[OF ii] have IH: "γ (i - 1) = 1" .
define L where "L = [3..< i]"
have id: "[3..<Suc (i - 1)] = L" "[3..<Suc i] = L @ [i]" "Suc (Suc (i - 2)) = i"
unfolding L_def using i by auto
define B where "B = (λ l. β l / α l)"
define A where "A = (λ l. α l / β l)"
define Q where "Q = (λ l. f (l - 1) ^ (δ (l - 2) + δ (l - 1)))"
define R where "R = (λ i l. B l ^ (n (l - 1) - n (i - 1) + 1))"
define P where "P = (λ i l. R i l * Q l)"
have fi0: "f (i - 1) ≠ 0" using f0[of "i - 1"] less(2-) by auto
have fi0': "f (i - 2) ≠ 0" using f0[of "i - 2"] less(2-) by auto
{
fix j
assume "j ∈ set L"
hence "j ≥ 3" "j < i" unfolding L_def by auto
with less(3) have j: "j - 1 ≠ 0" "j - 1 < k" by auto
hence Q: "Q j ≠ 0" unfolding Q_def using f0[of "j - 1"] by auto
from j α0 β0[of j] have 0: "α j ≠ 0" "β j ≠ 0" by auto
hence "B j ≠ 0" "A j ≠ 0" unfolding B_def A_def by auto
note Q this
} note L0 = this
let ?exp = "δ (i - 2)"
have "γ i = γ i / γ (i - 1)" unfolding IH by simp
also have "… = (- 1) ^ σ i * pow_int (f (i - 1)) (1 - int (δ (i - 1))) *
(∏l←L. P i l) * P i i /
((- 1) ^ σ (i - 1) * pow_int (f (i - 2)) (1 - int (δ (i - 2))) *
(∏l←L. P (i - 1) l))" (is "_ = ?a1 * ?f1 * ?L1 * P i i / (?a2 * ?f2 * ?L2)")
unfolding γ_def id P_def Q_def R_def B_def by (simp add: numeral_2_eq_2)
also have "… = (?a1 * ?a2) * (?f1 * P i i) / ?f2 * (?L1 / ?L2)" unfolding divide_prod_assoc by simp
also have "?a1 * ?a2 = (-1)^(σ i + σ (i - 1))" (is "_ = ?a") unfolding power_add by simp
also have "?L1 / ?L2 = (∏l←L. R i l) / (∏l←L. R (i - 1) l) * ((∏l←L. Q l) / (∏l←L. Q l))"
unfolding P_def prod_list_multf divide_prod_assoc by simp
also have "… = (∏l←L. R i l) / (∏l←L. R (i - 1) l)" (is "_ = ?L1 / ?L2")
proof -
have "(∏l←L. Q l) ≠ 0" unfolding prod_list_zero_iff using L0 by auto
thus ?thesis by simp
qed
also have "?f1 * P i i = (?f1 * pow_int (f (i - 1)) (int ?exp + int (δ (i - 1)))) * R i i" unfolding P_def Q_def
exp_pow_int by simp
also have "?f1 * pow_int (f (i - 1)) (int ?exp + δ (i - 1)) = pow_int (f (i - 1))
(1 + int ?exp)" (is "_ = ?f1")
unfolding pow_int_add[OF fi0, symmetric] by simp
also have "R i i = β i / α i" unfolding B_def R_def Fract_conv_to_fract by simp
also have "α i = f (i - 1) ^ Suc ?exp" unfolding α_def by simp
also have "β i / … = β i * pow_int (f (i - 1)) (- 1 - ?exp)"
(is "_ = ?β * ?f12")
unfolding exp_pow_int pow_int_divide by simp
finally have "γ i = (?a * (?f1 * ?f12)) * ?β / ?f2 * (?L1 / ?L2)"
by simp
also have "?a * (?f1 * ?f12) = ?a" unfolding pow_int_add[OF fi0, symmetric] by simp
also have "?L1 / ?L2 = pow_int (∏l←L. A l) (- ?exp)"
proof -
have id: "i - 1 - 1 = i - 2" by simp
have "set L ⊆ {l. 3 ≤ l ∧ l ≤ k ∧ l < i}" unfolding L_def using less(3) by auto
thus ?thesis unfolding R_def id
proof (induct L)
case (Cons l L)
from Cons(2) have l: "3 ≤ l" "l ≤ k" "l < i" and L: "set L ⊆ {l. 3 ≤ l ∧ l ≤ k ∧ l < i}" by auto
note IH = Cons(1)[OF L]
from l α0 β0[of l] have 0: "α l ≠ 0" "β l ≠ 0" by auto
hence B0: "B l ≠ 0" unfolding B_def by auto
have "(∏l←l # L. B l ^ (n (l - 1) - n (i - 1) + 1)) / (∏l←l # L. B l ^ (n (l - 1) - n (i - 2) + 1))
= (B l ^ (n (l - 1) - n (i - 1) + 1) * (∏l←L. B l ^ (n (l - 1) - n (i - 1) + 1))) /
(B l ^ (n (l - 1) - n (i - 2) + 1) * (∏l←L. B l ^ (n (l - 1) - n (i - 2) + 1)))"
(is "_ = (?l1 * ?L1) / (?l2 * ?L2)") by simp
also have "… = (?l1 / ?l2) * (?L1 / ?L2)" by simp
also have "?L1 / ?L2 = pow_int (prod_list (map A L)) (- int (δ (i - 2)))" by (rule IH)
also have "?l1 / ?l2 = pow_int (B l) (int (n (l - 1) - n (i - 1)) - int (n (l - 1) - n (i - 2)))" unfolding exp_pow_int pow_int_divide pow_int_add[OF B0, symmetric]
by simp
also have "int (n (l - 1) - n (i - 1)) - int (n (l - 1) - n (i - 2)) = int ?exp"
proof -
have "n (l - 1) ≥ n (i - 2)" "n (l - 1) ≥ n (i - 1)" "n (i - 2) ≥ n (i - 1)"
using i l less(3)
by (intro n_ge_trans, auto)+
hence id: "int (n (l - 1) - n (i - 1)) = int (n (l - 1)) - int (n (i - 1))"
"int (n (l - 1) - n (i - 2)) = int (n (l - 1)) - int (n (i - 2))"
"int (n (i - 2) - n (i - 1)) = int (n (i - 2)) - int (n (i - 1))"
by simp_all
have id2: "int ?exp = int (n (i - 2) - n (i - 1))"
unfolding δ using i by (cases i; cases "i - 1", auto)
show ?thesis unfolding id2 unfolding id by simp
qed
also have "pow_int (B l) … = pow_int (inverse (B l)) (- …)" unfolding pow_int_def
by (cases "int (δ (i - 2))" rule: linorder_cases, auto simp: field_simps)
also have "inverse (B l) = A l" unfolding B_def A_def by simp
also have "pow_int (A l) (- int ?exp) * pow_int (prod_list (map A L)) (- int ?exp)
= pow_int (prod_list (map A (l # L))) (- int ?exp)"
by (simp add: pow_int_mult)
finally show ?case .
qed simp
qed
also have "β i = (- 1) ^ (?exp + 1) * f (i - 2) * h (i - 2) ^ ?exp"
unfolding βi[OF iii] ..
finally have "γ i = (((- 1) ^ (σ i + σ (i - 1)) * (- 1) ^ (?exp + 1))) *
(pow_int (f (i - 2)) 1 *
pow_int (f (i - 2)) (int ?exp - 1)) *
h (i - 2) ^ ?exp /
(∏l←L. A l) ^ ?exp" (is "_ = ?a * ?f1 * ?H / ?L") unfolding pow_int_divide exp_pow_int by simp
also have "?f1 = pow_int (f (i - 2)) (int ?exp)" (is "_ = ?f1") unfolding pow_int_add[OF fi0', symmetric]
by simp
also have "h (i - 2) = (- 1) ^ (n 1 + n (i - 2) + (i - 2) + 1) / f (i - 2) *
(∏l←L. A l)" (is "_ = ?a2 / ?f2 * ?L") unfolding B_eq_17[OF iv] A_def id L_def by simp
also have "((- (1 :: 'a fract)) ^ (σ i + σ (i - 1)) * (- 1) ^ (?exp + 1)) =
((- 1) ^ (σ i + σ (i - 1) + ?exp + 1))" (is "_ = ?a1") by (simp add: power_add)
finally have "γ i = ?a1 * ?f1 * (?a2 / ?f2 * ?L) ^ ?exp / ?L ^ ?exp" by simp
also have "… = (?a1 * ?a2^?exp) * (?f1 / ?f2 ^ ?exp) * (?L^?exp / ?L ^ ?exp)"
unfolding power_mult_distrib power_divide by auto
also have " ?L ^ ?exp / ?L ^ ?exp = 1"
proof -
have "?L ≠ 0" unfolding prod_list_zero_iff using L0 by auto
thus ?thesis by simp
qed
also have "?f1 / ?f2 ^ ?exp = 1" unfolding exp_pow_int pow_int_divide
pow_int_add[OF fi0', symmetric] by simp
also have "?a2^?exp = (- 1) ^ ((n 1 + n (i - 2) + (i - 2) + 1) * ?exp)"
by (rule semiring_normalization_rules)
also have "?a1 * … = (- 1) ^ (σ i + σ (i - 1) + ?exp + 1 + (n 1 + n (i - 2) + (i - 2) + 1) * ?exp)"
(is "_ = _ ^ ?e")
by (simp add: power_add)
also have "… = (-1)^0"
proof -
define e where "e = ?e"
have *: "?e = (2 * ?exp + σ i + σ (i - 1) + 1 + (n 1 + n (i - 2) + (i - 2)) * ?exp)" by simp
define A where "A = (λ i l. (n (l - 2) + n (i - 1) + 1) * (n (l - 1) + n (i - 1) + 1))"
define B where "B = (λ i. (n (i - 1) + 1) * (n (i - 1) + 1))"
define C where "C = (λ l. (n (l - 1) + n (l - 2) + n (l - 1) * n (l - 2)))"
define D where "D = (λ l. n (l - 1) + n (l - 2))"
define m2 where "m2 = n (i - 2)"
define m1 where "m1 = n (i - 1)"
define m0 where "m0 = n 1"
define i3 where "i3 = i - 3"
have m12: "m2 ≥ m1" unfolding m2_def m1_def using n_ge[of "i - 2"] i less(3)
by (cases i, auto)
have idd: "Suc (i - 2) = i - 1" "i - 1 - 1 = i - 2" using i by auto
have id4: "i - 2 = Suc i3" unfolding i3_def using i by auto
from i have "3 < i" by auto
hence "∃ k. sum_list (map D L) = n 1 + n (i - 2) + 2 * k" unfolding L_def
proof (induct i rule: less_induct)
case (less i)
show ?case
proof (cases "i = 4")
case True
thus ?thesis by (simp add: D_def)
next
case False
obtain ii where i: "i = Suc ii" and ii: "ii < i" "3 < ii" using False less(2) by (cases i, auto)
from less(1)[OF ii] obtain k where IH: "sum_list (map D [3 ..< ii]) = n 1 + n (ii - 2) + 2 * k" by auto
have "map D [3 ..< i] = map D [3 ..< ii] @ [D ii]" unfolding i using ii by auto
hence "sum_list (map D [3..<i]) = n 1 + n (ii - 2) + 2 * k + D ii" using IH by simp
also have "… = n 1 + n (ii - 1) + 2 * (n (ii - 2) + k)" unfolding D_def by simp
also have "n (ii - 1) = n (i - 2)" unfolding i by simp
finally show ?thesis by blast
qed
qed
then obtain kk where DL: "sum_list (map D L) = n 1 + n (i - 2) + 2 * kk" ..
let ?l = "i - 3"
have len: "length L = i - 3" unfolding L_def using i by auto
have A: "A i l = B i + D l * n (i - 1) + C l" for i l
unfolding A_def B_def C_def D_def ring_distribs by simp
have id2: "[3..<Suc i] = 3 # [Suc 3 ..< Suc i]"
unfolding L_def using i by (auto simp: upt_rec[of 3])
have "even e = even ?e" unfolding e_def by simp
also have "… = even ((1 + (n 1 + n (i - 2) + (i - 2)) * ?exp) + (σ i + σ (i - 1)))"
(is "_ = even (?g + ?j)")
unfolding * by (simp add: ac_simps)
also have "?j = (∑l←L @ [i]. A i l) + (∑l←L. A (i - 1) l)"
unfolding σ_def id A_def by simp
also have "… = 2 * (∑l←L. C l) + (Suc ?l) * B i + (∑l←L @ [i]. D l * n (i - 1)) + C i +
?l * B (i - 1) + (∑l←L. D l * n (i - 1 - 1))"
unfolding A sum_list_addf by (simp add: sum_list_triv len)
also have "… = ((Suc ?l * B i + C i +
?l * B (i - 1) + D i * n (i - 1)) + ((∑l←L. D l) * (n (i - 1) + n (i - 2)) + 2 * (∑l←L. C l)))"
(is "_ = ?i + ?j")
unfolding sum_list_mult_const by (simp add: ring_distribs numeral_2_eq_2)
also have "?j =
(n 1 + n (i - 2)) * (n (i - 1) + n (i - 2)) + 2 * (kk * (n (i - 1) + n (i - 2)) + (∑l←L. C l))"
(is "_ = ?h + 2 * ?f")
unfolding DL by (simp add: ring_distribs)
finally have "even e = even (?g + ?i + ?h + 2 * ?f)" by presburger
also have "… = even (?g + ?i + ?h)" by presburger
also have "?g + ?i + ?h =
i3 * (m2 - m1 + m1 * m1 + m2 * m2)
+ (m2 - m1 + m1 + m2) * (m0 + m2)
+ (m1 + m2 + (m2 - m1))
+ 2 * (m1 * m2 + m1 * m1 + 1 + i3 + m1 * Suc i3 + m2 * i3)" unfolding idd B_def D_def C_def δ
m1_def[symmetric] m2_def[symmetric] m0_def[symmetric]
unfolding i3_def[symmetric] id4
by (simp add: ring_distribs)
also have "(m1 + m2 + (m2 - m1)) = 2 * m2" using m12 by simp
also have "(m2 - m1 + m1 + m2) * (m0 + m2) = 2 * (m2 * (m0 + m2))" using m12 by simp
finally obtain l1 l2 l3 where
"even e = even (i3 * (m2 - m1 + m1 * m1 + m2 * m2) + 2 * l1 + 2 * l2 + 2 * l3)"
by blast
also have "… = even (i3 * (m2 - m1 + m1 * m1 + m2 * m2))" by simp
also have "… = even (i3 * (2 * m1 + (m2 - m1 + m1 * m1 + m2 * m2)))" by simp
also have "2 * m1 + (m2 - m1 + m1 * m1 + m2 * m2) = m1 + m2 + m1 * m1 + m2 * m2"
using m12 by simp
also have "even (i3 * …)" by auto
finally have "even e" .
thus ?thesis unfolding e_def
by (intro minus_1_even_eqI, auto)
qed
finally show "γ i = 1" by simp
qed
qed
context
fixes i :: nat
assumes i: "3 ≤ i" "i ≤ k"
begin
lemma B_theorem_3_b: "Θ i * f i = ff (lead_coeff (H i))"
using arg_cong[OF fundamental_theorem_eq_6[folded H_def, OF i], of lead_coeff] unfolding f[of i]
lead_coeff_smult by simp
lemma B_theorem_3_main: "Θ i * f i / γ (i + 1) = (-1)^(n 1 + n i + i + 1) / f i * (∏l←[3..< Suc (Suc i)]. (α l / β l))"
proof (cases "f i = 0")
case True
thus ?thesis by simp
next
case False note ff0 = this
from i(1) have "Suc (Suc i) > 3" by auto
hence id: "[3 ..< Suc (i + 1)] = [3 ..< Suc i] @ [Suc i]" "[3 ..< Suc (Suc i)] = [3 ..< Suc i] @ [Suc i]" by auto
have cong: "⋀ a b c d. a = c ⟹ b = d ⟹ a * b = c * (d :: 'a fract)" by auto
define AB where "AB = (λ l. β l / α l)"
define ABP where "ABP = (λ l. AB l ^ (n (l - 1) - n i) * f (l - 1) ^ (δ (l - 2) + δ (l - 1)))"
define PR where "PR = (∏l←[3..<Suc i]. ABP l)"
define PR2 where "PR2 = (∏l←[3..<Suc i]. AB l)"
from F0[of i]
have "Θ i * f i / γ (i + 1) = (
((- 1) ^ τ i * (- 1) ^ σ (i + 1)) * (pow_int (f i) (int (δ (i - 1)) - 1) *
PR * f i /
pow_int (f i) (1 - int (δ i)) / ((∏l←[3..<Suc i]. ABP l * AB l) *
AB (Suc i) * f i ^ (δ (i - 1) + δ i))))"
unfolding id prod_list.append map_append Θ_def γ_def divide_prod_assoc
by (simp add: field_simps power_add AB_def ABP_def PR_def)
also have "(- 1 :: 'a fract) ^ τ i * (- 1) ^ σ (i + 1) = (- 1) ^ (τ i + σ (i + 1))"
unfolding power_add by (auto simp: field_simps)
also have "… = (- 1) ^ (n 1 + n i + i + 1)"
proof (cases "i = 2")
case True
show ?thesis unfolding τ_def σ_def True by (auto, rule minus_1_even_eqI, auto)
next
case False
define a where "a = (λ l. n (l - 2) + n i)"
define b where "b = (λ l. n (l - 1) + n i)"
define c where "c = (∑l←[3..<Suc i]. (a l * b l + n i))"
define d where "d = c + (∑l←[3..<i]. n (l - 1))"
define e where "e = (n (i - 1) + n i + 1) * n i"
have "(τ i + σ (i + 1)) =
((∑l←[3..<Suc i]. (a l * b l) + (a l + 1) * (b l + 1))) + (a (Suc i) + 1) * (b (Suc i) + 1)"
unfolding σ_def τ_def id a_def b_def sum_list_addf by simp
also have "(∑l←[3..<Suc i]. (a l * b l) + (a l + 1) * (b l + 1)) =
(∑l←[3..<Suc i]. 2 * a l * b l + (a l + b l) + 1)"
by (rule arg_cong, rule map_cong, auto)
also have "… = (∑l←[3..<Suc i]. 2 * (a l * b l + n i) + (n (l - 1) + n (l - 2)) + 1)"
by (simp add: field_simps a_def b_def)
also have "… = 2 * c + (∑l←[3..<Suc i]. (n (l - 1) + n (l - 2))) + length [3 ..< Suc i]"
unfolding sum_list_addf c_def sum_list_const_mult sum_list_triv by simp
also have "(∑l←[3..<Suc i]. (n (l - 1) + n (l - 2)))
= (∑l←[3..<Suc i]. n (l - 1)) + (∑l←[3..<Suc i]. n (l - 2))"
by (simp add: sum_list_addf)
also have "(∑l←[3..<Suc i]. n (l - 2)) = (∑l←3 # [4..<Suc i]. n (l - 2))"
by (rule arg_cong, rule map_cong, insert i False, auto simp: upt_rec[of 3])
also have "… = n 1 + (∑l←[(Suc 3)..<Suc i]. n (l - 2))" by auto
also have "(∑l←[(Suc 3)..<Suc i]. n (l - 2)) = (∑l←[3..<i]. n (l - 1))"
proof (rule arg_cong[of _ _ sum_list], rule nth_equalityI, force, auto simp: nth_append, goal_cases)
case (1 j)
hence "i - 2 = Suc (Suc j)" by simp
thus ?case by simp
qed
also have "(∑l←[3..<Suc i]. n (l - 1)) = (∑l←[3..<i] @ [i]. n (l - 1))"
by (rule arg_cong, rule map_cong, insert i False, auto)
finally have "τ i + σ (i + 1) =
2 * d + n (i - 1) + n 1 + length [3..<Suc i] + (a (Suc i) + 1) * (b (Suc i) + 1)"
by (simp add: d_def)
also have "length [3 ..< Suc i] = i - 2" using i by auto
also have "(a (Suc i) + 1) * (b (Suc i) + 1) = 2 * e + n (i - 1) + n i + 1" unfolding a_def b_def e_def
by simp
finally have id: "τ i + σ (i + 1) = 2 * (d + n (i - 1) + e) + n 1 + (i - 2) + n i + 1"
by simp
show ?thesis
by (rule minus_1_even_eqI, unfold id, insert i, auto)
qed
also have "(∏l←[3..<Suc i]. ABP l * AB l) = PR * PR2"
unfolding PR_def prod_list_multf PR2_def by simp
also have "(pow_int (f i) (int (δ (i - 1)) - 1) * PR * f i / pow_int (f i) (1 - int (δ i))
/ (PR * PR2 * AB (Suc i) * f i ^ (δ (i - 1) + δ i))) =
((pow_int (f i) (int (δ (i - 1)) - 1) * pow_int (f i) 1 * pow_int (f i) (int (δ i) - 1)
/ pow_int (f i) (int (δ (i - 1) + δ i)))) * (PR / PR / (PR2 * AB (Suc i)))"
(is "… = ?x * ?y")
unfolding exp_pow_int[symmetric] by (simp add: pow_int_divide ac_simps)
also have "?x = pow_int (f i) (-1)"
unfolding pow_int_divide pow_int_add[OF ff0, symmetric] by simp
also have "… = 1 / (f i)"
unfolding pow_int_def by simp
also have "PR / PR = 1"
proof -
have "PR ≠ 0" unfolding PR_def prod_list_zero_iff set_map
proof
assume "0 ∈ ABP ` set [3 ..< Suc i]"
then obtain j where j: "3 ≤ j" "j < Suc i" and 0: "ABP j = 0" by auto
with i have jk: "j ≤ k" and j1: "j - 1 ≠ 0" "j - 1 < k" by auto
hence 1: "α j ≠ 0" "f (j - 1) ≠ 0" using α0 f0 by auto
with 0 have "AB j = 0" unfolding ABP_def by simp
from this[unfolded AB_def] 1(1) β0[of j] show False by auto
qed
thus ?thesis by simp
qed
also have "PR2 * AB (Suc i) = (∏l←[3..<Suc (Suc i)]. AB l)" unfolding id PR2_def by auto
also have "1 / … = inverse …" by (simp add: inverse_eq_divide)
also have "… = (∏l←[3..<Suc (Suc i)]. α l / β l)" unfolding AB_def
inverse_prod_list map_map o_def
by (auto cong: map_cong)
finally show ?thesis by simp
qed
lemma B_theorem_3: "h i = Θ i * f i" "h i = ff (lead_coeff (H i))"
proof -
have "Θ i * f i = Θ i * f i / γ (i + 1)"
using B_theorem_2[of "i + 1"] i by auto
also have "… = (- 1) ^ (n 1 + n i + i + 1) / f i *
(∏l←[3..<Suc (Suc i)]. α l / β l)" by (rule B_theorem_3_main)
also have "… = h i" using B_eq_17[of i] i by simp
finally show "h i = Θ i * f i" ..
thus "h i = ff (lead_coeff (H i))" using B_theorem_3_b by auto
qed
end
lemma h0: "i ≤ k ⟹ h i ≠ 0"
proof (induct i)
case (Suc i)
thus ?case unfolding h.simps[of "Suc i"] using f0 by (auto simp del: h.simps)
qed auto
lemma deg_G12: "degree G1 ≥ degree G2" using n12
unfolding n F1 F2 by auto
lemma R0: shows "R 0 = [: resultant G1 G2 :]"
proof(cases "n 2 = 0")
case True
hence d:"degree G2 = 0" unfolding n F2 by auto
from degree0_coeffs[OF d] F2 F12 obtain a where
G2: "G2 = [:a:]" and a: "a ≠ 0" by auto
have "sdiv_poly [:a * a ^ degree G1:] a = [:a ^ degree G1:]" using a
unfolding sdiv_poly_def by auto
note dp = this
show ?thesis using G2 F12
unfolding R_def δ n F1 F2 Suc_1 by (auto split:if_splits simp:mult.commute dp)
next
case False
from False n12 have d:"degree G2 ≠ 0" "degree G2 ≤ degree G1" unfolding n F2 F1 by auto
from False have "R 0 = subresultant 0 G1 G2" unfolding R_def by simp
also have "… = [: resultant G1 G2 :]" unfolding subresultant_resultant by simp
finally show ?thesis .
qed
context
fixes div_exp :: "'a ⇒ 'a ⇒ nat ⇒ 'a"
assumes div_exp_sound: "div_exp_sound div_exp"
begin
interpretation div_exp_sound div_exp by (rule div_exp_sound)
lemma subresultant_prs_main: assumes "subresultant_prs_main Gi_1 Gi hi_1 = (Gk, hk)"
and "F i = ffp Gi"
and "F (i - 1) = ffp Gi_1"
and "h (i - 1) = ff hi_1"
and "i ≥ 3" "i ≤ k"
shows "F k = ffp Gk ∧ h k = ff hk ∧ (∀ j. i ≤ j ⟶ j ≤ k ⟶ F j ∈ range ffp ∧ β (Suc j) ∈ range ff)"
proof -
obtain m where m: "m = k - i" by auto
show ?thesis using m assms
proof (induct m arbitrary: Gi_1 Gi hi_1 i rule: less_induct)
case (less m Gi_1 Gi hi_1 i)
note IH = less(1)
note m = less(2)
note res = less(3)
note id = less(4-6)
note i = less(7-8)
let ?pmod = "pseudo_mod Gi_1 Gi"
let ?ni = "degree Gi"
let ?ni_1 = "degree Gi_1"
let ?gi = "lead_coeff Gi"
let ?gi_1 = "lead_coeff Gi_1"
let ?d1 = "?ni_1 - ?ni"
obtain hi where hi: "hi = div_exp ?gi hi_1 ?d1" by auto
obtain divisor where div: "divisor = (-1) ^ (?d1 + 1) * ?gi_1 * (hi_1 ^ ?d1)" by auto
obtain G1_p1 where G1_p1: "G1_p1 = sdiv_poly ?pmod divisor" by auto
note res = res[unfolded subresultant_prs_main.simps[of Gi_1] Let_def,
folded hi, folded div, folded G1_p1]
have h_i: "h i = f i ^ δ (i - 1) / h (i - 1) ^ (δ (i - 1) - 1)" unfolding h.simps[of i] using i by simp
have hi_ff: "h i ∈ range ff" using B_theorem_3[OF _ i(2)] i by auto
have d1: "δ (i - 1) = ?d1" unfolding δ n using id(1,2) using i by simp
have fi: "f i = ff ?gi" unfolding f id by simp
have fi1: "f (i - 1) = ff ?gi_1" unfolding f id by simp
have eq': "h i = ff (lead_coeff Gi) ^ δ (i - 1) / ff hi_1 ^ (δ (i - 1) - 1)"
unfolding h_i fi id ..
have idh: "h i = ff hi" using hi_ff h_i fi id
unfolding hi d1[symmetric]
by (subst div_exp[of ?gi "δ (i - 1)" hi_1], unfold eq'[symmetric], insert assms, blast+)
have "β (Suc i) = (- 1) ^ (δ (i - 1) + 1) * f (i - 1) * h (i - 1) ^ δ (i - 1)"
using βi[of "Suc i"] i by auto
also have "… = ff ((- 1) ^ (δ (i - 1) + 1) * lead_coeff Gi_1 * hi_1 ^ δ (i - 1))"
unfolding id f by (simp add: hom_distribs)
also have "… ∈ range ff" by blast
finally have beta: "β (Suc i) ∈ range ff" .
have pm: "pseudo_mod (F (i - 1)) (F i) = ffp ?pmod"
unfolding to_fract_hom.pseudo_mod_hom[symmetric] id by simp
have eq: "(?pmod = 0) = (i = k)"
using pm i pmod[of "Suc i"] F0[of "Suc i"] i β0[of "Suc i"] by auto
show ?case
proof (cases "i = k")
case True
with res eq have res: "Gk = Gi" "hk = hi" by auto
with pmod
have "F k = ffp Gk ∧ h k = ff hk" unfolding res idh[symmetric] id[symmetric] True by auto
thus ?thesis using beta unfolding True by auto
next
case False
with res eq have res:
"subresultant_prs_main Gi G1_p1 hi = (Gk, hk)" by auto
from m False i have m: "m - 1 < m" "m - 1 = k - Suc i" by auto
have si: "Suc i - 1 = i" and ii: "3 ≤ Suc i" "Suc i ≤ k" and iii: "3 ≤ Suc i" "Suc i ≤ Suc k"
using False i by auto
have *: "(∀j≥Suc i. j ≤ k ⟶ F j ∈ range ffp ∧ β (Suc j) ∈ range ff) = (∀j≥i. j ≤ k ⟶ F j ∈ range ffp ∧ β (Suc j) ∈ range ff)"
by (rule for_all_Suc, insert id(1) beta, auto)
show ?thesis
proof (rule IH[OF m res, unfolded si, OF _ id(1) idh ii, unfolded *])
have F_ffp: "F (Suc i) ∈ range ffp" using fundamental_theorem_eq_4[OF ii, symmetric] B_theorem_2[OF iii] by auto
from pmod[OF iii] have "smult (β (Suc i)) (F (Suc i)) = pseudo_mod (F (i - 1)) (F i)"
by simp
from arg_cong[OF this, of "λ x. smult (inverse (β (Suc i))) x"]
have Fsi: "F (Suc i) = smult (inverse (β (Suc i))) (pseudo_mod (F (i - 1)) (F i))"
using β0[of "Suc i"] by auto
show "F (Suc i) = ffp G1_p1"
proof (rule smult_inverse_sdiv_poly[OF F_ffp Fsi G1_p1 _ pm])
from i ii have iv: "4 ≤ Suc i" "Suc i ≤ Suc k" by auto
have *: "Suc i - 2 = i - 1" by auto
show "β (Suc i) = ff divisor" unfolding βi[OF iv] div d1 * fi1
using id by (simp add: hom_distribs)
qed
qed
qed
qed
qed
lemma subresultant_prs: assumes res: "subresultant_prs G1 G2 = (Gk, hk)"
shows "F k = ffp Gk ∧ h k = ff hk ∧ (i ≠ 0 ⟶ F i ∈ range ffp) ∧ (3 ≤ i ⟶ i ≤ Suc k ⟶ β i ∈ range ff)"
proof -
let ?pmod = "pseudo_mod G1 G2"
have pm: "pseudo_mod (F 1) (F 2) = ffp ?pmod"
unfolding to_fract_hom.pseudo_mod_hom[symmetric] F1 F2 by simp
let ?g2 = "lead_coeff G2"
let ?n2 = "degree G2"
obtain d1 where d1: "d1 = degree G1 - ?n2" by auto
obtain h2 where h2: "h2 = ?g2 ^ d1" by auto
have "(?pmod = 0) = (pseudo_mod (F 1) (F 2) = 0)" using pm by auto
also have "… = (k < 3)" using k2 pmod[of 3] F0[of 3] β0[of 3] by auto
finally have eq: "?pmod = 0 ⟷ k = 2" using k2 by linarith
note res = res[unfolded subresultant_prs_def Let_def eq, folded d1, folded h2]
have idh2: "h 2 = ff h2" unfolding h2 d1 h.simps[of 2] δ n F1
using F2 by (simp add: numeral_2_eq_2 f hom_distribs)
have main: "F k = ffp Gk ∧ h k = ff hk ∧ (i ≥ 3 ⟶ i ≤ k ⟶ F i ∈ range ffp ∧ β (Suc i) ∈ range ff)" for i
proof (cases "k = 2")
case True
with res have "Gk = G2" "hk = h2" by auto
thus ?thesis using True idh2 F2 by auto
next
case False
hence "(k = 2) = False" by simp
note res = res[unfolded this if_False]
have F_2: "F (3 - 1) = ffp G2" using F2 by simp
have h2: "h (3 - 1) = ff h2" using idh2 by simp
have n2: "degree G2 = n (3 - 1)" unfolding n using F2 by simp
from False k2 have k3: "3 ≤ k" by auto
have "F k = ffp Gk ∧ h k = ff hk ∧ (∀j≥3. j ≤ k ⟶ F j ∈ range ffp ∧ β (Suc j) ∈ range ff)"
proof (rule subresultant_prs_main[OF res _ F_2 h2 le_refl k3])
let ?pow = "(- 1) ^ (δ 1 + 1) :: 'a fract"
from pmod[of 3] k3
have "smult (β 3) (F 3) = pseudo_mod (F 1) (F 2)" by simp
also have "… = pseudo_mod (ffp G1) (ffp G2)" using F1 F2 by auto
also have "… = ffp (pseudo_mod G1 G2)" unfolding to_fract_hom.pseudo_mod_hom by simp
also have "β 3 = (- 1) ^ (δ 1 + 1)" unfolding β3 by simp
finally have "smult ((- 1) ^ (δ 1 + 1)) (F 3) = ffp (pseudo_mod G1 G2)" by simp
also have "smult ((- 1) ^ (δ 1 + 1)) (F 3) = [: ?pow :] * F 3"
by simp
also have "[: ?pow :] = (- 1) ^ (δ 1 + 1)" by (unfold hom_distribs, simp)
finally have "(- 1) ^ (δ 1 + 1) * F 3 = ffp (pseudo_mod G1 G2)" by simp
from arg_cong[OF this, of "λ i. (- 1) ^ (δ 1 + 1) * i"]
have "F 3 = (- 1) ^ (δ 1 + 1) * ffp (pseudo_mod G1 G2)" by simp
also have "δ 1 = d1" unfolding δ n d1 using F1 F2 by (simp add: numeral_2_eq_2)
finally show F3: "F 3 = ffp ((- 1) ^ (d1 + 1) * pseudo_mod G1 G2)" by (simp add: hom_distribs)
qed
thus ?thesis by auto
qed
show ?thesis
proof (intro conjI impI)
assume "i ≠ 0"
then consider (12) "i = 1 ∨ i = 2" | (i3) "i ≥ 3 ∧ i ≤ k" | (ik) "i > k" by linarith
thus "F i ∈ range ffp"
proof cases
case 12
thus ?thesis using F1 F2 by auto
next
case i3
thus ?thesis using main by auto
next
case ik
hence "F i = 0" using F0 by auto
thus ?thesis by simp
qed
next
assume "3 ≤ i" and "i ≤ Suc k"
then consider (3) "i = 3" | (4) "3 ≤ i - 1" "i - 1 ≤ k" by linarith
thus "β i ∈ range ff"
proof (cases)
case 3
have "β i = ff ((- 1) ^ (δ 1 + 1))" unfolding 3 β3 by (auto simp: hom_distribs)
thus ?thesis by blast
next
case 4
with main[of "i - 1"] show ?thesis by auto
qed
qed (insert main, auto)
qed
lemma resultant_impl_main: "resultant_impl_main G1 G2 = resultant G1 G2"
proof -
from F0[of 2] F12(2) have k2: "k ≥ 2" by auto
obtain Gk hk where sub: "subresultant_prs G1 G2 = (Gk, hk)" by force
from subresultant_prs[OF this] have *: "F k = ffp Gk" "h k = ff hk" by auto
have "resultant_impl_main G1 G2 = (if degree (F k) = 0 then hk else 0)"
unfolding resultant_impl_main_def sub split * using F2 F12 by auto
also have "… = resultant G1 G2"
proof (cases "n k = 0")
case False
with fundamental_theorem_eq_7[of 0] show ?thesis unfolding n[of k] * R0 by auto
next
case True
from H_def[of k, unfolded True] have R: "R 0 = H k" by simp
show ?thesis
proof (cases "k = 2")
case False
with k2 have k3: "k ≥ 3" by auto
from B_theorem_3[OF k3] R0 R have "h k = ff (resultant G1 G2)" by simp
from this[folded *] * have "hk = resultant G1 G2" by simp
with True show ?thesis unfolding n by auto
next
case 2: True
have id: "(if degree (F k) = 0 then hk else 0) = hk" using True unfolding n by simp
from F0[of 3, unfolded 2] have "F 3 = 0" by simp
with pmod[of 3, unfolded 2] β0[of 3] have "pseudo_mod (F 1) (F 2) = 0" by auto
hence pm: "pseudo_mod G1 G2 = 0" unfolding F1 F2 to_fract_hom.pseudo_mod_hom by simp
from subresultant_prs_def[of G1 G2, unfolded sub Let_def this]
have id: "Gk = G2" "hk = lead_coeff G2 ^ (degree G1 - degree G2)" by auto
from F12 F1 F2 have "G1 ≠ 0" "G2 ≠ 0" by auto
from resultant_pseudo_mod_0[OF pm deg_G12 this]
have res: "resultant G1 G2 = (if degree G2 = 0 then lead_coeff G2 ^ degree G1 else 0)"
by simp
from True[unfolded 2 n F2] have "degree G2 = 0" by simp
thus ?thesis unfolding res 2 F2 id by simp
qed
qed
finally show ?thesis .
qed
end
end
text ‹At this point, we have soundness of the resultant-implementation, provided that we can
instantiate the locale by constructing suitable values of F, b, h, etc. Now we show the
existence of suitable locale parameters by constructively computing them.›
context
fixes G1 G2 :: "'a :: idom_divide poly"
begin
private function F and b and h where "F i = (if i = (0 :: nat) then 1
else if i = 1 then map_poly to_fract G1 else if i = 2 then map_poly to_fract G2
else (let G = pseudo_mod (F (i - 2)) (F (i - 1))
in if F (i - 1) = 0 ∨ G = 0 then 0 else smult (inverse (b i)) G))"
| "b i = (if i ≤ 2 then 1 else
if i = 3 then (- 1) ^ (degree (F 1) - degree (F 2) + 1)
else if F (i - 2) = 0 then 1 else (- 1) ^ (degree (F (i - 2)) - degree (F (i - 1)) + 1) * lead_coeff (F (i - 2)) *
h (i - 2) ^ (degree (F (i - 2)) - degree (F (i - 1))))"
| "h i = (if (i ≤ 1) then 1 else if i = 2 then (lead_coeff (F 2) ^ (degree (F 1) - degree (F 2))) else
if F i = 0 then 1 else (lead_coeff (F i) ^ (degree (F (i - 1)) - degree (F i)) / (h (i - 1) ^ ((degree (F (i - 1)) - degree (F i)) - 1))))"
by pat_completeness auto
termination
proof
show "wf (measure (case_sum (λ fi. 3 * fi +1) (case_sum (λ bi. 3 * bi) (λ hi. 3 * hi + 2))))" by simp
qed (auto simp: termination_simp)
declare h.simps[simp del] b.simps[simp del] F.simps[simp del]
private lemma Fb0: assumes base: "G1 ≠ 0" "G2 ≠ 0"
shows "(F i = 0 ⟶ F (Suc i) = 0) ∧ b i ≠ 0 ∧ h i ≠ 0"
proof (induct i rule: less_induct)
case (less i)
note * [simp] = F.simps[of i] b.simps[of i] h.simps[of i]
consider (0) "i = 0" | (1) "i = 1" | (2) "i ≥ 2" by linarith
thus ?case
proof cases
case 0
show ?thesis unfolding * unfolding 0 by simp
next
case 1
show ?thesis unfolding * unfolding 1 using assms by simp
next
case 2
have F: "F i = 0 ⟹ F (Suc i) = 0" unfolding F.simps[of "Suc i"] using 2 by simp
from assms have F2: "F 2 ≠ 0" unfolding F.simps[of 2] by simp
from 2 have "i - 1 < i" "i - 2 < i" by auto
note IH = less[OF this(1)] less[OF this(2)]
hence b: "b (i - 1) ≠ 0" and h: "h (i - 1) ≠ 0" "h (i - 2) ≠ 0" by auto
from h have hi: "h i ≠ 0" unfolding h.simps[of i] using 2 F2 by auto
have bi: "b i ≠ 0" unfolding b.simps[of i] using h(2) by auto
show ?thesis using hi bi F by blast
qed
qed
private definition "k = (LEAST i. F (Suc i) = 0)"
private lemma k_exists: "∃ i. F (Suc i) = 0"
proof -
obtain n i where "i ≥ 3" "length (coeffs (F (Suc i))) = n" by blast
thus ?thesis
proof (induct n arbitrary: i rule: less_induct)
case (less n i)
let ?ii = "Suc (Suc i)"
let ?i = "Suc i"
from less(2) have i: "?i ≥ 3" by auto
let ?mod = "pseudo_mod (F (?ii - 2)) (F ?i)"
have Fi: "F ?ii = (if F ?i = 0 ∨ ?mod = 0 then 0 else smult (inverse (b ?ii)) ?mod)"
unfolding F.simps[of ?ii] using i by auto
show ?case
proof (cases "F ?ii = 0")
case False
hence Fi: "F ?ii = smult (inverse (b ?ii)) ?mod" and mod: "?mod ≠ 0" and Fi1: "F ?i ≠ 0"
unfolding Fi by auto
from pseudo_mod[OF Fi1, of "F (?ii - 2)"] mod have "degree ?mod < degree (F ?i)" by simp
hence deg: "degree (F ?ii) < degree (F ?i)" unfolding Fi by auto
hence "length (coeffs (F ?ii)) < length (coeffs (F ?i))" unfolding degree_eq_length_coeffs by auto
from less(1)[OF _ i refl, folded less(3), OF this] show ?thesis by auto
qed blast
qed
qed
private lemma k: "F (Suc k) = 0" "i < k ⟹ F (Suc i) ≠ 0"
proof -
show "F (Suc k) = 0" unfolding k_def using k_exists by (rule LeastI2_ex)
assume "i < k" from not_less_Least[OF this[unfolded k_def]] show "F (Suc i) ≠ 0" .
qed
lemma enter_subresultant_prs: assumes len: "length (coeffs G1) ≥ length (coeffs G2)"
and G2: "G2 ≠ 0"
shows "∃ F n d f k b. subresultant_prs_locale2 F n d f k b G1 G2"
proof (intro exI)
from G2 len have G1: "G1 ≠ 0" by auto
from len have deg_le: "degree (F 2) ≤ degree (F 1)"
by (simp add: F.simps degree_eq_length_coeffs)
from G2 G1 have F1: "F 1 ≠ 0" and F2: "F 2 ≠ 0" by (auto simp: F.simps)
note Fb0 = Fb0[OF G1 G2]
interpret s: subresultant_prs_locale F "λ i. degree (F i)" "λ i. degree (F i) - degree (F (Suc i))"
"λ i. lead_coeff (F i)" k b G1 G2
proof (unfold_locales, rule refl, rule refl, rule refl, rule deg_le, rule F1, rule F2)
from k(1) F1 have k0: "k ≠ 0" by (cases k, auto)
show Fk: "(F i = 0) = (k < i)" for i
proof
assume "F i = 0" with k(2)[of "i - 1"]
have "¬ (i - 1 < k)" by (cases i, auto simp: F.simps)
thus "i > k" using k0 by auto
next
assume "i > k"
then obtain j l where i: "i = j + l" and "j = Suc k" and "l = i - Suc k" and Fj: "F j = 0" using k(1)
by auto
with F1 F2 k0 have j2: "j ≥ 2" by auto
show "F i = 0" unfolding i
proof (induct l)
case (Suc l)
thus ?case unfolding F.simps[of "j + Suc l"] using j2 by auto
qed (auto simp: Fj)
qed
show b: "b i ≠ 0" for i using Fb0 by blast
show "F 1 = map_poly to_fract G1" unfolding F.simps[of 1] by simp
show "F 2 = map_poly to_fract G2" unfolding F.simps[of 2] by simp
fix i
let ?mod = "pseudo_mod (F (i - 2)) (F (i - 1))"
assume i: "3 ≤ i" "i ≤ Suc k"
from Fk[of "i - 1"] i have "F (i - 1) ≠ 0" by auto
with i have Fi: "F i = (if ?mod = 0 then 0 else smult (inverse (b i)) ?mod)" unfolding F.simps[of i]
Let_def by simp
show "smult (b i) (F i) = ?mod"
proof (cases "?mod = 0")
case True
thus ?thesis unfolding Fi by simp
next
case False
with Fi have Fi: "F i = smult (inverse (b i)) ?mod" by simp
from arg_cong[OF this, of "smult (b i)"] b[of i] show ?thesis by simp
qed
qed
note s.h.simps[simp del]
show "subresultant_prs_locale2 F (λ i. degree (F i)) (λ i. degree (F i) - degree (F (Suc i)))
(λ i. lead_coeff (F i)) k b G1 G2"
proof
show "b 3 = (- 1) ^ (degree (F 1) - degree (F (Suc 1)) + 1)" unfolding b.simps numeral_2_eq_2 by simp
fix i
assume i: "4 ≤ i" "i ≤ Suc k"
with s.F0[of "i - 2"] have "F (i - 2) ≠ 0" by auto
hence bi: "b i = (- 1) ^ (degree (F (i - 2)) - degree (F (i - 1)) + 1) * lead_coeff (F (i - 2)) *
h (i - 2) ^ (degree (F (i - 2)) - degree (F (i - 1)))" unfolding b.simps
using i by auto
have "i < k ⟹ s.h i = h i" for i
proof (induct i)
case 0
thus ?case by (simp add: h.simps s.h.simps)
next
case (Suc i)
from Suc(2) s.F0[of "Suc i"] have "F (Suc i) ≠ 0" by auto
with Suc show ?case unfolding h.simps[of "Suc i"] s.h.simps[of "Suc i"] numeral_2_eq_2 by simp
qed
hence sh: "s.h (i - 2) = h (i - 2)" using i by simp
from i have *: "Suc (i - 2) = i - 1" by auto
show "b i = (- 1) ^ (degree (F (i - 2)) - degree (F (Suc (i - 2))) + 1) * lead_coeff (F (i - 2)) *
s.h (i - 2) ^ (degree (F (i - 2)) - degree (F (Suc (i - 2))))"
unfolding sh bi * ..
qed
qed
end
text ‹Now we obtain the soundness lemma outside the locale.›
context div_exp_sound
begin
lemma resultant_impl_main: assumes len: "length (coeffs G1) ≥ length (coeffs G2)"
shows "resultant_impl_main G1 G2 = resultant G1 G2"
proof (cases "G2 = 0")
case G2: False
from enter_subresultant_prs[OF len G2] obtain F n d f k b
where "subresultant_prs_locale2 F n d f k b G1 G2" by auto
interpret subresultant_prs_locale2 F n d f k b G1 G2 by fact
show ?thesis by (rule resultant_impl_main, standard)
next
case G2: True
show ?thesis unfolding G2
resultant_impl_main_def using resultant_const(2)[of G1 0] by simp
qed
theorem resultant_impl: "resultant_impl = resultant"
proof (intro ext)
fix f g :: "'a poly"
show "resultant_impl f g = resultant f g"
proof (cases "length (coeffs f) ≥ length (coeffs g)")
case True
thus ?thesis unfolding resultant_impl_def resultant_impl_main[OF True] by auto
next
case False
hence "length (coeffs g) ≥ length (coeffs f)" by auto
from resultant_impl_main[OF this]
show ?thesis unfolding resultant_impl_def resultant_swap[of f g] using False
by (auto simp: Let_def)
qed
qed
end
subsection ‹Code Equations›
text ‹In the following code-equations, we only compute the required values, e.g., $h_k$
is not required if $n_k > 0$, we compute $(-1)^{\ldots} * \ldots$ via a case-analysis,
and we perform special cases for $\delta_i = 1$, which is the most frequent case.›
context div_exp_param
begin
partial_function(tailrec) subresultant_prs_main_impl where
"subresultant_prs_main_impl f Gi_1 Gi ni_1 d1_1 hi_2 = (let
gi_1 = lead_coeff Gi_1;
ni = degree Gi;
hi_1 = (if d1_1 = 1 then gi_1 else div_exp gi_1 hi_2 d1_1);
d1 = ni_1 - ni;
pmod = pseudo_mod Gi_1 Gi
in (if pmod = 0 then f (Gi, (if d1 = 1 then lead_coeff Gi
else div_exp (lead_coeff Gi) hi_1 d1)) else
let
gi = lead_coeff Gi;
divisor = (-1) ^ (d1 + 1) * gi_1 * (hi_1 ^ d1) ;
Gi_p1 = sdiv_poly pmod divisor
in subresultant_prs_main_impl f Gi Gi_p1 ni d1 hi_1))"
definition subresultant_prs_impl where
"subresultant_prs_impl f G1 G2 = (let
pmod = pseudo_mod G1 G2;
n2 = degree G2;
delta_1 = (degree G1 - n2);
g2 = lead_coeff G2;
h2 = g2 ^ delta_1
in if pmod = 0 then f (G2,h2) else let
G3 = (- 1) ^ (delta_1 + 1) * pmod;
g3 = lead_coeff G3;
n3 = degree G3;
d2 = n2 - n3;
pmod = pseudo_mod G2 G3
in if pmod = 0 then f (G3, if d2 = 1 then g3 else div_exp g3 h2 d2)
else let divisor = (- 1) ^ (d2 + 1) * g2 * h2 ^ d2; G4 = sdiv_poly pmod divisor
in subresultant_prs_main_impl f G3 G4 n3 d2 h2
)"
end
context div_exp_sound
begin
lemma div_exp_1: "div_exp g h (Suc 0) = g"
using div_exp[of g "Suc 0" h] by simp
lemma subresultant_prs_impl: "subresultant_prs_impl f G1 G2 = f (subresultant_prs G1 G2)"
proof -
define h2 where "h2 = lead_coeff G2 ^ (degree G1 - degree G2)"
define G3 where "G3 = ((- 1) ^ (degree G1 - degree G2 + 1) * pseudo_mod G1 G2)"
define G4 where "G4 = sdiv_poly (pseudo_mod G2 G3)
((- 1) ^ (degree G2 - degree G3 + 1) * lead_coeff G2 *
h2 ^ (degree G2 - degree G3))"
define d2 where "d2 = degree G2 - degree G3"
have dl1: "(if d = 1 then (g :: 'a) else div_exp g h d) = div_exp g h d" for d g h
by (cases "d = 1", auto simp: div_exp_1)
show ?thesis
unfolding subresultant_prs_impl_def subresultant_prs_def Let_def
subresultant_prs_main.simps[of G2]
if_distrib[of f] dl1
proof (rule if_cong[OF refl refl if_cong[OF refl refl]], unfold h2_def[symmetric],
unfold G3_def[symmetric], unfold G4_def[symmetric], unfold d2_def[symmetric])
note simp = subresultant_prs_main_impl.simps[of f] subresultant_prs_main.simps
show "subresultant_prs_main_impl f G3 G4 (degree G3) d2 h2 =
f (subresultant_prs_main G3 G4 (div_exp (lead_coeff G3) h2 d2))"
proof (induct G4 arbitrary: G3 d2 h2 rule: wf_induct[OF wf_measure[of degree]])
case (1 G4 G3 d2 h2)
let ?M = "pseudo_mod G3 G4"
show ?case
proof (cases "?M = 0")
case True
thus ?thesis unfolding simp[of G3] Let_def dl1 by simp
next
case False
hence id: "(?M = 0) = False" by auto
let ?c = "((- 1) ^ (degree G3 - degree G4 + 1) * lead_coeff G3 *
(div_exp (lead_coeff G3) h2 d2) ^ (degree G3 - degree G4))"
let ?N = "sdiv_poly ?M ?c"
show ?thesis
proof (cases "G4 = 0")
case G4: False
have "degree ?N ≤ degree ?M" unfolding sdiv_poly_def by (rule degree_map_poly_le)
also have "… < degree G4" using pseudo_mod[OF G4, of G3] False by auto
finally show ?thesis unfolding simp[of G3] Let_def id if_False dl1
by (intro 1(1)[rule_format], auto)
next
case 0: True
with False have "G3 ≠ 0" by auto
show ?thesis unfolding 0 unfolding simp[of G3] Let_def unfolding dl1 simp[of 0] by simp
qed
qed
qed
qed
qed
definition
"resultant_impl_rec = subresultant_prs_main_impl (λ (Gk,hk). if degree Gk = 0 then hk else 0)"
definition
"resultant_impl_start = subresultant_prs_impl (λ (Gk,hk). if degree Gk = 0 then hk else 0)"
lemma resultant_impl_start_code:
"resultant_impl_start G1 G2 =
(let pmod = pseudo_mod G1 G2;
n2 = degree G2;
n1 = degree G1;
g2 = lead_coeff G2;
d1 = n1 - n2
in if pmod = 0 then if n2 = 0 then if d1 = 0 then 1 else if d1 = 1 then g2 else g2 ^ d1 else 0
else let
G3 = if even d1 then - pmod else pmod;
n3 = degree G3;
pmod = pseudo_mod G2 G3
in if pmod = 0
then if n3 = 0 then
let d2 = n2 - n3;
g3 = lead_coeff G3
in (if d2 = 1 then g3 else
div_exp g3 (if d1 = 1 then g2 else g2 ^ d1) d2) else 0
else let
h2 = (if d1 = 1 then g2 else g2 ^ d1);
d2 = n2 - n3;
divisor = (if d2 = 1 then g2 * h2 else if even d2 then - g2 * h2 ^ d2 else g2 * h2 ^ d2);
G4 = sdiv_poly pmod divisor
in resultant_impl_rec G3 G4 n3 d2 h2)"
proof -
obtain d1 where d1: "degree G1 - degree G2 = d1" by auto
have id1: "(if even d1 then - pmod else pmod) = (-1)^ (d1 + 1) * (pmod :: 'a poly)" for pmod by simp
have id3: "(if d2 = 1 then g2 * h2 else if even d2 then - g2 * h2 ^ d2 else g2 * h2 ^ d2) =
((- 1) ^ (d2 + 1) * g2 * h2 ^ d2)"
for d2 and g2 h2 :: 'a by auto
show ?thesis
unfolding resultant_impl_start_def subresultant_prs_impl_def resultant_impl_rec_def[symmetric] Let_def split
unfolding d1
unfolding id1
unfolding id3
by (rule if_cong[OF refl if_cong if_cong], auto simp: power2_eq_square)
qed
lemma resultant_impl_rec_code:
"resultant_impl_rec Gi_1 Gi ni_1 d1_1 hi_2 = (
let ni = degree Gi;
pmod = pseudo_mod Gi_1 Gi
in
if pmod = 0
then if ni = 0
then
let
d1 = ni_1 - ni;
gi = lead_coeff Gi
in if d1 = 1 then gi else
let gi_1 = lead_coeff Gi_1;
hi_1 = (if d1_1 = 1 then gi_1 else div_exp gi_1 hi_2 d1_1) in
div_exp gi hi_1 d1
else 0
else let
d1 = ni_1 - ni;
gi_1 = lead_coeff Gi_1;
hi_1 = (if d1_1 = 1 then gi_1 else div_exp gi_1 hi_2 d1_1);
divisor = if d1 = 1 then gi_1 * hi_1 else if even d1 then - gi_1 * hi_1 ^ d1 else gi_1 * hi_1 ^ d1;
Gi_p1 = sdiv_poly pmod divisor
in resultant_impl_rec Gi Gi_p1 ni d1 hi_1)"
unfolding resultant_impl_rec_def subresultant_prs_main_impl.simps[of _ Gi_1] split Let_def
unfolding resultant_impl_rec_def[symmetric]
by (rule if_cong[OF _ if_cong _], auto)
lemma resultant_impl_main_code: "resultant_impl_main G1 G2 =
(if G2 = 0 then if degree G1 = 0 then 1 else 0
else resultant_impl_start G1 G2)"
unfolding resultant_impl_main_def
resultant_impl_start_def subresultant_prs_impl by simp
lemma resultant_impl_code: "resultant_impl f g =
(if length (coeffs f) ≥ length (coeffs g) then resultant_impl_main f g
else let res = resultant_impl_main g f in
if even (degree f) ∨ even (degree g) then res else - res)"
unfolding resultant_impl_def resultant_impl_def ..
lemma resultant_code: "resultant = resultant_impl"
using resultant_impl by fastforce
lemmas resultant_code_lemmas =
resultant_impl_code
resultant_impl_main_code
resultant_impl_start_code
resultant_impl_rec_code
end
global_interpretation div_exp_Lazard: div_exp_sound "dichotomous_Lazard :: 'a :: factorial_ring_gcd ⇒ _"
defines
resultant_impl_Lazard = div_exp_Lazard.resultant_impl and
resultant_impl_main_Lazard = div_exp_Lazard.resultant_impl_main and
resultant_impl_start_Lazard = div_exp_Lazard.resultant_impl_start and
resultant_impl_rec_Lazard = div_exp_Lazard.resultant_impl_rec
by (rule dichotomous_Lazard)
declare div_exp_Lazard.resultant_code_lemmas[code]
text ‹As default use Lazard-implementation, which implements
resultants on factorial rings.›
declare div_exp_Lazard.resultant_code[code]
text ‹We also provide a second implementation without Lazard's
optimization, which works on integral domains.›
global_interpretation div_exp_basic: div_exp_sound basic_div_exp
defines
resultant_impl_basic = div_exp_basic.resultant_impl and
resultant_impl_main_basic = div_exp_basic.resultant_impl_main and
resultant_impl_start_basic = div_exp_basic.resultant_impl_start and
resultant_impl_rec_basic = div_exp_basic.resultant_impl_rec
by (rule basic_div_exp)
declare div_exp_basic.resultant_code_lemmas[code]
thm div_exp_basic.resultant_code
end