Theory Jordan_Normal_Form_Existence

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

theory Jordan_Normal_Form_Existence
imports 
  Jordan_Normal_Form
  Column_Operations
  Schur_Decomposition
begin

hide_const (open) Coset.order

text‹We prove existence of Jordan normal forms by means of first applying Schur's algorithm
 to convert a matrix into upper-triangular form, and then applying the following algorithm 
 to convert a upper-triangular matrix into a Jordan normal form. It only consists of 
 basic row- and column-operations.›

subsection ‹Pseudo Code Algorithm›

text ‹The following algorithm is used to compute JNFs from upper-triangular matrices.
It was generalized from cite‹Sect.~11.1.4› in "PO07" where this algorithm was not
explicitly specified but only applied on an example. We further introduced step 2
which does not occur in the textbook description.

\begin{enumerate} 
\item[1.] Eliminate entries within blocks besides EV $a$ and above EV $b$ for $a \neq b$:
    for $A_{ij}$ with EV $a$ left of it, and EV $b$ below of it, perform
      @{term "add_col_sub_row (Aij / (b - a)) i j"}.
    The iteration should be by first increasing $j$ and the inner loop by decreasing $i$.
      
\item[2.] Move rows of same EV together, can only be done after 1., 
    otherwise triangular-property is lost.
    Say both rows $i$ and $j$ ($i < j$) contain EV $a$, but all rows between $i$ and $j$ have different EV.
    Then perform
      @{term "swap_cols_rows (i+1) j"}, 
      @{term "swap_cols_rows (i+2) j"}, 
      \ldots
      @{term "swap_cols_rows (j-1) j"}.
    Afterwards row $j$ will be at row $i+1$, and rows $i+1,\ldots,j-1$ will be moved to $i+2,\ldots,j$.
    The global iteration works by increasing $j$.

\item[3.] Transform each EV-block into JNF, do this for increasing upper $n \times k$ matrices, 
    where each new column $k$ will be treated as follows.
\begin{enumerate}
 \item[a)] Eliminate entries $A_{ik}$ in rows of form $0 \ldots 0\ ev\ 1\ 0 \ldots 0\ A_{ik}$:
         @{term "add_col_sub_row (-Aik) (i+1) k"}.
       Perform elimination by increasing $i$.
\item[b)] Figure out largest JB (of $n-1 \times n-1$ sub-matrix) with lowest row of form $0 \ldots 0\ ev\ 0 \ldots 0\ A_{lk}$
       where $A_{lk} \neq 0$, and set $x := A_{lk}$.
\item[c)] If such a JB does not exist, continue with next column.
  Otherwise, eliminate all other non-zero-entries $y := A_{ik}$ via row $l$:
         @{term "add_col_sub_row (y/x) i l"},
         @{term "add_col_sub_row (y/x) (i-1) (l-1)"},
         @{term "add_col_sub_row (y/x) (i-2) (l-2)"}, \ldots
         where the number of steps is determined by the size of the JB left-above of $A_{ik}$.
         Perform an iteration over $i$.
\item[d)] Normalize value in row $l$ to 1:
         @{term "mult_col_div_row (1/x) k"}.
\item[e)] Move the 1 down from row $l$ to row $k-1$:
         @{term "swap_cols_rows (l+1) k"},
         @{term "swap_cols_rows (l+2) k"},
         \ldots,
         @{term "swap_cols_rows (k-1) k"}.
\end{enumerate}
\end{enumerate}
›


subsection ‹Real Algorithm›

fun lookup_ev :: "'a  nat  'a mat  nat option" where
  "lookup_ev ev 0 A = None"
| "lookup_ev ev (Suc i) A = (if A $$ (i,i) = ev then Some i else lookup_ev ev i A)"

function swap_cols_rows_block :: "nat  nat  'a mat  'a mat" where
  "swap_cols_rows_block i j A = (if i < j then
    swap_cols_rows_block (Suc i) j (swap_cols_rows i j A) else A)"
  by pat_completeness auto
termination by (relation "measure (λ (i,j,A). j - i)") auto

fun identify_block :: "'a :: one mat  nat  nat" where
  "identify_block A 0 = 0"
| "identify_block A (Suc i) = (if A $$ (i,Suc i) = 1 then
    identify_block A i else (Suc i))"

function identify_blocks_main :: "'a :: ring_1 mat  nat  (nat × nat) list  (nat × nat) list" where
  "identify_blocks_main A 0 list = list"
| "identify_blocks_main A (Suc i_end) list = (
    let i_begin = identify_block A i_end
    in identify_blocks_main A i_begin ((i_begin, i_end) # list)
    )"
  by pat_completeness auto

definition identify_blocks :: "'a :: ring_1 mat  nat  (nat × nat)list" where
  "identify_blocks A i = identify_blocks_main A i []"

fun find_largest_block :: "nat × nat  (nat × nat)list  nat × nat" where
  "find_largest_block block [] = block"
| "find_largest_block (m_start,m_end) ((i_start,i_end) # blocks) = 
    (if i_end - i_start  m_end - m_start then
      find_largest_block (i_start,i_end) blocks else
      find_largest_block (m_start,m_end) blocks)"

fun lookup_other_ev :: "'a  nat  'a mat  nat option" where
  "lookup_other_ev ev 0 A = None"
| "lookup_other_ev ev (Suc i) A = (if A $$ (i,i)  ev then Some i else lookup_other_ev ev i A)"

partial_function (tailrec) partition_ev_blocks :: "'a mat  'a mat list  'a mat list" where
  [code]: "partition_ev_blocks A bs = (let n = dim_row A in
    if n = 0 then bs
    else (case lookup_other_ev (A $$ (n-1, n-1)) (n-1) A of 
      None  A # bs 
    | Some i  case split_block A (Suc i) (Suc i) of (UL,_,_,LR)  partition_ev_blocks UL (LR # bs)))"

context 
  fixes n :: nat
  and ty :: "'a :: field itself"
begin

function step_1_main :: "nat  nat  'a mat  'a mat" where
  "step_1_main i j A = (if j  n then A else if i = 0 then step_1_main (j+1) (j+1) A
    else let 
      i' = i - 1;
      ev_left = A $$ (i',i');
      ev_below = A $$ (j,j);
      aij = A $$ (i',j);
      B = if (ev_left  ev_below  aij  0) then add_col_sub_row (aij / (ev_below - ev_left)) i' j A else A
     in step_1_main i' j B)"
  by pat_completeness auto
termination by (relation "measures [λ (i,j,A). n - j, λ (i,j,A). i]") auto

function step_2_main :: "nat  'a mat  'a mat" where
  "step_2_main j A = (if j  n then A 
    else 
      let ev = A $$ (j,j);
        B = (case lookup_ev ev j A of 
          None  A
        | Some i  swap_cols_rows_block (Suc i) j A
          )
     in step_2_main (Suc j) B)"
  by pat_completeness auto
termination by (relation "measure (λ (j,A). n - j)") auto

fun step_3_a :: "nat  nat  'a mat  'a mat" where 
  "step_3_a 0 j A = A"
| "step_3_a (Suc i) j A = (let 
    aij = A $$ (i,j);
    B = (if A $$ (i,i+1) = 1  aij  0 
    then add_col_sub_row (- aij) (Suc i) j A else A)
    in step_3_a i j B)"

fun step_3_c_inner_loop :: "'a  nat  nat  nat  'a mat  'a mat" where
  "step_3_c_inner_loop val l i 0 A = A"
| "step_3_c_inner_loop val l i (Suc k) A = step_3_c_inner_loop val (l - 1) (i - 1) k 
     (add_col_sub_row val i l A)"

fun step_3_c :: "'a  nat  nat  (nat × nat)list  'a mat  'a mat" where
  "step_3_c x l k [] A = A"
| "step_3_c x l k ((i_begin,i_end) # blocks) A = (
    let 
      B = (if i_end = l then A else 
        step_3_c_inner_loop (A $$ (i_end,k) / x) l i_end (Suc i_end - i_begin) A)
      in step_3_c x l k blocks B)"

function step_3_main :: "nat  'a mat  'a mat" where
  "step_3_main k A = (if k  n then A 
    else let 
      B = step_3_a (k-1) k A; ― ‹3_a›
      all_blocks = identify_blocks B k;
      blocks = filter (λ block. B $$ (snd block,k)  0) all_blocks;
      F = (if blocks = [] ― ‹column k› has only 0s›
        then B
        else let          
          (l_start,l) = find_largest_block (hd blocks) (tl blocks); ― ‹3_b›
          x = B $$ (l,k); 
          C = step_3_c x l k blocks B; ― ‹3_c›
          D = mult_col_div_row (inverse x) k C; ― ‹3_d›
          E = swap_cols_rows_block (Suc l) k D ― ‹3_e›
        in E)
      in step_3_main (Suc k) F)"
  by pat_completeness auto
termination by (relation "measure (λ (k,A). n - k)") auto

end

definition step_1 :: "'a :: field mat  'a mat" where
  "step_1 A = step_1_main (dim_row A) 0 0 A"

definition step_2 :: "'a :: field mat  'a mat" where
  "step_2 A = step_2_main (dim_row A) 0 A"

definition step_3 :: "'a :: field mat  'a mat" where
  "step_3 A = step_3_main (dim_row A) 1 A"

declare swap_cols_rows_block.simps[simp del]
declare step_1_main.simps[simp del]
declare step_2_main.simps[simp del]
declare step_3_main.simps[simp del]

function jnf_vector_main :: "nat  'a :: one mat  (nat × 'a) list" where
  "jnf_vector_main 0 A = []"
| "jnf_vector_main (Suc i_end) A = (let 
    i_start = identify_block A i_end 
    in jnf_vector_main i_start A @ [(Suc i_end - i_start, A $$ (i_start,i_start))])"
  by pat_completeness auto

definition jnf_vector :: "'a :: one mat  (nat × 'a) list" where
  "jnf_vector A = jnf_vector_main (dim_row A) A"

definition triangular_to_jnf_vector :: "'a :: field mat  (nat × 'a) list" where
  "triangular_to_jnf_vector A  let B = step_2 (step_1 A)
    in concat (map (jnf_vector o step_3) (partition_ev_blocks B []))"

subsection ‹Preservation of Dimensions›

lemma swap_cols_rows_block_dims_main: 
  "dim_row (swap_cols_rows_block i j A) = dim_row A  dim_col (swap_cols_rows_block i j A) = dim_col A"
proof (induct i j A rule: swap_cols_rows_block.induct)
  case (1 i j A)
  thus ?case unfolding swap_cols_rows_block.simps[of i j A]
    by (auto split: if_splits)
qed

lemma swap_cols_rows_block_dims[simp]: 
  "dim_row (swap_cols_rows_block i j A) = dim_row A"
  "dim_col (swap_cols_rows_block i j A) = dim_col A"
  "A  carrier_mat n n  swap_cols_rows_block i j A  carrier_mat n n"
  using swap_cols_rows_block_dims_main unfolding carrier_mat_def by auto

lemma step_1_main_dims_main: 
  "dim_row (step_1_main n i j A) = dim_row A  dim_col (step_1_main n i j A) = dim_col A"
proof (induct i j A taking: n rule: step_1_main.induct)
  case (1 i j A)
  thus ?case unfolding step_1_main.simps[of n i j A]
    by (auto split: if_splits simp: Let_def)
qed

lemma step_1_main_dims[simp]: 
  "dim_row (step_1_main n i j A) = dim_row A"
  "dim_col (step_1_main n i j A) = dim_col A"
  using step_1_main_dims_main by blast+

lemma step_2_main_dims_main: 
  "dim_row (step_2_main n j A) = dim_row A  dim_col (step_2_main n j A) = dim_col A"
proof (induct j A taking: n rule: step_2_main.induct)
  case (1 j A)
  thus ?case unfolding step_2_main.simps[of n j A]
    by (auto split: if_splits option.splits simp: Let_def)
qed

lemma step_2_main_dims[simp]: 
  "dim_row (step_2_main n j A) = dim_row A"
  "dim_col (step_2_main n j A) = dim_col A"
  using step_2_main_dims_main by blast+

lemma step_3_a_dims_main: 
  "dim_row (step_3_a i j A) = dim_row A  dim_col (step_3_a i j A) = dim_col A"
  by (induct i j A rule: step_3_a.induct, auto split: if_splits simp: Let_def)

lemma step_3_a_dims[simp]: 
  "dim_row (step_3_a i j A) = dim_row A"
  "dim_col (step_3_a i j A) = dim_col A"
  using step_3_a_dims_main by blast+

lemma step_3_c_inner_loop_dims_main: 
  "dim_row (step_3_c_inner_loop val l i j A) = dim_row A  dim_col (step_3_c_inner_loop val l i j A) = dim_col A"
  by (induct val l i j A rule: step_3_c_inner_loop.induct, auto)

lemma step_3_c_inner_loop_dims[simp]: 
  "dim_row (step_3_c_inner_loop val l i j A) = dim_row A"
  "dim_col (step_3_c_inner_loop val l i j A) = dim_col A"
  using step_3_c_inner_loop_dims_main by blast+

lemma step_3_c_dims_main: 
  "dim_row (step_3_c x l k i A) = dim_row A  dim_col (step_3_c x l k i A) = dim_col A"
  by (induct x l k i A rule: step_3_c.induct, auto simp: Let_def)

lemma step_3_c_dims[simp]: 
  "dim_row (step_3_c x l k i A) = dim_row A"
  "dim_col (step_3_c x l k i A) = dim_col A"
  using step_3_c_dims_main by blast+
  
lemma step_3_main_dims_main: 
  "dim_row (step_3_main n k A) = dim_row A  dim_col (step_3_main n k A) = dim_col A"
proof (induct k A taking: n rule: step_3_main.induct)
  case (1 k A)
  thus ?case unfolding step_3_main.simps[of n k A]
    by (auto split: if_splits prod.splits option.splits simp: Let_def)
qed

lemma step_3_main_dims[simp]: 
  "dim_row (step_3_main n j A) = dim_row A"
  "dim_col (step_3_main n j A) = dim_col A"
  using step_3_main_dims_main by blast+

lemma triangular_to_jnf_steps_dims[simp]: 
  "dim_row (step_1 A) = dim_row A"
  "dim_col (step_1 A) = dim_col A"
  "dim_row (step_2 A) = dim_row A"
  "dim_col (step_2 A) = dim_col A"
  "dim_row (step_3 A) = dim_row A"
  "dim_col (step_3 A) = dim_col A"
  unfolding step_1_def step_2_def step_3_def o_def by auto

subsection ‹Properties of Auxiliary Algorithms›

lemma lookup_ev_Some: 
  "lookup_ev ev j A = Some i  
  i < j  A $$ (i,i) = ev  ( k. i < k  k < j  A $$ (k,k)  ev)"
  by (induct j, auto split: if_splits, case_tac "k = j", auto)

lemma lookup_ev_None: "lookup_ev ev j A = None  i < j  A $$ (i,i)  ev"
  by (induct j, auto split: if_splits, case_tac "i = j", auto)

lemma swap_cols_rows_block_index[simp]: 
  "i < dim_row A  i < dim_col A  j < dim_row A  j < dim_col A 
   low  high  high < dim_row A  high < dim_col A  
   swap_cols_rows_block low high A $$ (i,j) = A $$ 
    (if i = low then high else if i > low  i  high then i - 1 else i,
     if j = low then high else if j > low  j  high then j - 1 else j)"
proof (induct low high A rule: swap_cols_rows_block.induct)
  case (1 low high A)
  let ?r = "dim_row A" let ?c = "dim_col A"
  let ?A = "swap_cols_rows_block low high A"
  let ?B = "swap_cols_rows low high A"
  let ?C = "swap_cols_rows_block (Suc low) high ?B"
  note [simp] = swap_cols_rows_block.simps[of low high A]
  from 1(2-) have lh: "low  high" by simp
  show ?case
  proof (cases "low < high")
    case False
    with lh have lh: "low = high" by auto
    from False have "?A = A" by simp
    with lh show ?thesis by auto
  next
    case True
    hence A: "?A = ?C" by simp
    from True lh have "Suc low  high" by simp
    note IH = 1(1)[unfolded swap_cols_rows_carrier, OF True 1(2-5) this 1(7-)]
    note * = 1(2-)
    let ?i = "if i = Suc low then high else if Suc low < i  i  high then i - 1 else i"
    let ?j = "if j = Suc low then high else if Suc low < j  j  high then j - 1 else j"
    have cong: " i j i' j'. i = i'  j = j'  A $$ (i,j) = A $$ (i',j')" by auto
    have ij: "?i < ?r" "?i < ?c" "?j < ?r" "?j < ?c" "low < ?r" "high < ?r" using * True by auto
    show ?thesis unfolding A IH
      by (subst swap_cols_rows_index[OF ij], rule cong, insert * True, auto)
  qed
qed
  
lemma find_largest_block_main: assumes "find_largest_block block blocks = (m_b, m_e)"
  shows "(m_b, m_e)  insert block (set blocks)
   ( b  insert block (set blocks). m_e - m_b  snd b - fst b)"
  using assms
proof (induct block blocks rule: find_largest_block.induct)
  case (2 m_start m_end i_start i_end blocks)  
  let ?res = "find_largest_block (m_start,m_end) ((i_start,i_end) # blocks)"
  show ?case
  proof (cases "i_end - i_start  m_end - m_start")
    case True
    with 2(3-) have "find_largest_block (i_start,i_end) blocks = (m_b, m_e)" by auto
    note IH = 2(1)[OF True this]
    thus ?thesis using True by auto
  next
    case False
    with 2(3-) have "find_largest_block (m_start,m_end) blocks = (m_b, m_e)" by auto
    note IH = 2(2)[OF False this]
    thus ?thesis using False by auto
  qed
qed simp

lemma find_largest_block: assumes bl: "blocks  []"
  and find: "find_largest_block (hd blocks) (tl blocks) = (m_begin, m_end)"
  shows "(m_begin,m_end)  set blocks"
  " i_begin i_end. (i_begin,i_end)  set blocks  m_end - m_begin  i_end - i_begin"
proof -
  from bl have id: "insert (hd blocks) (set (tl blocks)) = set blocks" by (cases blocks, auto)
  from find_largest_block_main[OF find, unfolded id] 
  show "(m_begin,m_end)  set blocks"
    " i_begin i_end. (i_begin,i_end)  set blocks  m_end - m_begin  i_end - i_begin" by auto
qed

context 
  fixes ev :: "'a :: one"
  and A :: "'a mat"
begin

lemma identify_block_main: assumes "identify_block A j = i"
  shows "i  j  (i = 0  A $$ (i - 1, i)  1)  ( k. i  k  k < j  A $$ (k, Suc k) = 1)"
    (is "?P j")
  using assms
proof (induct j)
  case (Suc j)
  show ?case
  proof (cases "A $$ (j, Suc j) = 1")
    case False
    with Suc(2) show ?thesis by auto
  next
    case True
    with Suc(2) have "identify_block A j = i" by simp
    note IH = Suc(1)[OF this] 
    {
      fix k
      assume "k  i" and "k < Suc j"      
      hence "A $$ (k, Suc k) = 1"
      proof (cases "k < j")
        case True
        with IH k  i show ?thesis by auto
      next
        case False
        with k < Suc j have "k = j" by auto
        with True show ?thesis by auto
      qed
    }
    with IH show ?thesis by auto
  qed
qed simp


lemma identify_block_le: "identify_block A i  i"
  using identify_block_main[OF refl] by blast
end


lemma identify_block: assumes "identify_block A j = i"
  shows "i  j"
  "i = 0  A $$ (i - 1, i)  1"
  "i  k  k < j  A $$ (k, Suc k) = 1"
proof -
  let ?ev = "A $$ (j,j)"
  note main = identify_block_main[OF assms]
  from main show "i  j" by blast
  from main show "i = 0  A $$ (i - 1, i)  1" by blast
  assume "i  k"
  with main have main: "k < j  A $$ (k, Suc k) = 1" by blast
  thus "k < j  A $$ (k, Suc k) = 1" by blast
qed
    
lemmas identify_block_le' = identify_block(1)

lemma identify_block_le_rev: "j = identify_block A i  j  i"
  using identify_block_le'[of A i j] by auto
  
termination identify_blocks_main by (relation "measure (λ (_,i,list). i)", 
  auto simp: identify_block_le le_imp_less_Suc)

termination jnf_vector_main by (relation "measure (λ (i,A). i)", 
  auto simp: identify_block_le le_imp_less_Suc)

lemma identify_blocks_main: assumes "(i_start,i_end)  set (identify_blocks_main A i list)" 
  and " i_s i_e. (i_s, i_e)  set list  i_s  i_e  i_e < k"
  and "i  k"
  shows "i_start  i_end  i_end < k" using assms
proof (induct A i list rule: identify_blocks_main.induct)
  case (2 A i_e list)
  obtain i_b where id: "identify_block A i_e = i_b" by force
  note IH = 2(1)[OF id[symmetric]]
  let ?res = "identify_blocks_main A (Suc i_e) list"  
  let ?rec = "identify_blocks_main A i_b ((i_b, i_e) # list)"
  have res: "?res = ?rec" using id by simp
  from 2(2)[unfolded res] have "(i_start, i_end)  set ?rec" .
  note IH = IH[OF this]
  from 2(3-) have iek: "i_e < k" by simp
  from identify_block_le'[OF id] have ibe: "i_b  i_e" .
  from ibe iek have "i_b  k" by simp
  note IH = IH[OF _ this]
  show ?thesis
    by (rule IH, insert ibe iek 2(3-), auto)
qed simp

lemma identify_blocks: assumes "(i_start,i_end)  set (identify_blocks B k)" 
  shows "i_start  i_end  i_end < k"
  using identify_blocks_main[OF assms[unfolded identify_blocks_def]] by auto

subsection ‹Proving Similarity›

context
begin
private lemma swap_cols_rows_block_similar: assumes "A  carrier_mat n n"
  and "j < n" and "i  j"
  shows "similar_mat (swap_cols_rows_block i j A) A"
  using assms
proof (induct i j A rule: swap_cols_rows_block.induct)
  case (1 i j A)
  hence A: "A  carrier_mat n n"
    and jn: "j < n" and ij: "i  j" by auto
  note [simp] = swap_cols_rows_block.simps[of i j]
  show ?case
  proof (cases "i < j")
    case False
    thus ?thesis using similar_mat_refl[OF A] by simp
  next
    case True note ij = this
    let ?B = "swap_cols_rows i j A"
    let ?C = "swap_cols_rows_block (Suc i) j ?B"
    from swap_cols_rows_similar[OF A _ jn, of i] ij jn
    have BA: "similar_mat ?B A" by auto
    have CB: "similar_mat ?C ?B"
      by (rule 1(1)[OF ij _ jn], insert ij A, auto)
    from similar_mat_trans[OF CB BA] show ?thesis using ij by simp
  qed
qed

private lemma step_1_main_similar: "i  j  A  carrier_mat n n  similar_mat (step_1_main n i j A) A"
proof (induct i j A taking: n rule: step_1_main.induct)
  case (1 i j A)
  note [simp] = step_1_main.simps[of n i j A]
  from 1(3-) have ij: "i  j" and A: "A  carrier_mat n n" by auto
  show ?case
  proof (cases "j  n")
    case True
    thus ?thesis using similar_mat_refl[OF A] by simp
  next
    case False 
    hence jn: "j < n" by simp
    note IH = 1(1-2)[OF False]
    show ?thesis
    proof (cases "i = 0")
      case True
      from IH(1)[OF this _ A]
      show ?thesis using jn by (simp, simp add: True)
    next
      case False
      let ?evi = "A $$ (i - 1, i - 1)"
      let ?evj = "A $$ (j,j)"
      let ?B = "if ?evi  ?evj  A $$ (i - 1, j)  0 then 
        add_col_sub_row (A $$ (i - 1, j) / (?evj - ?evi)) (i - 1) j A else A"
      obtain B where B: "B = ?B" by auto
      have Bn: "B  carrier_mat n n" unfolding B using A by simp
      from False ij jn have *: "i - 1 < n" "j < n" "i - 1  j" by auto
      have BA: "similar_mat B A" unfolding B using similar_mat_refl[OF A]
        add_col_sub_row_similar[OF A *] by auto
      from ij have "i - 1  j" by simp
      note IH = IH(2)[OF False refl refl refl refl B this Bn]
      from False jn have id: "step_1_main n i j A = step_1_main n (i - 1) j B"
        unfolding B by (simp add: Let_def)
      show ?thesis unfolding id
        by (rule similar_mat_trans[OF IH BA])
    qed
  qed
qed

private lemma step_2_main_similar: "A  carrier_mat n n  similar_mat (step_2_main n j A) A"
proof (induct j A taking: n rule: step_2_main.induct)
  case (1 j A)
  note [simp] = step_2_main.simps[of n j A]
  from 1(2) have A: "A  carrier_mat n n" .
  show ?case
  proof (cases "j  n")
    case True
    thus ?thesis using similar_mat_refl[OF A] by simp
  next
    case False 
    hence jn: "j < n" by simp
    note IH = 1(1)[OF False]
    let ?look = "lookup_ev (A $$ (j,j)) j A"
    let ?B = "case ?look of 
          None  A
        | Some i  swap_cols_rows_block (Suc i) j A"
    obtain B where B: "B = ?B" by auto
    have id: "step_2_main n j A = step_2_main n (Suc j) B" unfolding B using False by simp
    have Bn: "B  carrier_mat n n" unfolding B using A by (auto split: option.splits)
    have BA: "similar_mat B A" 
    proof (cases ?look)
      case None
      thus ?thesis unfolding B using similar_mat_refl[OF A] by simp
    next
      case (Some i)
      from lookup_ev_Some[OF this]
      show ?thesis unfolding B Some
        by (auto intro: swap_cols_rows_block_similar[OF A jn])
    qed
    show ?thesis unfolding id
      by (rule similar_mat_trans[OF IH[OF refl B Bn] BA])
  qed
qed

private lemma step_3_a_similar: "A  carrier_mat n n  i < j  j < n  similar_mat (step_3_a i j A) A"
proof (induct i j A rule: step_3_a.induct)
  case (1 j A)
  thus ?case by (simp add: similar_mat_refl)
next
  case (2 i j A)
  from 2(2-) have A: "A  carrier_mat n n" and ij: "Suc i < j" and j: "j < n" by auto
  let ?B = "if A $$ (i, i + 1) = 1  A $$ (i, j)  0 
    then add_col_sub_row (- A $$ (i, j)) (Suc i) j A else A"
  obtain B where B: "B = ?B" by auto
  from A have Bn: "B  carrier_mat n n" unfolding B by simp
  note IH = 2(1)[OF refl B Bn _ j]
  have id: "step_3_a (Suc i) j A = step_3_a i j B" unfolding B by (simp add: Let_def)
  from ij j have *: "Suc i < n" "j < n" "Suc i  j" by auto
  have BA: "similar_mat B A" unfolding B
    using similar_mat_refl[OF A] add_col_sub_row_similar[OF A *] by auto
  show ?case unfolding id
    by (rule similar_mat_trans[OF IH BA], insert ij, auto)
qed

private lemma step_3_c_inner_loop_similar: 
  "A  carrier_mat n n  l  i  k - 1  l  k - 1  i  l < n  i < n  
  similar_mat (step_3_c_inner_loop val l i k A) A"
proof (induct val l i k A rule: step_3_c_inner_loop.induct)
  case (1 val l i A)
  thus ?case using similar_mat_refl[of A] by simp
next
  case (2 val l i k A)
  let ?res = "step_3_c_inner_loop val l i (Suc k) A"
  from 2(2-) have A: "A  carrier_mat n n" and 
    *: "l  i" "k  l" "k  i" "l < n" "i < n" by auto
  let ?B = "add_col_sub_row val i l A"
  have BA: "similar_mat ?B A"
    by (rule add_col_sub_row_similar[OF A], insert *, auto)
  have B: "?B  carrier_mat n n" using A unfolding carrier_mat_def by simp
  show ?case
  proof (cases k)
    case 0
    hence id: "?res = ?B" by simp
    thus ?thesis using BA by simp
  next
    case (Suc kk)
    with * have "l - 1  i - 1" "k - 1  l - 1" "k - 1  i - 1" "l - 1 < n" "i - 1 < n" by auto
    note IH = 2(1)[OF B this]
    show ?thesis unfolding step_3_c_inner_loop.simps
      by (rule similar_mat_trans[OF IH BA])
  qed
qed

private lemma step_3_c_similar: 
  "A  carrier_mat n n  l < k  k < n 
   ( i_begin i_end. (i_begin, i_end)  set blocks   i_end  k  i_end - i_begin  l)
   similar_mat (step_3_c x l k blocks A) A"
proof (induct x l k blocks A rule: step_3_c.induct)
  case (1 x l k A)
  thus ?case using similar_mat_refl[OF 1(1)] by simp
next
  case (2 x l k i_begin i_end blocks A)
  let ?res = "step_3_c x l k ((i_begin,i_end) # blocks) A"
  from 2(2-4) have A: "A  carrier_mat n n" and lki: "l < k" "k < n" by auto
  from 2(5) have i: "i_end  k" "i_end - i_begin  l" by auto
  let ?y = "A $$ (i_end,k)"
  let ?inner = "step_3_c_inner_loop (?y / x) l i_end (Suc i_end - i_begin) A"
  obtain B where B: 
    "B = (if i_end = l then A else ?inner)" by auto    
  hence id: "?res = step_3_c x l k blocks B"
    by simp
  have BA: "similar_mat B A" 
  proof (cases "i_end = l")
    case True
    thus ?thesis unfolding B using similar_mat_refl[OF A] by simp
  next
    case False
    hence B: "B = ?inner" and li: "l  i_end" by (auto simp: B)      
    show ?thesis unfolding B 
      by (rule step_3_c_inner_loop_similar[OF A li], insert lki i, auto)
  qed
  have Bn: "B  carrier_mat n n" using A unfolding B carrier_mat_def by simp
  note IH = 2(1)[OF B Bn lki(1-2) 2(5)]
  show ?case unfolding id
    by (rule similar_mat_trans[OF IH BA], auto)
qed

private lemma step_3_main_similar: "A  carrier_mat n n  k > 0  similar_mat (step_3_main n k A) A"
proof (induct k A taking: n rule: step_3_main.induct)
  case (1 k A)
  from 1(2-) have A: "A  carrier_mat n n" and k: "k > 0" by auto
  note [simp] = step_3_main.simps[of n k A]
  show ?case
  proof (cases "k  n")
    case True
    thus ?thesis using similar_mat_refl[OF A] by simp
  next
    case False
    hence kn: "k < n" by simp
    obtain B where B: "B = step_3_a (k - 1) k A" by auto
    note IH = 1(1)[OF False B]
    from A have Bn: "B  carrier_mat n n" unfolding B carrier_mat_def by simp
    from k have "k - 1 < k" by simp
    from step_3_a_similar[OF A this kn] have BA: "similar_mat B A" unfolding B .
    obtain all_blocks where ab: "all_blocks = identify_blocks B k" by simp
    obtain blocks where blocks: "blocks = filter (λ block. B $$ (snd block, k)  0) all_blocks" by simp
    obtain F where F: "F = (if blocks = [] then B
       else let (l_begin,l) = find_largest_block (hd blocks) (tl blocks); x = B $$ (l, k); C = step_3_c x l k blocks B;
            D = mult_col_div_row (inverse x) k C; E = swap_cols_rows_block (Suc l) k D
        in E)" by simp
    note IH = IH[OF ab blocks F]
    have Fn: "F  carrier_mat n n" unfolding F Let_def carrier_mat_def using Bn 
      by (simp split: prod.splits)
    have FB: "similar_mat F B" 
    proof (cases "blocks = []")
      case True
      thus ?thesis unfolding F using similar_mat_refl[OF Bn] by simp
    next
      case False
      obtain l_start l where l: "find_largest_block (hd blocks) (tl blocks) = (l_start, l)" by force
      obtain x where x: "x = B $$ (l,k)" by simp
      obtain C where C: "C = step_3_c x l k blocks B" by simp
      obtain D where D: "D = mult_col_div_row (inverse x) k C" by auto
      obtain E where E: "E = swap_cols_rows_block (Suc l) k D" by auto
      from find_largest_block[OF False l] have lb: "(l_start,l)  set blocks"
        and llarge: " i_begin i_end. (i_begin,i_end)  set blocks  l - l_start  i_end - i_begin" by auto
      from lb have x0: "x  0" unfolding blocks x by simp
      {
        fix i_start i_end
        assume "(i_start,i_end)  set blocks"
        hence "(i_start,i_end)  set (identify_blocks B k)" unfolding blocks ab by simp
        from identify_blocks[OF this]
        have "i_end < k" by simp
      } note block_bound = this
      from block_bound[OF lb]
      have lk: "l < k" .
      from False have F: "F = E" unfolding E D C x F l Let_def by simp
      from Bn have Cn: "C  carrier_mat n n" unfolding C carrier_mat_def by simp
      have CB: "similar_mat C B" unfolding C
      proof (rule step_3_c_similar[OF Bn lk kn])
        fix i_begin i_end
        assume i: "(i_begin, i_end)  set blocks"
        from llarge[OF i] block_bound[OF i] 
        show "i_end  k  i_end - i_begin  l" by auto
      qed
      from x0 have "inverse x  0" by simp
      from mult_col_div_row_similar[OF Cn kn this] 
      have DC: "similar_mat D C" unfolding D .
      from Cn have Dn: "D  carrier_mat n n" unfolding D carrier_mat_def by simp
      from lk have "Suc l  k" by auto
      from swap_cols_rows_block_similar[OF Dn kn this] 
      have ED: "similar_mat E D" unfolding E .
      from similar_mat_trans[OF ED similar_mat_trans[OF DC 
        similar_mat_trans[OF CB similar_mat_refl[OF Bn]]]]
      show ?thesis unfolding F .
    qed
    have "0 < Suc k" by simp
    note IH = IH[OF Fn this]
    have id: "step_3_main n k A = step_3_main n (Suc k) F" using kn 
      by (simp add: F Let_def blocks ab B)
    show ?thesis unfolding id
      by (rule similar_mat_trans[OF IH similar_mat_trans[OF FB BA]])
  qed
qed

lemma step_1_similar: "A  carrier_mat n n  similar_mat (step_1 A) A"
  unfolding step_1_def by (rule step_1_main_similar, auto)

lemma step_2_similar: "A  carrier_mat n n  similar_mat (step_2 A) A"
  unfolding step_2_def by (rule step_2_main_similar, auto)

lemma step_3_similar: "A  carrier_mat n n  similar_mat (step_3 A) A"
  unfolding step_3_def by (rule step_3_main_similar, auto)

end


subsection ‹Invariants for Proving that Result is in JNF›
context 
  fixes n :: nat
  and ty :: "'a :: field itself"
begin

(* upper triangular, ensured by precondition and then maintained *)
definition uppert :: "'a mat  nat  nat  bool" where
  "uppert A i j  j < i  A $$ (i,j) = 0" 

(* zeros at places where EVs differ, ensured by step 1 and then maintained *)
definition diff_ev :: "'a mat  nat  nat  bool" where
  "diff_ev A i j  i < j  A $$ (i,i)  A $$ (j,j)  A $$ (i,j) = 0"

(* same EVs are grouped together, ensured by step 2 and then maintained *)
definition ev_blocks_part :: "nat  'a mat  bool" where
  "ev_blocks_part m A   i j k. i < j  j < k  k < m  A $$ (k,k) = A $$ (i,i)  A $$ (j,j) = A $$ (i,i)"

definition ev_blocks :: "'a mat  bool" where
  "ev_blocks  ev_blocks_part n"

text ‹In step 3, there is a separation at which iteration we are.
  The columns left of $k$ will be in JNF, the columns right of $k$ or equal to $k$ will satisfy @{const uppert}, @{const diff_ev}, 
  and @{const ev_blocks}, and the column at $k$ will have one of the following properties, which are ensured in the
  different phases of step 3.›

(* Left of a one is a zero: In a row of shape 0 ... 0 ev 1 0 ... 0 entry, the entry will be 0,
   ensured by step 3 a and then maintained *)
private definition one_zero :: "'a mat  nat  nat  bool" where
  "one_zero A i j  
    (Suc i < j  A $$ (i,Suc i) = 1  A $$ (i,j) = 0)  
    (j < i  A $$ (i,j) = 0) 
    (i < j  A $$ (i,i)  A $$ (j,j)  A $$ (i,j) = 0)"

(* There is exactly one row   0 ... 0 ev 0 ... 0 entry with entry != 0 and that entry is x,
   ensured by step 3 c and then maintained *)
private definition single_non_zero :: "nat  nat  'a  'a mat  nat  nat  bool" where
  "single_non_zero  λ l k x. (λ A i j. (i  {k,l}  A $$ (i,k) = 0)  A $$ (l,k) = x)"

(* There is at exactly one row   0 ... 0 ev 0 ... 0 entry with entry != 0 and that entry is 1,
   ensured by step 3 d and then maintained *)
private definition single_one :: "nat  nat  'a mat  nat  nat  bool" where
  "single_one  λ l k. (λ A i j. (i  {k,l}  A $$ (i,k) = 0)  A $$ (l,k) = 1)"

(* there is at most one row   0 ... 0 ev 0 ... 0 entry with entry != 0 and that entry is 1 and there
   are no zeros between ev and the entry, ensured by step 3 e *)
private definition lower_one :: "nat  'a mat  nat  nat  bool" where
  "lower_one k A i j = (j = k  
    (A $$ (i,j) = 0  i = j  (A $$ (i,j) = 1  j = Suc i  A $$ (i,i) = A $$ (j,j))))"

(* the desired property: we have a jordan block matrix *)
definition jb :: "'a mat  nat  nat  bool" where
  "jb A i j  (Suc i = j  A $$ (i,j)  {0,1}) 
   (i  j  (Suc i  j  A $$ (i,i)  A $$ (j,j))  A $$ (i,j) = 0)"

text ‹The following properties are useful to easily ensure the above invariants 
  just from invariants of other matrices. The properties are essential in showing
  that the blocks identified in step 3b are the same as one would identify for
  the matrices in the upcoming steps 3c and 3d.›
 
definition same_diag :: "'a mat  'a mat  bool" where
  "same_diag A B   i < n. A $$ (i,i) = B $$ (i,i)"

private definition same_upto :: "nat  'a mat  'a mat  bool" where
  "same_upto j A B   i' j'. i' < n  j' < j  A $$ (i',j') = B $$ (i',j')"

text ‹Definitions stating where the properties hold›

definition inv_all :: "('a mat  nat  nat  bool)  'a mat  bool" where
  "inv_all p A   i j. i < n  j < n  p A i j"

private definition inv_part :: "('a mat  nat  nat  bool)  'a mat  nat  nat  bool" where
  "inv_part p A m_i m_j   i j. i < n  j < n  j < m_j  j = m_j  i  m_i  p A i j"

private definition inv_upto :: "('a mat  nat  nat  bool)  'a mat  nat  bool" where
  "inv_upto p A m   i j. i < n  j < n  j < m  p A i j"

private definition inv_from :: "('a mat  nat  nat  bool)  'a mat  nat  bool" where
  "inv_from p A m   i j. i < n  j < n  j > m  p A i j"

private definition inv_at :: "('a mat  nat  nat  bool)  'a mat  nat  bool" where
  "inv_at p A m   i. i < n  p A i m"

private definition inv_from_bot :: "('a mat  nat  bool)  'a mat  nat  bool" where
  "inv_from_bot p A mi   i. i  mi  i < n  p A i"

text ‹Auxiliary Lemmas on Handling, Comparing, and Accessing Invariants›

lemma jb_imp_uppert: "jb A i j  uppert A i j"
  unfolding jb_def uppert_def by auto

private lemma ev_blocks_partD:
  "ev_blocks_part m A  i < j  j < k  k < m  A $$ (k,k) = A $$ (i,i)  A $$ (j,j) = A $$ (i,i)"
  unfolding ev_blocks_part_def by auto

private lemma ev_blocks_part_leD:
  assumes "ev_blocks_part m A"
  "i  j" "j  k" "k < m" "A $$ (k,k) = A $$ (i,i)" 
  shows "A $$ (j,j) = A $$ (i,i)"
proof -  
  show ?thesis
  proof (cases "i = j  j = k")
    case False
    with assms(2-3) have "i < j" "j < k" by auto
    from ev_blocks_partD[OF assms(1) this assms(4-)] show ?thesis .
  next
    case True
    thus ?thesis using assms(5) by auto
  qed
qed

private lemma ev_blocks_partI:
  assumes " i j k. i < j  j < k  k < m  A $$ (k,k) = A $$ (i,i)  A $$ (j,j) = A $$ (i,i)"
  shows "ev_blocks_part m A"
  using assms unfolding ev_blocks_part_def by blast

private lemma ev_blocksD:
  "ev_blocks A  i < j  j < k  k < n  A $$ (k,k) = A $$ (i,i)  A $$ (j,j) = A $$ (i,i)"
  unfolding ev_blocks_def by (rule ev_blocks_partD)

private lemma ev_blocks_leD:
  "ev_blocks A  i  j  j  k  k < n  A $$ (k,k) = A $$ (i,i)  A $$ (j,j) = A $$ (i,i)"
  unfolding ev_blocks_def by (rule ev_blocks_part_leD)

lemma inv_allD: "inv_all p A  i < n  j < n  p A i j"
  unfolding inv_all_def by auto

private lemma inv_allI: assumes " i j. i < n  j < n  p A i j"
  shows "inv_all p A" 
  using assms unfolding inv_all_def by blast

private lemma inv_partI: assumes " i j. i < n  j < n  j < m_j  j = m_j  i  m_i  p A i j"
  shows "inv_part p A m_i m_j"
  using assms unfolding inv_part_def by auto

private lemma inv_partD: assumes "inv_part p A m_i m_j" "i < n" "j < n" 
  shows "j < m_j  p A i j"
  and "j = m_j  i  m_i  p A i j"
  and "j < m_j  j = m_j  i  m_i  p A i j"
  using assms unfolding inv_part_def by auto

private lemma inv_uptoI: assumes " i j. i < n  j < n  j < m  p A i j"
  shows "inv_upto p A m"
  using assms unfolding inv_upto_def by auto

private lemma inv_uptoD: assumes "inv_upto p A m" "i < n" "j < n" "j < m"
  shows "p A i j"
  using assms unfolding inv_upto_def by auto

private lemma inv_upto_Suc: assumes "inv_upto p A m"
  and " i. i < n  p A i m"
  shows "inv_upto p A (Suc m)"
proof (intro inv_uptoI)
  fix i j
  assume "i < n" "j < n" "j < Suc m"
  thus "p A i j" using inv_uptoD[OF assms(1), of i j] assms(2)[of i] by (cases "j = m", auto)
qed

private lemma inv_upto_mono: assumes " i j. i < n  j < k  p A i j  q A i j"
  shows "inv_upto p A k  inv_upto q A k"
  using assms unfolding inv_upto_def by auto

private lemma inv_fromI: assumes " i j. i < n  j < n  j > m  p A i j"
  shows "inv_from p A m"
  using assms unfolding inv_from_def by auto

private lemma inv_fromD: assumes "inv_from p A m" "i < n" "j < n" "j > m"
  shows "p A i j"
  using assms unfolding inv_from_def by auto

private lemma inv_atI[intro]: assumes " i. i < n  p A i m"
  shows "inv_at p A m"
  using assms unfolding inv_at_def by auto

private lemma inv_atD: assumes "inv_at p A m" "i < n"
  shows "p A i m"
  using assms unfolding inv_at_def by auto
  
private lemma inv_all_imp_inv_part: "m i  n  m_j  n  inv_all p A  inv_part p A m_i m_j"
  unfolding inv_all_def inv_part_def by auto

private lemma inv_all_eq_inv_part: "inv_all p A = inv_part p A n n"
  unfolding inv_all_def inv_part_def by auto

private lemma inv_part_0_Suc: "m_j < n  inv_part p A 0 m_j = inv_part p A n (Suc m_j)"
  unfolding inv_part_def by (auto, case_tac "j = m_j", auto)

private lemma inv_all_uppertD: "inv_all uppert A  j < i  i < n  A $$ (i,j) = 0"
  unfolding inv_all_def uppert_def by auto

private lemma inv_all_diff_evD: "inv_all diff_ev A  i < j  j < n 
   A $$ (i, i)  A $$ (j, j)  A $$ (i,j) = 0"
  unfolding inv_all_def diff_ev_def by auto

private lemma inv_all_diff_ev_uppertD: assumes "inv_all diff_ev A"
  "inv_all uppert A"
  "i < n" "j < n"
  and neg: "A $$ (i, i)  A $$ (j, j)"
  shows "A $$ (i,j) = 0"
proof -
  from neg have "i  j" by auto
  hence "i < j  j < i" by arith
  thus ?thesis
  proof 
    assume "i < j"
    from inv_all_diff_evD[OF assms(1) this j < n neg] show ?thesis .
  next
    assume "j < i"
    from inv_all_uppertD[OF assms(2) this i < n] show ?thesis .
  qed
qed

private lemma inv_from_bot_step: "p A i  inv_from_bot p A (Suc i)  inv_from_bot p A i"
  unfolding inv_from_bot_def by (auto, case_tac "ia = i", auto)

private lemma same_diag_refl[simp]: "same_diag A A" unfolding same_diag_def by auto
private lemma same_diag_trans: "same_diag A B  same_diag B C  same_diag A C" 
  unfolding same_diag_def by auto

private lemma same_diag_ev_blocks: "same_diag A B  ev_blocks A  ev_blocks B"
  unfolding same_diag_def ev_blocks_def ev_blocks_part_def by auto

private lemma same_uptoI[intro]: assumes " i' j'. i' < n  j' < j  A $$ (i',j') = B $$ (i',j')"
  shows "same_upto j A B"
  using assms unfolding same_upto_def by blast

private lemma same_uptoD[dest]: assumes "same_upto j A B" "i' < n" "j' < j" 
  shows "A $$ (i',j') = B $$ (i',j')"
  using assms unfolding same_upto_def by blast

private lemma same_upto_refl[simp]: "same_upto j A A" unfolding same_upto_def by auto

private lemma same_upto_trans: "same_upto j A B  same_upto j B C  same_upto j A C" 
  unfolding same_upto_def by auto

private lemma same_upto_inv_upto_jb: "same_upto j A B  inv_upto jb A j  inv_upto jb B j"
  unfolding inv_upto_def same_upto_def jb_def by auto

lemma jb_imp_diff_ev: "jb A i j  diff_ev A i j"
  unfolding jb_def diff_ev_def by auto

private lemma ev_blocks_diag: 
  "same_diag A B  ev_blocks B  ev_blocks A"
  unfolding ev_blocks_def ev_blocks_part_def same_diag_def by auto

private lemma inv_all_imp_inv_from: "inv_all p A  inv_from p A k"
  unfolding inv_all_def inv_from_def by auto

private lemma inv_all_imp_inv_at: "inv_all p A  k < n  inv_at p A k"
  unfolding inv_all_def inv_at_def by auto

private lemma inv_from_upto_at_all: 
  assumes "inv_upto jb A k" "inv_from diff_ev A k" "inv_from uppert A k" "inv_at p A k"
  and " i. i < n  p A i k  diff_ev A i k  uppert A i k"
  shows "inv_all diff_ev A" "inv_all uppert A"
proof -
  {
    fix i j
    assume ij: "i < n" "j < n"
    have "diff_ev A i j  uppert A i j"
    proof (cases "j < k")
      case True
      with assms(1) ij have "jb A i j" unfolding inv_upto_def by auto
      thus ?thesis using ij unfolding jb_def diff_ev_def uppert_def by auto
    next
      case False note ge = this
      show ?thesis
      proof (cases "j = k")
        case True
        with assms(4-) ij show ?thesis unfolding inv_at_def by auto
      next
        case False
        with ge have "j > k" by auto
        with assms(2-3) ij show ?thesis unfolding inv_from_def by auto
      qed
    qed
  }
  thus "inv_all diff_ev A" "inv_all uppert A" unfolding inv_all_def by auto
qed

private lemma lower_one_diff_uppert:
  "i < n  lower_one k B i k  diff_ev B i k  uppert B i k"
   unfolding lower_one_def diff_ev_def uppert_def by auto

definition ev_block :: "nat  'a mat  bool" where 
  " n. ev_block n A = ( i j. i < n  j < n  A $$ (i,i) = A $$ (j,j))"

lemma ev_blockD: " n. ev_block n A  i < n  j < n  A $$ (i,i) = A $$ (j,j)"
  unfolding ev_block_def carrier_mat_def by blast

lemma same_diag_ev_block: "same_diag A B  ev_block n A  ev_block n B"
  unfolding ev_block_def carrier_mat_def same_diag_def by metis


subsection ‹Alternative Characterization of @{const identify_blocks} in Presence of @{const ev_block}

private lemma identify_blocks_main_iff: assumes *: "k  k'" 
  "k'  k  k > 0  A $$ (k - 1, k)  1" and "k' < n"
  shows "set (identify_blocks_main A k list) = 
  set list  {(i,j) | i j. i  j  j < k  ( l. i  l  l < j  A $$ (l, Suc l) = 1)
   (Suc j  k'  A $$ (j, Suc j)  1)  (i > 0  A $$ (i - 1, i)  1)}" (is "_ = _  ?ss A k")
  using *
proof (induct A k list rule: identify_blocks_main.induct)
  case 1
  show ?case unfolding identify_blocks_main.simps by auto
next
  case (2 A i_e list)
  let ?s = "?ss A"
  obtain i_b where id: "identify_block A i_e = i_b" by force
  note IH = 2(1)[OF id[symmetric]]
  let ?res = "identify_blocks_main A (Suc i_e) list"  
  let ?rec = "identify_blocks_main A i_b ((i_b, i_e) # list)"
  note idb = identify_block[OF id]
  hence res: "?res = ?rec" using id by simp
  from 2(2-) have iek: "i_e < k'" by simp
  from identify_block_le'[OF id] have ibe: "i_b  i_e" .
  from ibe iek have "i_b  k'" by simp
  have "k'  i_b  0 < i_b   A $$ (i_b - 1, i_b)  1" 
    using idb(2) by auto
  note IH = IH[OF i_b  k' this]
  have cong: " a b c d. insert a c = d  set (a # b)  c = set b  d" by auto
  show ?case unfolding res IH
  proof (rule cong)
    from ibe have "?s i_b  ?s (Suc i_e)" by auto
    moreover 
    have inter: "l. i_b  l  l < i_e  A $$ (l, Suc l) = 1" using idb by blast
    have last: "Suc i_e  k'  A $$ (i_e, Suc i_e)  1" using 2(3) by auto
    have "(i_b, i_e)  ?s (Suc i_e)" using ibe idb(2) inter last by blast
    ultimately have "insert (i_b, i_e) (?s i_b)  ?s (Suc i_e)" by auto
    moreover  
    {
      fix i j
      assume ij: "(i,j)  ?s (Suc i_e)"
      hence "(i,j)  insert (i_b, i_e) (?s i_b)"
      proof (cases "j < i_b")
        case True
        with ij show ?thesis by blast
      next
        case False
        with ij have "i_b  j" "j  i_e" by auto
        {
          assume j: "j < i_e"
          from idb(3)[OF i_b  j this] have 1: "A $$ (j, Suc j) = 1" .
          from j Suc i_e  k' have "Suc j  k'" by auto
          with ij 1 have False by auto
        }
        with j  i_e have j: "j = i_e" by (cases "j = i_e", auto)
        {
          assume i: "i < i_b  i > i_b"
          hence False
          proof
            assume "i < i_b"
            hence "i_b > 0" by auto
            with idb(2) have *: "A $$ (i_b - 1, i_b)  1" by auto
            from i < i_b i_b  i_e i_e < k' have "i  i_b - 1" "i_b - 1  k'" by auto
            from i < i_b i_b  i_e j have **: "i  i_b - 1" "i_b - 1 < j" "Suc (i_b - 1) = i_b" by auto
            from ij have " l. li  l < j  A $$ (l, Suc l) = 1" by auto
            from this[OF **(1-2)] **(3) * show False by auto
          next
            assume "i > i_b"
            with ij j have "A $$ (i - 1, i)  1" and 
              *: "i - 1  i_b" "i - 1  i_e" "i - 1 < i_e" "Suc (i - 1) = i" by auto
            with idb(3)[OF *(1,3)] show False by auto
          qed
        }
        hence i: "i = i_b" by arith
        show ?thesis unfolding i j by simp
      qed
    }
    hence "?s (Suc i_e)  insert (i_b, i_e) (?s i_b)" by blast
    ultimately 
    show "insert (i_b, i_e) (?s i_b) = ?s (Suc i_e)" by blast
  qed
qed
  

private lemma identify_blocks_iff: assumes "k < n"
  shows "set (identify_blocks A k) = 
  {(i,j) | i j. i  j  j < k  ( l. i  l  l < j  A $$ (l, Suc l) = 1)
   (Suc j  k  A $$ (j, Suc j)  1)  (i > 0  A $$ (i - 1, i)  1)}" 
  unfolding identify_blocks_def using identify_blocks_main_iff[OF le_refl _ k < n] by auto

private lemma identify_blocksD: assumes "k < n" and "(i,j)  set (identify_blocks A k)"
  shows "i  j" "j < k" 
  " l. i  l  l < j  A $$ (l, Suc l) = 1"
  "Suc j  k  A $$ (j, Suc j)  1"
  "i > 0  A $$ (i - 1, i - 1)  A $$ (k,k)  A $$ (i - 1, i)  1" 
  using assms unfolding identify_blocks_iff[OF assms(1)] by auto

private lemma identify_blocksI: assumes inv: "k < n"
  "i  j" "j < k" " l. i  l  l < j  A $$ (l, Suc l) = 1"
  "Suc j  k  A $$ (j, Suc j)  1" "i > 0  A $$ (i - 1, i)  1"
  shows "(i,j)  set (identify_blocks A k)"
  unfolding identify_blocks_iff[OF inv(1)] using inv by blast

private lemma identify_blocks_rev: assumes "A $$ (i, Suc i) = 0  Suc i < k  Suc i = k"
  and inv: "k < n"
  shows "(identify_block A i, i)  set (identify_blocks A k)"
proof -
  obtain j where id: "identify_block A i = j" by force
  note idb = identify_block[OF this]
  show ?thesis unfolding id 
    by (rule identify_blocksI[OF inv], insert idb assms, auto)
qed


subsection ‹Proving the Invariants›

private lemma add_col_sub_row_diag: assumes A: "A  carrier_mat n n"
  and ut: "inv_all uppert A"
  and ijk: "i < j" "j < n" "k < n"
  shows "add_col_sub_row a i j A $$ (k,k) = A $$ (k,k)"
proof -
  from inv_all_uppertD[OF ut]
  show ?thesis
    by (subst add_col_sub_index_row, insert A ijk, auto)
qed

private lemma add_col_sub_row_diff_ev_part_old: assumes A: "A  carrier_mat n n"
  and ij: "i  j" "i  0" "i < n" "j < n" "i' < n" "j' < n"
  and choice: "j' < j  j' = j  i'  i"
  and old: "inv_part diff_ev A i j"
  and ut: "inv_all uppert A"
  shows "diff_ev (add_col_sub_row a (i - 1) j A) i' j'"
  unfolding diff_ev_def
proof (intro impI)
  assume ij': "i' < j'"
  let ?A = "add_col_sub_row a (i - 1) j A"
  assume neq: "?A $$ (i',i')  ?A $$ (j',j')"
  from A have dim: "dim_row A = n" "dim_col A = n" by auto
  note utd = inv_all_uppertD[OF ut]
  let ?i = "i - 1"
  have "?i < j" using i  j i  0 i < n by auto
  from utd[OF this j < n] have Aji: "A $$ (j,?i) = 0" by simp
  from add_col_sub_row_diag[OF A ut ?i < j j < n]
  have diag: " k. k < n  ?A $$ (k,k) = A $$ (k,k)" .
  from neq[unfolded diag[OF i' < n] diag[OF j' < n]] 
  have neq: "A $$ (i',i')  A $$ (j',j')" by auto
  {
    from inv_partD(3)[OF old i' < n j' < n choice]
    have "diff_ev A i' j'" by auto
    with neq ij' have "A $$ (i',j') = 0" unfolding diff_ev_def by auto
  } note zero = this
  {
    assume "i'  ?i" "j' = j"
    with ij' ij(1) choice have "i' > ?i" by auto
    from utd[OF this] ij
    have "A $$ (i', ?i) = 0" by auto
  } note 1 = this
  {
    assume "j'  j" "i' = ?i"
    with ij' ij(1) choice have "j > j'" by auto
    from utd[OF this] ij
    have "A $$ (j, j') = 0" by auto
  } note 2 = this
  from ij' ij choice have "(i' = ?i  j' = j) = False" by arith
  note id = add_col_sub_index_row[of i' A j' j a ?i, unfolded dim this if_False, 
    OF i' < n i' < n j' < n j' < n j < n]
  show "?A $$ (i',j') = 0" unfolding id zero using 1 2 by auto
qed

private lemma add_col_sub_row_uppert: assumes "A  carrier_mat n n"
  and "i < j"
  and "j < n"
  and inv: "inv_all uppert (A :: 'a mat)"
  shows "inv_all uppert (add_col_sub_row a i j A)"
  unfolding inv_all_def uppert_def
proof (intro allI impI)
  fix i' j'
  assume *: "i' < n" "j' < n" "j' < i'"
  note inv = inv_allD[OF inv, unfolded uppert_def]
  show "add_col_sub_row a i j A $$ (i', j') = 0"
    by (subst add_col_sub_index_row, insert assms * inv, auto)
qed

private lemma step_1_main_inv: "i  j 
   A  carrier_mat n n 
   inv_all uppert A 
   inv_part diff_ev A i j 
   inv_all uppert (step_1_main n i j A)  inv_all diff_ev (step_1_main n i j A)"
proof (induct i j A taking: n rule: step_1_main.induct)
  case (1 i j A)
  let ?i = "i - 1"
  note [simp] = step_1_main.simps[of n i j A]
  from 1(3-) have ij: "i  j" and A: "A  carrier_mat n n" and inv: "inv_all uppert A" 
    "inv_part diff_ev A i j" by auto
  show ?case
  proof (cases "j  n")
    case True
    thus ?thesis using inv by (simp add: inv_all_eq_inv_part, auto simp: inv_part_def)
  next
    case False 
    hence jn: "j < n" by simp
    note IH = 1(1-2)[OF False]
    show ?thesis
    proof (cases "i = 0")
      case True
      from inv[unfolded True inv_part_0_Suc[OF jn]]
      have inv2: "inv_part diff_ev A n (j + 1)" by simp
      have "inv_part diff_ev A (j + 1) (j + 1)" 
      proof (intro inv_partI)
        fix i' j'
        assume ij: "i' < n" "j' < n" and choice: "j' < j + 1  j' = j + 1  j + 1  i'"
        from inv_partD[OF inv2 ij] choice
        show "diff_ev A i' j'" using jn unfolding diff_ev_def by auto
      qed
      from IH(1)[OF True _ A inv(1) this]
      show ?thesis using jn by (simp, simp add: True)
    next
      case False 
      let ?evi = "A $$ (?i,?i)"
      let ?evj = "A $$ (j,j)"
      let ?choice = "?evi  ?evj  A $$ (?i, j)  0"
      let ?A = "add_col_sub_row (A $$ (?i, j) /