Theory Jordan_Normal_Form.Jordan_Normal_Form

(*  
    Author:      René Thiemann 
                 Akihisa Yamada
    License:     BSD
*)
section ‹Jordan Normal Form›

text ‹This theory defines Jordan normal forms (JNFs) in a sparse representation, i.e., 
  as block-diagonal matrices. We also provide a closed formula for powers of JNFs, 
  which allows to estimate the growth rates of JNFs.›

theory Jordan_Normal_Form
imports 
  Matrix
  Char_Poly
  Polynomial_Interpolation.Missing_Unsorted
begin

definition jordan_block :: "nat  'a :: {zero,one}  'a mat" where 
  "jordan_block n a = mat n n (λ (i,j). if i = j then a else if Suc i = j then 1 else 0)"

lemma jordan_block_index[simp]: "i < n  j < n  
  jordan_block n a $$ (i,j) = (if i = j then a else if Suc i = j then 1 else 0)"
  "dim_row (jordan_block n k) = n"
  "dim_col (jordan_block n k) = n"
  unfolding jordan_block_def by auto

lemma jordan_block_carrier[simp]: "jordan_block n k  carrier_mat n n" 
  unfolding carrier_mat_def by auto

lemma jordan_block_char_poly: "char_poly (jordan_block n a) = [: -a, 1:]^n"
  unfolding char_poly_defs by (subst det_upper_triangular[of _ n], auto simp: prod_list_diag_prod)

lemma jordan_block_pow_carrier[simp]:
  "jordan_block n a ^m r  carrier_mat n n" by auto
lemma jordan_block_pow_dim[simp]:
  "dim_row (jordan_block n a ^m r) = n" "dim_col (jordan_block n a ^m r) = n" by auto

lemma jordan_block_pow: "(jordan_block n (a :: 'a :: comm_ring_1)) ^m r = 
  mat n n (λ (i,j). if i  j then of_nat (r choose (j - i)) * a ^ (r + i - j) else 0)"
proof (induct r)
  case 0
  {
    fix i j :: nat
    assume "i  j" "i  j"
    hence "j - i > 0" by auto
    hence "0 choose (j - i) = 0" by simp
  } note [simp] = this
  show ?case
    by (simp, rule eq_matI, auto)
next
  case (Suc r)
  let ?jb = "jordan_block n a"
  let ?rij = "λ r i j. of_nat (r choose (j - i)) * a ^ (r + i - j)"
  let ?v = "λ i j. if i  j then of_nat (r choose (j - i)) * a ^ (r + i - j) else 0"
  have "?jb ^m Suc r = mat n n (λ (i,j). if i  j then ?rij r i j else 0) * ?jb" by (simp add: Suc)
  also have " = mat n n (λ (i,j). if i  j then ?rij (Suc r) i j else 0)"
  proof -
    {
      fix j
      assume j: "j < n"
      hence col: "col (jordan_block n a) j = vec n (λi. if i = j then a else if Suc i = j then 1 else 0)"
        unfolding jordan_block_def col_mat[OF j] by simp
      fix f
      have "vec n f  col (jordan_block n a) j = (f j * a + (if j = 0 then 0 else f (j - 1)))"
      proof -
        define p where "p = (λ i. vec n f $ i * col (jordan_block n a) j $ i)"
        have "vec n f  col (jordan_block n a) j = (i = 0 ..< n. p i)"
          unfolding scalar_prod_def p_def by simp
        also have " = p j + sum p ({0 ..< n} - {j})" using j
          by (subst sum.remove[of _ j], auto)
        also have "p j = f j * a" unfolding p_def col using j by auto
        also have "sum p ({0 ..< n} - {j}) = (if j = 0 then 0 else f (j - 1))"
        proof (cases j)
          case 0
          have "sum p ({0 ..< n} - {j}) = 0"
            by (rule sum.neutral, auto simp: p_def col 0)
          thus ?thesis using 0 by simp
        next
          case (Suc jj)
          with j have jj: "jj  {0 ..< n} - {j}" by auto
          have "sum p ({0 ..< n} - {j}) = p jj + sum p ({0 ..< n} - {j} - {jj})"
            by (subst sum.remove[OF _ jj], auto)
          also have "p jj = f (j - 1)" unfolding p_def col using jj
            by (auto simp: Suc)
          also have "sum p ({0 ..< n} - {j} - {jj}) = 0"
            by (rule sum.neutral, auto simp: p_def col, auto simp: Suc)
          finally show ?thesis unfolding Suc by simp
        qed
        finally show ?thesis .
      qed
    } note scalar_to_sum = this
    {
      fix i j
      assume i: "i < n" and ij: "i > j"
      hence j: "j < n" by auto
      have "vec n (?v i)  col (jordan_block n a) j = 0"
        unfolding scalar_to_sum[OF j] using ij i j by auto
    } note easy_case = this
    {
      fix i j
      assume j: "j < n" and ij: "i  j"
      hence i: "i < n" and id: " p q. (if i  j then p else q) = p" by auto
      have "vec n (?v i)  col (jordan_block n a) j =
        (of_nat (r choose (j - i)) * (a ^ (Suc (r + i - j)))) +
          (if j = 0 then 0
         else if i  j - 1 then of_nat (r choose (j - 1 - i)) * a ^ (r + i - (j - 1)) else 0)"
      unfolding scalar_to_sum[OF j]
      using ij by simp
      also have " = of_nat (Suc r choose (j - i)) * a ^ (Suc (r + i) - j)"
      proof (cases j)
        case (Suc jj)
        {
          assume "i  Suc jj" and "¬ i  jj"
          hence "i = Suc jj" by auto 
          hence "a * a ^ (r + i - Suc jj) = a ^ (r + i - jj)" by simp
        } 
        moreover
        {
          assume ijj: "i  jj"
          have "of_nat (r choose (Suc jj - i)) * (a * a ^ (r + i - Suc jj)) 
          + of_nat (r choose (jj - i)) * a ^ (r + i - jj) =
            of_nat (Suc r choose (Suc jj - i)) * a ^ (r + i - jj)"
          proof (cases "r + i < jj")
            case True
            hence gt: "jj - i > r" "Suc jj - i > r" "Suc jj - i > Suc r" by auto
            show ?thesis 
              unfolding binomial_eq_0[OF gt(1)] binomial_eq_0[OF gt(2)] binomial_eq_0[OF gt(3)]
              by simp
          next
            case False 
            hence ge: "r + i  jj" by simp
            show ?thesis
            proof (cases "jj = r + i")
              case True
              have gt: "r < Suc r" by simp
              show ?thesis unfolding True by (simp add: binomial_eq_0[OF gt])
            next
              case False
              with ge have lt: "jj < r + i" by auto
              hence "r + i - jj = Suc (r + i - Suc jj)" by simp 
              hence prod: "a * a ^ (r + i - Suc jj) = a ^ (r + i - jj)" by simp
              from ijj have id: "Suc jj - i = Suc (jj - i)" by simp
              have binom: "Suc r choose (Suc jj - i) = 
                r choose (Suc jj - i) + (r choose (jj - i))"
                unfolding id
                by (subst binomial_Suc_Suc, simp)
              show ?thesis unfolding prod binom  
                by (simp add: field_simps)
            qed
          qed
        }
        ultimately show ?thesis using ij unfolding Suc by auto
      qed auto
      finally have "vec n (?v i)  col (jordan_block n a) j 
        = of_nat (Suc r choose (j - i)) * a ^ (Suc (r + i) - j)" .
    } note main_case = this
    show ?thesis
      by (rule eq_matI, insert easy_case main_case, auto)
  qed
  finally show ?case by simp
qed

definition jordan_matrix :: "(nat × 'a :: {zero,one})list  'a mat" where
  "jordan_matrix n_as = diag_block_mat (map (λ (n,a). jordan_block n a) n_as)"

lemma jordan_matrix_dim[simp]: 
  "dim_row (jordan_matrix n_as) = sum_list (map fst n_as)"
  "dim_col (jordan_matrix n_as) = sum_list (map fst n_as)"
  unfolding jordan_matrix_def
  by (subst dim_diag_block_mat, auto, (induct n_as, auto simp: Let_def)+)

lemma jordan_matrix_carrier[simp]: 
  "jordan_matrix n_as  carrier_mat (sum_list (map fst n_as)) (sum_list (map fst n_as))"
  unfolding carrier_mat_def by auto

lemma jordan_matrix_upper_triangular: "i < sum_list (map fst n_as)
   j < i  jordan_matrix n_as $$ (i,j) = 0"
  unfolding jordan_matrix_def
  by (rule diag_block_upper_triangular, auto simp: jordan_matrix_def[symmetric])

lemma jordan_matrix_pow: "(jordan_matrix n_as) ^m r = 
  diag_block_mat (map (λ (n,a). (jordan_block n a) ^m r) n_as)"
  unfolding jordan_matrix_def
  by (subst diag_block_pow_mat, force, rule arg_cong[of _ _ diag_block_mat], auto)

lemma jordan_matrix_char_poly: 
  "char_poly (jordan_matrix n_as) = ((n, a)n_as. [:- a, 1:] ^ n)"
proof -
  let ?n = "sum_list (map fst n_as)"
  have "diag_mat
     ([:0, 1:] m 1m (sum_list (map fst n_as)) + map_mat (λa. [:- a:]) (jordan_matrix n_as)) =
    concat (map (λ(n, a). replicate n [:- a, 1:]) n_as)" unfolding jordan_matrix_def
  proof (induct n_as)
    case (Cons na n_as)
    obtain n a where na: "na = (n,a)" by force
    let ?n2 = "sum_list (map fst n_as)"
    note fbo = four_block_one_mat
    note mz = zero_carrier_mat
    note mo = one_carrier_mat
    have mA: " A. A  carrier_mat (dim_row A) (dim_col A)" unfolding carrier_mat_def by auto
    let ?Bs = "map (λ(x, y). jordan_block x y) n_as"
    let ?B = "diag_block_mat ?Bs"
    from jordan_matrix_dim[of n_as, unfolded jordan_matrix_def]
    have dimB: "dim_row ?B = ?n2" "dim_col ?B = ?n2" by auto
    hence B: "?B  carrier_mat ?n2 ?n2" unfolding carrier_mat_def by simp
    show ?case unfolding na fbo
    apply (simp add: Let_def fbo[symmetric] del: fbo)
    apply (subst smult_four_block_mat[OF mo mz mz mo])
    apply (subst map_four_block_mat[OF jordan_block_carrier mz mz mA])
    apply (subst add_four_block_mat[of _ n n _ ?n2 _ ?n2], auto simp: dimB B)
    apply (subst diag_four_block_mat[of _ n _ ?n2], auto simp: dimB B)
    apply (subst Cons, auto simp: jordan_block_def diag_mat_def, 
      intro nth_equalityI, auto)
    done
  qed (force simp: diag_mat_def)
  also have "prod_list ... = ((n, a)n_as. [:- a, 1:] ^ n)"
    by (induct n_as, auto)
  finally
  show ?thesis unfolding char_poly_defs
    by (subst det_upper_triangular[of _ ?n], auto simp: jordan_matrix_upper_triangular)
qed

definition jordan_nf :: "'a :: semiring_1 mat  (nat × 'a)list  bool" where
  "jordan_nf A n_as  (0  fst ` set n_as  similar_mat A (jordan_matrix n_as))"

lemma jordan_nf_powE: assumes A: "A  carrier_mat n n" and jnf: "jordan_nf A n_as" 
  obtains P Q where "P  carrier_mat n n" "Q  carrier_mat n n" and 
  "char_poly A = ((na, a)n_as. [:- a, 1:] ^ na)"
  " k. A ^m k = P * (jordan_matrix n_as)^m k * Q"
proof -
  from A have dim: "dim_row A = n" by auto
  assume obt: "P Q. P  carrier_mat n n  Q  carrier_mat n n  
    char_poly A = ((na, a)n_as. [:- a, 1:] ^ na)  
    (k. A ^m k = P * jordan_matrix n_as ^m k * Q)  thesis"
  from jnf[unfolded jordan_nf_def] obtain P Q where
    simw: "similar_mat_wit A (jordan_matrix n_as) P Q"
    and sim: "similar_mat A (jordan_matrix n_as)" unfolding similar_mat_def by blast
  show thesis
  proof (rule obt)
    show " k. A ^m k = P * jordan_matrix n_as ^m k * Q"
      by (rule similar_mat_wit_pow_id[OF simw])
    show "char_poly A = ((na, a)n_as. [:- a, 1:] ^ na)"
      unfolding char_poly_similar[OF sim] jordan_matrix_char_poly ..    
  qed (insert simw[unfolded similar_mat_wit_def Let_def dim], auto)
qed

lemma choose_poly_bound: assumes "i  d"
  shows "r choose i  max 1 (r^d)"
proof (cases "i  r")
  case False
  hence "r choose i = 0" by simp
  thus ?thesis by arith
next
  case True
  show ?thesis
  proof (cases r)
    case (Suc rr)
    from binomial_le_pow[OF True] have "r choose i  r ^ i" by simp
    also have "  r^d" using power_increasing[OF i  d, of r] Suc by auto
    finally show ?thesis by simp
  qed (insert True, simp)
qed  

context
  fixes b :: "'a :: archimedean_field"
  assumes b: "0 < b" "b < 1"
begin
      
lemma poly_exp_constant_bound: " p.  x. c * b ^ x * of_nat x ^ deg  p" 
proof (cases "c  0")
  case True
  show ?thesis
    by (rule exI[of _ 0], intro allI, 
    rule mult_nonpos_nonneg[OF mult_nonpos_nonneg[OF True]], insert b, auto)
next
  case False
  hence c: "c  0" by simp
  from poly_exp_bound[OF b, of deg] obtain p where " x. b ^ x * of_nat x ^ deg  p" by auto
  from mult_left_mono[OF this c]
  show ?thesis by (intro exI[of _ "c * p"], auto simp: ac_simps)
qed

lemma poly_exp_max_constant_bound: " p.  x. c * b ^ x * max 1 (of_nat x ^ deg)  p" 
proof -
  from poly_exp_constant_bound[of c deg] obtain p where
    p: " x. c * b ^ x * of_nat x ^ deg  p" by auto
  show ?thesis
  proof (rule exI[of _ "max p c"], intro allI)
    fix x
    let ?exp = "of_nat x ^ deg :: 'a"
    show "c * b ^ x * max 1 ?exp  max p c"
    proof (cases "x = 0")
      case False
      hence "?exp  of_nat 0" by simp
      hence "?exp  1" by (metis less_one not_less of_nat_1 of_nat_less_iff of_nat_power)
      hence "max 1 ?exp = ?exp" by simp
      thus ?thesis using p[of x] by simp
    qed (cases deg, auto)
  qed
qed
end

context
  fixes a :: "'a :: real_normed_field"
begin
lemma jordan_block_bound: 
  assumes i: "i < n" and j: "j < n"
  shows "norm ((jordan_block n a ^m k) $$ (i,j)) 
     norm a ^ (k + i - j) * max 1 (of_nat k ^ (n - 1))"
    (is "?lhs  ?rhs")
proof -
  have id: "(jordan_block n a ^m k) $$ (i,j) = (if i  j then of_nat (k choose (j - i)) * a ^ (k + i - j) else 0)"
    unfolding jordan_block_pow using i j by auto
  from i j have diff: "j - i  n - 1" by auto
  show ?thesis
  proof (cases "i  j")
    case False
    thus ?thesis unfolding id by simp
  next
    case True
    hence "?lhs = norm (of_nat (k choose (j - i)) * a ^ (k + i - j))" unfolding id by simp
    also have "  norm (of_nat (k choose (j - i)) :: 'a) * norm (a ^ (k + i - j))"
      by (rule norm_mult_ineq)
    also have "  (max 1 (of_nat k ^ (n - 1))) * norm a ^ (k + i - j)"
    proof (rule mult_mono[OF _ norm_power_ineq _ norm_ge_zero])
      have "k choose (j - i)  max 1 (k ^ (n - 1))" 
        by (rule choose_poly_bound[OF diff])
      hence "norm (of_nat (k choose (j - i)) :: 'a)  of_nat (max 1 (k ^ (n - 1)))"
        unfolding norm_of_nat of_nat_le_iff .
      also have " = max 1 (of_nat k ^ (n - 1))" by (metis max_def of_nat_1 of_nat_le_iff of_nat_power)
      finally show "norm (of_nat (k choose (j - i)) :: 'a)  max 1 (real_of_nat k ^ (n - 1))" .
    qed simp
    also have " = ?rhs" by simp
    finally show ?thesis .
  qed
qed

lemma jordan_block_poly_bound: 
  assumes i: "i < n" and j: "j < n" and a: "norm a = 1"
  shows "norm ((jordan_block n a ^m k) $$ (i,j))  max 1 (of_nat k ^ (n - 1))"
    (is "?lhs  ?rhs")
proof -
  from jordan_block_bound[OF i j, of k, unfolded a]
  show ?thesis by simp
qed


theorem jordan_block_constant_bound: assumes a: "norm a < 1" 
  shows " p.  i j k. i < n  j < n  norm ((jordan_block n a ^m k) $$ (i,j))  p"
proof (cases "a = 0") 
  case True
  show ?thesis
  proof (rule exI[of _ 1], intro allI impI)
    fix i j k
    assume *: "i < n" "j < n"
    {
      assume ij: "i  j"
      have "norm ((of_nat (k choose (j - i)) :: 'a) * 0 ^ (k + i - j))  1" (is "norm ?lhs  1")
      proof (cases "k + i > j")
        case True
        hence "?lhs = 0" by simp
        also have "norm ()  1" by simp
        finally show ?thesis .
      next
        case False
        hence id: "?lhs = (of_nat (k choose (j - i)) :: 'a)" and j: "j - i  k" by auto
        from j have "k choose (j - i) = 0  k choose (j - i) = 1" by (simp add: nat_less_le)
        thus "norm ?lhs  1"
        proof
          assume k: "k choose (j - i) = 0"
          show ?thesis unfolding id k by simp
        next
          assume k: "k choose (j - i) = 1"
          show ?thesis unfolding id unfolding k by simp
        qed
      qed
    }    
    thus "norm ((jordan_block n a ^m k) $$ (i,j))  1" unfolding True
      unfolding jordan_block_pow using * by auto
  qed
next
  case False
  hence na: "norm a > 0" by auto
  define c where "c = inverse (norm a ^ n)"
  define deg where "deg = n - 1"
  have c: "c > 0" unfolding c_def using na by auto
  define b where "b = norm a"
  from a na have "0 < b" "b < 1" unfolding b_def by auto
  from poly_exp_max_constant_bound[OF this, of c deg]
  obtain p where " k. c * b ^ k * max 1 (of_nat k ^ deg)  p" by auto
  show ?thesis
  proof (intro exI[of _ p], intro allI impI)
    fix i j k
    assume ij: "i < n" "j < n"
    from jordan_block_bound[OF this]
    have "norm ((jordan_block n a ^m k) $$ (i, j))
       norm a ^ (k + i - j) * max 1 (real_of_nat k ^ (n - 1))" .
    also have "  c * norm a ^ k * max 1 (real_of_nat k ^ (n - 1))"
    proof (rule mult_right_mono)
      from ij have "i - j  n" by auto
      show "norm a ^ (k + i - j)  c * norm a ^ k"
      proof (rule mult_left_le_imp_le)
        show "0 < norm a ^ n" using na by auto
        let ?lhs = "norm a ^ n * norm a ^ (k + i - j)"
        let ?rhs = "norm a ^ n * (c * norm a ^ k)"
        from ij have ge: "n + (k + i - j)  k" by arith
        have "?lhs = norm a ^ (n + (k + i - j))" by (simp add: power_add)
        also have "  norm a ^ k" using ge a na using less_imp_le power_decreasing by blast
        also have " = ?rhs" unfolding c_def using na by simp
        finally show "?lhs  ?rhs" .
      qed
    qed simp
    also have " = c * b ^ k * max 1 (real_of_nat k ^ deg)" unfolding b_def deg_def ..
    also have "  p" by fact
    finally show "norm ((jordan_block n a ^m k) $$ (i, j))  p" .
  qed
qed

definition norm_bound :: "'a mat  real  bool" where
  "norm_bound A b   i j. i < dim_row A  j < dim_col A  norm (A $$ (i,j))  b"

lemma norm_boundI[intro]:
  assumes " i j. i < dim_row A  j < dim_col A  norm (A $$ (i,j))  b"
  shows "norm_bound A b"
  unfolding norm_bound_def using assms by blast

lemma  jordan_block_constant_bound2:
"p. norm (a :: 'a :: real_normed_field) < 1 
    (i j k. i < n  j < n  norm ((jordan_block n a ^m k) $$ (i, j))  p)"
using jordan_block_constant_bound by auto

lemma jordan_matrix_poly_bound2:
  fixes n_as :: "(nat × 'a) list"
  assumes n_as: " n a. (n,a)  set n_as  n > 0  norm a  1"
  and N: " n a. (n,a)  set n_as  norm a = 1  n  N"
  shows "c1. k. e  elements_mat (jordan_matrix n_as ^m k).
    norm e  c1 + of_nat k ^ (N - 1)"
proof -
  from jordan_matrix_carrier[of n_as] obtain d where
    jm: "jordan_matrix n_as  carrier_mat d d" by blast
  define f where "f = (λn (a::'a) i j k. norm ((jordan_block n a ^m k) $$ (i,j)))"
  let ?g = "λk c1. c1 + of_nat k ^ (N-1)"
  let ?P = "λn (a::'a) i j k c1. f n a i j k  ?g k c1"
  define Q where "Q = (λn (a::'a) k c1. i j. i<n  j<n  ?P n a i j k c1)"
  have " c c' k n a i j. c  c'  ?P n a i j k c  ?P n a i j k c'" by auto  
  hence Q_mono: "n a c c'. c  c'  k. Q n a k c  k. Q n a k c'"
    unfolding Q_def by arith
  { fix n a assume na: "(n,a)  set n_as"
    obtain c where c: "norm a < 1  (i j k. i < n  j < n  f n a i j k  c)"
      apply (rule exE[OF jordan_block_constant_bound2])
      unfolding f_def using Jordan_Normal_Form.jordan_block_constant_bound2
      by metis
    define c1 where "c1 = max 1 c"
    then have "c1  1" "c1  c" by auto
    have "c1. k i j. i < n  j < n  ?P n a i j k c1"
    proof rule+
      fix i j k assume "i < n" "j < n"
      then have "0<n" by auto
      let ?jbs = "map (λ(n,a). jordan_block n a) n_as"
      have sq_jbs: "Ball (set ?jbs) square_mat" by auto
      have "jordan_matrix n_as ^m k = diag_block_mat (map (λA. A ^m k) ?jbs)"
        unfolding jordan_matrix_def using diag_block_pow_mat[OF sq_jbs] by auto
      show "?P n a i j k c1"
      proof (cases "norm a = 1")
        case True {
          have nN:"n-1  N-1" using N[OF na] True by auto
          have "f n a i j k  max 1 (of_nat k ^ (n-1))"
            using Jordan_Normal_Form.jordan_block_poly_bound True i<n j<n
            unfolding f_def by auto
          also have "...  max 1 (of_nat k ^ (N-1))"
            proof (cases "k=0")
              case False then show ?thesis
                by (subst max.mono[OF _ power_increasing[OF nN]], auto)
            qed (simp add: power_eq_if)
          also have "...  max c1 (of_nat k ^ (N-1))" using c11 by auto
          also have "...  c1 + (of_nat k ^ (N-1))" using c11 by auto
          finally show ?thesis by simp
        } next
        case False {
          then have na1: "norm a < 1" using n_as[OF na] 0<n by auto
          hence "f n a i j k  c" using c i<n j<n by auto
          also have "...  c1" using cc1.
          also have "...  c1 + of_nat k ^ (N-1)" by auto
          finally show ?thesis by auto
        }
      qed
    qed
  }
  hence "na. c1. na  set n_as  (k. Q (fst na) (snd na) k c1)"
    unfolding Q_def by auto
  from choice[OF this] obtain c'
    where c': " na k. na  set n_as  Q (fst na) (snd na) k (c' na)" by blast
  define c where "c = max 0 (Max (set (map c' n_as)))"
  { fix n a assume na: "(n,a)  set n_as"
    then have Q: " k. Q n a k (c' (n,a))" using c'[OF na] by auto
    from na have "c' (n,a)  set (map c' n_as)" by auto
    from Max_ge[OF _ this] have "c' (n,a)  c" unfolding c_def by auto
    from Q_mono[OF this Q] have " k. Q n a k c" by blast
  }
  hence Q: "k n a. (n,a)  set n_as  Q n a k c" by auto
  have c0: "c  0" unfolding c_def by simp
  { fix k n a e
    assume na:"(n,a)  set n_as"
    let ?jbk = "jordan_block n a ^m k"
    assume "e  elements_mat ?jbk"
    from elements_matD[OF this] obtain i j
      where "i < n" "j < n" and [simp]: "e = ?jbk $$ (i,j)"
      by (simp only:pow_mat_dim_square[OF jordan_block_carrier],auto)
    hence "norm e  ?g k c" using Q[OF na] unfolding Q_def f_def by simp
  }
  hence norm_jordan:
    "k. (n,a)  set n_as. e  elements_mat (jordan_block n a ^m k).
     norm e  ?g k c" by auto
  { fix k
    let ?jmk = "jordan_matrix n_as ^m k"
    have "dim_row ?jmk = d" "dim_col ?jmk = d"
      using jm by (simp only:pow_mat_dim_square[OF jm])+
    let ?As = "(map (λ(n,a). jordan_block n a ^m k) n_as)"
    have "e. e  elements_mat ?jmk  norm e  ?g k c"
    proof -
      fix e assume e:"e  elements_mat ?jmk"
      obtain i j where ij: "i < d" "j < d" and "e = ?jmk $$ (i,j)"
        using elements_matD[OF e] by (simp only:pow_mat_dim_square[OF jm],auto)
      have "?jmk = diag_block_mat ?As"
        using jordan_matrix_pow[of n_as k] by auto
      hence "elements_mat ?jmk  {0}   (set (map elements_mat ?As))"
        using elements_diag_block_mat[of ?As] by auto
      hence e_mem: "e  {0}   (set (map elements_mat ?As))"
        using e by blast
      show "norm e  ?g k c"
      proof (cases "e = 0")
        case False
          then have "e   (set (map elements_mat ?As))" using e_mem by auto
          then obtain n a
            where "e  elements_mat (jordan_block n a ^m k)"
            and na: "(n,a)  set n_as" by force
          thus ?thesis using norm_jordan na by force
      qed (insert c0, auto)
    qed
  }
  thus ?thesis by auto
qed

lemma norm_bound_bridge:
  "e  elements_mat A. norm e  b  norm_bound A b"
  unfolding norm_bound_def by force

lemma norm_bound_mult: assumes A1: "A1  carrier_mat nr n"
  and A2: "A2  carrier_mat n nc"
  and b1: "norm_bound A1 b1"
  and b2: "norm_bound A2 b2"
  shows "norm_bound (A1 * A2) (b1 * b2 * of_nat n)"
proof 
  let ?A = "A1 * A2"
  let ?n = "of_nat n"
  fix i j
  assume i: "i < dim_row ?A" and j: "j < dim_col ?A"
  define v1 where "v1 = (λ k. row A1 i $ k)"
  define v2 where "v2 = (λ k. col A2 j $ k)"
  from assms(1-2) have dim: "dim_row A1 = nr" "dim_col A2 = nc" "dim_col A1 = n" "dim_row A2 = n" by auto
  {
    fix k
    assume k: "k < n"
    have n: "norm (v1 k)  b1" "norm (v2 k)  b2" 
      using i j k dim v1_def v2_def
      b1[unfolded norm_bound_def, rule_format, of i k] 
      b2[unfolded norm_bound_def, rule_format, of k j] by auto
    have "norm (v1 k * v2 k)  norm (v1 k) * norm (v2 k)" by (rule norm_mult_ineq)
    also have "  b1 * b2" by (rule mult_mono'[OF n], auto)
    finally have "norm (v1 k * v2 k)  b1 * b2" .
  } note bound = this
  have "?A $$ (i,j) = row A1 i  col A2 j" using dim i j by simp
  also have " = ( k = 0 ..< n. v1 k * v2 k)" unfolding scalar_prod_def 
    using dim i j v1_def v2_def by simp
  also have "norm ()  ( k = 0 ..< n. b1 * b2)" 
    by (rule sum_norm_le, insert bound, simp)
  also have " = b1 * b2 * ?n" by simp
  finally show "norm (?A $$ (i,j))  b1 * b2 * ?n" .
qed

lemma norm_bound_max: "norm_bound A (Max {norm (A $$ (i,j)) | i j. i < dim_row A  j < dim_col A})" 
  (is "norm_bound A (Max ?norms)")
proof 
  fix i j
  have fin: "finite ?norms" by (simp add: finite_image_set2)
  assume "i < dim_row A" and "j < dim_col A"     
  hence "norm (A $$ (i,j))  ?norms" by auto
  from Max_ge[OF fin this] show "norm (A $$ (i,j))  Max ?norms" .
qed

lemma jordan_matrix_poly_bound: fixes n_as :: "(nat × 'a)list"
  assumes n_as: " n a. (n,a)  set n_as  n > 0  norm a  1"
  and N: " n a. (n,a)  set n_as  norm a = 1  n  N"
  shows " c1.  k. norm_bound (jordan_matrix n_as ^m k) (c1 + of_nat k ^ (N - 1))" 
  using jordan_matrix_poly_bound2 norm_bound_bridge N n_as
  by metis

lemma jordan_nf_matrix_poly_bound: fixes n_as :: "(nat × 'a)list"
  assumes A: "A  carrier_mat n n"
  and n_as: " n a. (n,a)  set n_as  n > 0  norm a  1"
  and N: " n a. (n,a)  set n_as  norm a = 1  n  N"
  and jnf: "jordan_nf A n_as"
  shows " c1 c2.  k. norm_bound (A ^m k) (c1 + c2 * of_nat k ^ (N - 1))"
proof -
  let ?cp2 = "(n, a)n_as. [:- a, 1:] ^ n"
  let ?J = "jordan_matrix n_as"
  from jnf[unfolded jordan_nf_def]
  have sim: "similar_mat A ?J" by auto
  then obtain P Q where sim_wit: "similar_mat_wit A ?J P Q" unfolding similar_mat_def by auto
  from similar_mat_wit_pow_id[OF this] have pow: " k. A ^m k = P * ?J ^m k * Q" .
  from sim_wit[unfolded similar_mat_wit_def Let_def] A 
  have J: "?J  carrier_mat n n" and P: "P  carrier_mat n n" and Q: "Q  carrier_mat n n"
    unfolding carrier_mat_def by force+
  have "c1.  k. norm_bound (?J ^m k) (c1 + of_nat k ^ (N - 1))"
    by (rule jordan_matrix_poly_bound[OF n_as N])
  then obtain c1 where 
    bound_pow: " k. norm_bound ((?J ^m k)) (c1 + of_nat k ^ (N - 1))" by blast
  obtain bP where bP: "norm_bound P bP" using norm_bound_max[of P] by auto
  obtain bQ where bQ: "norm_bound Q bQ" using norm_bound_max[of Q] by auto
  let ?n = "of_nat n :: real"
  let ?c2 = "bP * ?n * bQ * ?n"
  let ?c1 = "?c2 * c1"
  {
    fix k
    have Jk: "?J ^m k  carrier_mat n n" using J by simp
    from norm_bound_mult[OF mult_carrier_mat[OF P Jk] Q 
      norm_bound_mult[OF P Jk bP bound_pow] bQ, folded pow] 
    have "norm_bound (A ^m k) (?c1 + ?c2 * of_nat k ^ (N - 1))"  (is "norm_bound _ ?exp") 
      by (simp add: field_simps)
  } note main = this
  show ?thesis 
    by (intro exI allI, rule main)
qed
end

context 
  fixes f_ty :: "'a :: field itself"
begin
lemma char_matrix_jordan_block: "char_matrix (jordan_block n a) b = (jordan_block n (a - b))"
  unfolding char_matrix_def jordan_block_def by auto

lemma diag_jordan_block_pow: "diag_mat (jordan_block n (a :: 'a) ^m k) = replicate n (a ^ k)"
  unfolding diag_mat_def jordan_block_pow
  by (intro nth_equalityI, auto)

lemma jordan_block_zero_pow: "(jordan_block n (0 :: 'a)) ^m k = 
  (mat n n (λ (i,j). if j  i  j - i = k then 1 else 0))"
proof -
  {
    fix i j
    assume  *: "j - i  k"
    have "of_nat (k choose (j - i)) * 0 ^ (k + i - j) = (0 :: 'a)"
    proof (cases "k + i - j > 0")
      case True thus ?thesis by (cases "k + i - j", auto)
    next
      case False
      with * have "j - i > k" by auto
      thus ?thesis by (simp add: binomial_eq_0)
    qed
  }
  thus ?thesis unfolding jordan_block_pow by (intro eq_matI, auto)
qed
end

lemma jordan_matrix_concat_diag_block_mat: "jordan_matrix (concat jbs) = diag_block_mat (map jordan_matrix jbs)"
  unfolding jordan_matrix_def[abs_def]
  by (induct jbs, auto simp: diag_block_mat_append Let_def)

lemma jordan_nf_diag_block_mat: assumes Ms: " A jbs. (A,jbs)  set Ms  jordan_nf A jbs"
  shows "jordan_nf (diag_block_mat (map fst Ms)) (concat (map snd Ms))"
proof -
  let ?Ms = "map (λ (A, jbs). (A, jordan_matrix jbs)) Ms"
  have id: "map fst ?Ms = map fst Ms" by auto
  have id2: "map snd ?Ms = map jordan_matrix (map snd Ms)" by auto
  {
    fix A B
    assume "(A,B)  set ?Ms"
    then obtain jbs where mem: "(A,jbs)  set Ms" and B: "B = jordan_matrix jbs" by auto
    from Ms[OF mem] have "similar_mat A B" unfolding B jordan_nf_def by auto
  }
  from similar_diag_mat_block_mat[of ?Ms, OF this, unfolded id id2] Ms
  show ?thesis
    unfolding jordan_nf_def jordan_matrix_concat_diag_block_mat by force
qed  


lemma jordan_nf_char_poly: assumes "jordan_nf A n_as"
  shows "char_poly A = ( (n,a)  n_as. [:- a, 1:] ^ n)"
  unfolding jordan_matrix_char_poly[symmetric]
  by (rule char_poly_similar, insert assms[unfolded jordan_nf_def], auto)

lemma jordan_nf_block_size_order_bound: assumes jnf: "jordan_nf A n_as"
  and mem: "(n,a)  set n_as"
  shows "n  order a (char_poly A)"
proof -
  from jnf[unfolded jordan_nf_def]
  have "similar_mat A (jordan_matrix n_as)" by auto
  from similar_matD[OF this] obtain m where "A  carrier_mat m m" by auto
  from degree_monic_char_poly[OF this] have A: "char_poly A  0" by auto
  from mem obtain as bs where nas: "n_as = as @ (n,a) # bs" 
    by (meson split_list)
  from jordan_nf_char_poly[OF jnf] 
  have cA: "char_poly A = ((n, a)n_as. [:- a, 1:] ^ n)" .
  also have " = [: -a, 1:] ^ n * ((n, a) as @ bs. [:- a, 1:] ^ n)" unfolding nas by auto
  also have "[: -a,1 :] ^ n dvd " unfolding dvd_def by blast
  finally have "[: -a,1 :] ^ n dvd char_poly A" by auto
  from order_max[OF this A] show ?thesis .
qed

lemma similar_mat_jordan_block_smult: fixes A :: "'a :: field mat" 
  assumes "similar_mat A (jordan_block n a)" 
   and k: "k  0" 
  shows "similar_mat (k m A) (jordan_block n (k * a))" 
proof -
  let ?J = "jordan_block n a" 
  let ?Jk = "jordan_block n (k * a)" 
  let ?kJ = "k m jordan_block n a" 
  from k have inv: "k ^ i  0" for i by auto
  let ?A = "mat_diag n (λ i. k^i)" 
  let ?B = "mat_diag n (λ i. inverse (k^i))"
  have "similar_mat_wit ?Jk ?kJ ?A ?B" 
  proof (rule similar_mat_witI)
    show "jordan_block n (k * a) = ?A * ?kJ * ?B"
      by (subst mat_diag_mult_left[of _ _ n], force, subst mat_diag_mult_right[of _ n],
       insert k inv, auto simp: jordan_block_def field_simps intro!: eq_matI)
  qed (auto simp: inv field_simps k)
  hence kJ: "similar_mat ?Jk ?kJ" 
    unfolding similar_mat_def by auto
  have "similar_mat A ?J" by fact
  hence "similar_mat (k m A) (k m ?J)" by (rule similar_mat_smult)
  with kJ show ?thesis
    using similar_mat_sym similar_mat_trans by blast
qed


lemma jordan_matrix_Cons:  "jordan_matrix (Cons (n,a) n_as) = four_block_mat 
  (jordan_block n a)                 (0m n (sum_list (map fst n_as))) 
  (0m (sum_list (map fst n_as)) n)   (jordan_matrix n_as)" 
  unfolding jordan_matrix_def by (simp, simp add: jordan_matrix_def[symmetric])

lemma similar_mat_jordan_matrix_smult:  fixes n_as :: "(nat × 'a :: field) list"
  assumes k: "k  0" 
  shows "similar_mat (k m jordan_matrix n_as) (jordan_matrix (map (λ (n,a). (n, k * a)) n_as))" 
proof (induct n_as)
  case Nil
  show ?case by (auto simp: jordan_matrix_def intro!: similar_mat_refl)
next
  case (Cons na n_as)
  obtain n a where na: "na = (n,a)" by force
  let ?l = "map (λ (n,a). (n, k * a))" 
  let ?n = "sum_list (map fst n_as)" 
  have "k m jordan_matrix (Cons na n_as) = k m four_block_mat 
     (jordan_block n a) (0m n ?n)
     (0m ?n n) (jordan_matrix n_as)" (is "?M = _ m four_block_mat ?A ?B ?C ?D")
    by (simp add: na jordan_matrix_Cons)
  also have " = four_block_mat (k m ?A) ?B ?C (k m ?D)" 
    by (subst smult_four_block_mat, auto)
  finally have jm: "?M = four_block_mat (k m ?A) ?B ?C (k m ?D)" .
  have [simp]: "fst (case x of (n :: nat, a)  (n, k * a)) = fst x" for x by (cases x, auto)
  have jmk: "jordan_matrix (?l (Cons na n_as)) = four_block_mat
     (jordan_block n (k * a)) ?B
     ?C (jordan_matrix (?l n_as))" (is "?kM = four_block_mat ?kA _ _ ?kD")
    by (simp add: na jordan_matrix_Cons o_def)
  show ?case unfolding jmk jm
    by (rule similar_mat_four_block_0_0[OF similar_mat_jordan_block_smult[OF _ k] Cons],
      auto intro!: similar_mat_refl)
qed

lemma jordan_nf_smult: fixes k :: "'a :: field" 
  assumes jn: "jordan_nf A n_as" 
  and k: "k  0" 
  shows "jordan_nf (k m A) (map (λ (n,a). (n, k * a)) n_as)" 
proof -
  let ?l = "map (λ (n,a). (n, k * a))" 
  from jn[unfolded jordan_nf_def] have sim: "similar_mat A (jordan_matrix n_as)" by auto
  from similar_mat_smult[OF this, of k] similar_mat_jordan_matrix_smult[OF k, of n_as]
  have "similar_mat (k m A) (jordan_matrix (map (λ(n, a). (n, k * a)) n_as))" 
    using similar_mat_trans by blast
  with jn show ?thesis unfolding jordan_nf_def by force
qed

lemma jordan_nf_order: assumes "jordan_nf A n_as" 
  shows "order a (char_poly A)  = sum_list (map fst (filter (λ na. snd na = a) n_as))" 
proof - 
  let ?p = "λ n_as. ((n, a)n_as. [:- a, 1:] ^ n)" 
  let ?s = "λ n_as. sum_list (map fst (filter (λ na. snd na = a) n_as))" 
  from jordan_nf_char_poly[OF assms]
  have "order a (char_poly A) = order a (?p n_as)" by simp
  also have " = ?s n_as" 
  proof (induct n_as)
    case (Cons nb n_as)
    obtain n b where nb: "nb = (n,b)" by force
    have "order a (?p (nb # n_as)) = order a ([: -b, 1:] ^ n * ?p n_as)" unfolding nb by simp
    also have " = order a ([: -b, 1:] ^ n) + order a (?p n_as)" 
      by (rule order_mult, auto simp: prod_list_zero_iff)
    also have " = (if a = b then n else 0) + ?s n_as" unfolding Cons order_linear_power by simp
    also have " = ?s (nb # n_as)" unfolding nb by auto
    finally show ?case .
  qed simp
  finally show ?thesis .
qed

subsection ‹Application for Complexity›

lemma factored_char_poly_norm_bound: assumes A: "A  carrier_mat n n"
  and linear_factors: "char_poly A = ( (a :: 'a :: real_normed_field)  as. [:- a, 1:])"
  and jnf_exists: " n_as. jordan_nf A n_as" 
  and le_1: " a. a  set as  norm a  1"
  and le_N: " a. a  set as  norm a = 1  length (filter ((=) a) as)  N"
  shows " c1 c2.  k. norm_bound (A ^m k) (c1 + c2 * of_nat k ^ (N - 1))"
proof -
  from jnf_exists obtain n_as 
    where jnf: "jordan_nf A n_as" by auto
  let ?cp1 = "( a  as. [:- a, 1:])"
  let ?cp2 = "(n, a)n_as. [:- a, 1:] ^ n"
  let ?J = "jordan_matrix n_as"
  from jnf[unfolded jordan_nf_def]
  have sim: "similar_mat A ?J" by auto
  from char_poly_similar[OF sim, unfolded linear_factors jordan_matrix_char_poly]
  have cp: "?cp1 = ?cp2" .
  show ?thesis
  proof (rule jordan_nf_matrix_poly_bound[OF A _ _ jnf])
    fix n a
    assume na: "(n,a)  set n_as"
    then obtain na1 na2 where n_as: "n_as = na1 @ (n,a) # na2"
      unfolding in_set_conv_decomp by auto
    then obtain p where "?cp2 = [: -a, 1 :]^n * p" unfolding n_as by auto
    from cp[unfolded this] have dvd: "[: -a, 1 :] ^ n dvd ?cp1" by auto
    let ?as = "filter ((=) a) as"
    let ?pn = "λ as. aas. [:- a, 1:]"
    let ?p = "λ as. aas. [: a, 1:]"
    have "?pn as = ?p (map uminus as)" by (induct as, auto)
    from poly_linear_exp_linear_factors[OF dvd[unfolded this]] 
    have "n  length (filter ((=) (- a)) (map uminus as))" .
    also have " = length (filter ((=) a) as)" 
      by (induct as, auto)
    finally have filt: "n  length (filter ((=) a) as)" .
    {
      assume "0 < n"
      with filt obtain b bs where "?as = b # bs" by (cases ?as, auto)
      from arg_cong[OF this, of set]
      have "a  set as" by auto 
      from le_1[rule_format, OF this]
      show "norm a  1" .
      note a  set as
    } note mem = this
    {
      assume "norm a = 1" 
      from le_N[OF mem this] filt show "n  N" by (cases n, auto)
    }
  qed
qed

end