two ways to define inverse of a diagonal matrix

Stream: new members

[Z. X. Wang]

import Mathlib

open Matrix

abbrev Matrix.diagPart {n α : Type*} [DecidableEq n] (X : Matrix n n α) [Zero α] :=
  Matrix.diagonal (fun i => X i i)

-- def D⁻¹ by entries
noncomputable def invdiagPart
  (A : Matrix (Fin n) (Fin n) ℝ)
  :
  Matrix (Fin n) (Fin n) ℝ
  :=
  fun i j ↦ if i = j then 1 / (A i j) else 0

example
  (A : Matrix (Fin n) (Fin n) ℝ)
  (hdiag : ∀ i, A i i > 0)
  :
  (diagPart A)⁻¹ = invdiagPart A
  := sorry

how can i prove that they are equal? And here is my atempt:

import Mathlib

open Matrix

abbrev Matrix.diagPart {n α : Type*} [DecidableEq n] (X : Matrix n n α) [Zero α] :=
  Matrix.diagonal (fun i => X i i)

-- def D⁻¹ by entries
noncomputable def invdiagPart
  (A : Matrix (Fin n) (Fin n) ℝ)
  :
  Matrix (Fin n) (Fin n) ℝ
  :=
  fun i j ↦ if i = j then 1 / (A i j) else 0

example
  (A : Matrix (Fin n) (Fin n) ℝ)
  (hdiag : ∀ i, A i i > 0)
  :
  (diagPart A)⁻¹ = invdiagPart A
  := by
  refine ext ?_
  intro i j
  unfold diagPart invdiagPart diagonal
  simp only [one_div]
  by_cases h1 : i = j
  . split_ifs
    rw [inv_def]
    simp only [Ring.inverse_eq_inv', smul_apply, smul_eq_mul]
    refine (IsUnit.inv_mul_eq_iff_eq_mul ?_).mpr ?_
    . have h2 : (diagPart A).det ≠ 0 := by
        unfold diagPart
        simp only [det_diagonal, ne_eq]
        by_contra a
        have h3 : ∏ i, A i i > 0 := by
          exact Finset.prod_pos fun i a ↦ hdiag i
        have h4 : ∏ i, A i i ≠ 0 := by
          exact Ne.symm (ne_of_lt h3)
        contradiction
      exact Ne.isUnit h2
    .
      sorry
  . split_ifs
    sorry

[Kevin Buzzard]

You seem to have stumbled upon a hole in the library (the sorry in the below code):

import Mathlib

set_option autoImplicit true

open Matrix

abbrev Matrix.diagPart {n α : Type*} [DecidableEq n] (X : Matrix n n α) [Zero α] :=
  Matrix.diagonal (fun i => X i i)

-- def D⁻¹ by entries
noncomputable def invdiagPart
  (A : Matrix (Fin n) (Fin n) ℝ)
  :
  Matrix (Fin n) (Fin n) ℝ
  :=
  fun i j ↦ if i = j then 1 / (A i j) else 0

lemma Matrix.diagonal_inv {n α : Type*} [DecidableEq n] [Fintype n] [Field α]
    {x : n → α} (hx : ∀ i, x i ≠ 0) : (diagonal x)⁻¹ = diagonal (fun i ↦ (x i)⁻¹) :=
  sorry

example
  (A : Matrix (Fin n) (Fin n) ℝ)
  (hdiag : ∀ i, A i i > 0)
  :
  (diagPart A)⁻¹ = invdiagPart A
  := by
ext i j
unfold diagPart invdiagPart
rw [Matrix.diagonal_inv]
· simp [diagonal_apply]
  grind
· grind

[Eric Wieser]

docs#Matrix.inv_diagonal exists

[Eric Wieser]

Note that invdiagPart is ill-typed, as you forgot the Matrix.of