some linear algebra

Stream: Is there code for X?

[Alex Kontorovich]

Can someone please let me know if there is a slick way to prove this? I have an n by n matrix Q (which is really a bilinear form on an n-dimensional real vector space F), and a basis vs for F having gram matrix G, that is, the matrix of Q inner products of the vs. Then Q can be expressed as claimed below:

import linear_algebra.matrix.bilinear_form
import data.real.basic

open_locale big_operators

lemma linear_algebra_thing {F : Type*} [add_comm_group F] [module ℝ F] (Q : bilin_form ℝ F)
{n : ℕ} (hdim : finite_dimensional.finrank ℝ F = n) (vs : (fin n) → F) (G : matrix (fin n) (fin n) ℝ) (hG : invertible G)
(hvG : ∀ i j, G i j = Q (vs i) (vs j)) (w₁ w₂ : F) :
Q w₁ w₂ = ∑ i, ∑ j, Q w₁ (vs i) * Q w₂ (vs j) * (G⁻¹ i j) :=
begin
  sorry,
end

[Moritz Doll]

before looking at the code more closely (I just came home) : please avoid bilin_form, we have https://leanprover-community.github.io/mathlib_docs/linear_algebra/sesquilinear_form.html which is way better

[Moritz Doll]

in your case it would be Q : F →ₗ[ℝ] F →ₗ[ℝ] ℝ

[Moritz Doll]

my first point would be that what you wrote in the text does not match the statement, since vs is not necessarily a basis, so I am not sure whether your claim is actually correct as stated.

[Moritz Doll]

so I guess you want rather something like

import linear_algebra.sesquilinear_form
import linear_algebra.matrix.default
import data.real.basic

open_locale big_operators

variables {F : Type*} [add_comm_group F] [module ℝ F]

lemma linear_algebra_thing (Q : F →ₗ[ℝ] F →ₗ[ℝ] ℝ)
{n : ℕ}  (vs : basis (fin n) ℝ F) (G : matrix (fin n) (fin n) ℝ) (hG : invertible G)
(hvG : ∀ i j, G i j = Q (vs i) (vs j)) (w₁ w₂ : F) :
Q w₁ w₂ = ∑ i, ∑ j, Q w₁ (vs i) * Q w₂ (vs j) * (G⁻¹ i j) :=
begin
  sorry,
end

[Alex Kontorovich]

Thanks for your comments! Well, the fact that the gram matrix is invertible and the cardinality of vs is the rank of F means (to a human) that it's a basis. For my application, I'd rather not assume that it's a basis, but sure, I'd be delighted if there was a slick Lean proof of your version of the statement.

[Moritz Doll]

I think I have to rephrase your problem a little more..

[Moritz Doll]

noncomputable
def coord_bilin_form (Q : F →ₗ[ℝ] F →ₗ[ℝ] ℝ) {n : ℕ}
  (vs : basis (fin n) ℝ F) (G : matrix (fin n) (fin n) ℝ) : F →ₗ[ℝ] F →ₗ[ℝ] ℝ :=
linear_map.mk₂ ℝ (λ w₁ w₂, ∑ i, ∑ j, Q w₁ (vs i) * Q w₂ (vs j) * (G⁻¹ i j))
  sorry sorry sorry sorry

lemma coord_bilin_form_apply (Q : F →ₗ[ℝ] F →ₗ[ℝ] ℝ) {n : ℕ}
  (vs : basis (fin n) ℝ F) (G : matrix (fin n) (fin n) ℝ) (w₁ w₂ : F) :
  coord_bilin_form Q vs G w₁ w₂ = ∑ i, ∑ j, Q w₁ (vs i) * Q w₂ (vs j) * (G⁻¹ i j) := rfl

lemma linear_algebra_thing (Q : F →ₗ[ℝ] F →ₗ[ℝ] ℝ)
{n : ℕ} (vs : basis (fin n) ℝ F) (G : matrix (fin n) (fin n) ℝ) (hG : invertible G)
(hvG : ∀ i j, G i j = Q (vs i) (vs j)) :
Q = coord_bilin_form Q vs G :=
begin
  refine linear_map.ext_basis vs vs _,
  intros i j,
  rw ←hvG i j,
  rw coord_bilin_form_apply,
  simp_rw [←hvG],
  sorry,
end

So now it is at least reduced to a calculation of matrices.

[Alex Kontorovich]

Ah great, thanks! The human proof of this is completely trivial, by the way. If Q is a matrix and V is a matrix of basis vectors for F with gram matrix G=V^t*Q*V, then taking inverses and moving the V's to the other side, we get Q⁻¹ = V * G⁻¹ * V^t, which is equivalent to the claim...

[Moritz Doll]

Depending on what you want to use this for it might be easier to express everything in terms of matrices and matrix multiplication first.

[Moritz Doll]

inverting and transposing matrices works quite well in Lean

[Moritz Doll]

my idea was to go the exact opposite road and state your lemma in a rather weird invariant way

[Moritz Doll]

The proof is really not nice:

import linear_algebra.sesquilinear_form
import linear_algebra.matrix.default
import data.real.basic

open_locale big_operators

variables {F : Type*} [add_comm_group F] [module ℝ F]

noncomputable
def coord_bilin_form (Q : F →ₗ[ℝ] F →ₗ[ℝ] ℝ) {n : ℕ}
  (vs : basis (fin n) ℝ F) (G : matrix (fin n) (fin n) ℝ) : F →ₗ[ℝ] F →ₗ[ℝ] ℝ :=
linear_map.mk₂ ℝ (λ w₁ w₂, ∑ i, ∑ j, Q w₁ (vs i) * Q w₂ (vs j) * (G⁻¹ i j))
  sorry sorry sorry sorry

lemma coord_bilin_form_apply (Q : F →ₗ[ℝ] F →ₗ[ℝ] ℝ) {n : ℕ}
  (vs : basis (fin n) ℝ F) (G : matrix (fin n) (fin n) ℝ) (w₁ w₂ : F) :
  coord_bilin_form Q vs G w₁ w₂ = ∑ i, ∑ j, Q w₁ (vs i) * Q w₂ (vs j) * (G⁻¹ i j) := rfl

lemma linear_algebra_thing (Q : F →ₗ[ℝ] F →ₗ[ℝ] ℝ)
{n : ℕ} (vs : basis (fin n) ℝ F) (G : matrix (fin n) (fin n) ℝ) (hG : invertible G)
(hvG : ∀ i j, G i j = Q (vs i) (vs j)) :
Q.flip = coord_bilin_form Q vs G :=
begin
  refine linear_map.ext_basis vs vs _,
  intros i j,
  rw linear_map.flip_apply,
  rw ←hvG j i,
  rw coord_bilin_form_apply,
  simp_rw [←hvG],
  simp_rw mul_comm (G i _),
  simp_rw mul_assoc,
  rw finset.sum_comm,
  simp_rw ←finset.mul_sum,
  simp_rw ←matrix.mul_apply,
  simp_rw matrix.mul_inv_of_invertible,
  simp_rw matrix.one_apply,
  simp,
end

[Moritz Doll]

it would be ok if there is a lemma that calculates (A.mul (B.mul C)) x y as the finsum, but I haven't found anything in that direction

[Moritz Doll]

but note the flip, it seems that I forgot one transposition when first rewriting your lemma

[Alex Kontorovich]

Great, thanks!

[Heather Macbeth]

Moritz Doll said:

def coord_bilin_form (Q : F →ₗ[ℝ] F →ₗ[ℝ] ℝ) {n : ℕ}
  (vs : basis (fin n) ℝ F) (G : matrix (fin n) (fin n) ℝ) : F →ₗ[ℝ] F →ₗ[ℝ] ℝ :=
linear_map.mk₂ ℝ (λ w₁ w₂, ∑ i, ∑ j, Q w₁ (vs i) * Q w₂ (vs j) * (G⁻¹ i j))
  sorry sorry sorry sorry

It is rarely necessary to prove linearity conditions like this, there is usually some abstract nonsense. In this case you could first construct each individual function

λ w₁ w₂, Q w₁ (vs i) * Q w₂ (vs j)

as a bilinear map using docs#linear_map.compl₂ and docs#linear_map.mul (if we don't have that construction already) and then add the appropriate scalar multiples of these together.

[Moritz Doll]

This is a very good point, thanks Heather. I was just lazy and wanted to get the bundled version as fast as possible.

[Moritz Doll]

There might be also a point to commute the factors (and the summations) in the definition so that the proof gets less annoying

[Moritz Doll]

but depending on how you replace the naive definition with a bundled one this might actually happen automatically

[Alex Kontorovich]

May I ask another stupid question in a similar vein? How would one prove something like this? I have an n by n matrix with -1s on the diagonal and 1s off the diagonal. I claim that its inverse matrix has something on the diagonal and something else off the diagonal. What's an efficient (?) way to verify this with a calculation? Thanks!

import data.matrix.basic
import data.real.basic

example {n : ℕ} (hn : n ≠ 0) (G : matrix (fin n) (fin n) ℝ) (hG : ∀ i j, G i j = ite (i=j) (-1) 1)
  (Ginv : matrix (fin n) (fin n) ℝ) (hGinv : ∀ i j, Ginv i j = ite (i = j) ((1-n)/(2*n)) (1/(2*n)))
  : Ginv.mul G = 1 :=
begin
  simp only [matrix.mul, matrix.dot_product, hGinv, hG],
  sorry,
end

[Yakov Pechersky]

Some progress:

import data.fin.basic

import data.matrix.basic
import data.real.basic

@[simp] lemma finset.filter_eq_univ {α : Type*} [fintype α] [decidable_eq α] (i : α) :
  finset.univ.filter (eq i) = {i} :=
by ext; simp [eq_comm]

@[simp] lemma finset.filter_eq_univ' {α : Type*} [fintype α] [decidable_eq α] (i : α) :
  finset.univ.filter (λ x, x = i) = {i} :=
by ext; simp [eq_comm]

@[simp] lemma finset.filter_ne_univ {α : Type*} [fintype α] [decidable_eq α] (i : α) :
  finset.univ.filter (λ x, x ≠ i) = finset.univ \ {i} :=
by ext; simp [eq_comm]

@[simp] lemma finset.filter_ne_univ' {α : Type*} [fintype α] [decidable_eq α] (i : α) :
  finset.univ.filter (λ x, i ≠ x) = finset.univ \ {i} :=
by ext; simp [eq_comm]

example {n : ℕ} (hn : n ≠ 0) (G : matrix (fin n) (fin n) ℝ) (hG : ∀ i j, G i j = ite (i=j) (-1) 1)
  (Ginv : matrix (fin n) (fin n) ℝ) (hGinv : ∀ i j, Ginv i j = ite (i = j) ((1-n)/(2*n)) (1/(2*n)))
  : Ginv.mul G = 1 :=
begin
  simp only [matrix.mul, matrix.dot_product, hGinv, hG],
  simp [finset.sum_ite],
end

[Junyan Xu]

I don't think this is correct? The matrix product has (n-1)/n on diagonal and -1/n off diagonal according to my calculation. And when n=2, the matrix G has rank 1 so it doesn't have an inverse.

[Yakov Pechersky]

Here's "proof" of it being incorrect, given how it's phrased right now:

import data.fin.basic

import data.matrix.basic
import data.real.basic

@[simp] lemma finset.filter_eq_univ {α : Type*} [fintype α] [decidable_eq α] (i : α) :
  finset.univ.filter (eq i) = {i} :=
by ext; simp [eq_comm]

@[simp] lemma finset.filter_eq_univ' {α : Type*} [fintype α] [decidable_eq α] (i : α) :
  finset.univ.filter (λ x, x = i) = {i} :=
by ext; simp [eq_comm]

@[simp] lemma finset.filter_ne_univ {α : Type*} [fintype α] [decidable_eq α] (i : α) :
  finset.univ.filter (λ x, x ≠ i) = finset.univ \ {i} :=
by ext; simp [eq_comm]

@[simp] lemma finset.filter_ne_univ' {α : Type*} [fintype α] [decidable_eq α] (i : α) :
  finset.univ.filter (λ x, i ≠ x) = finset.univ \ {i} :=
by ext; simp [eq_comm]

@[simp] lemma finset.filter_eq_singleton {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (eq i) {k} = if i = k then {i} else ∅ :=
by split_ifs; ext; simp [h] {contextual := tt}

@[simp] lemma finset.filter_eq_singleton' {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (λ x, x = i) {k} = if i = k then {i} else ∅ :=
by split_ifs; ext; simp [h, eq_comm] {contextual := tt}

@[simp] lemma finset.filter_ne_singleton {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (ne i) {k} = if i = k then ∅ else {k} :=
by split_ifs; ext; simp [h] {contextual := tt}

@[simp] lemma finset.filter_ne_singleton' {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (λ x, x ≠ i) {k} = if i = k then ∅ else {k} :=
by split_ifs; ext; simp [h, eq_comm] {contextual := tt}

@[simp] lemma finset.filter_eq_univ_erase {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (eq i) (finset.univ \ {k}) = if i = k then ∅ else {i} :=
by { split_ifs; ext; simp [h, eq_comm] {contextual := tt}; rintro rfl; simp [ne.symm h]  }

@[simp] lemma finset.filter_eq_univ_erase' {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (λ x, x =  i) (finset.univ \ {k}) = if i = k then ∅ else {i} :=
by split_ifs; ext; simp [h, eq_comm] {contextual := tt}

@[simp] lemma finset.filter_ne_univ_erase {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (ne i) (finset.univ \ {k}) = finset.univ \ {i, k} :=
by ext; simp [not_or_distrib, eq_comm, and.comm]

@[simp] lemma finset.filter_ne_univ_erase' {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (λ x, x ≠ i) (finset.univ \ {k}) = finset.univ \ {i, k} :=
by ext; simp [not_or_distrib, eq_comm, and.comm]

example {n : ℕ} (hn : n ≠ 0) (G : matrix (fin n) (fin n) ℝ) (hG : ∀ i j, G i j = ite (i=j) (-1) 1)
  (Ginv : matrix (fin n) (fin n) ℝ) (hGinv : ∀ i j, Ginv i j = ite (i = j) ((1-n)/(2*n)) (1/(2*n)))
  : Ginv.mul G = 1 :=
begin
  simp only [matrix.mul, matrix.dot_product, hGinv, hG],
  ext i k,
  cases eq_or_ne k i with h h,
  { simp [finset.sum_ite, h, matrix.one_apply],
    field_simp,
    ring_nf,
    simp [mul_comm, mul_assoc, mul_left_comm, hn],
    -- false
    sorry },
  { simp [finset.sum_ite, h, matrix.one_apply, ne.symm h],
    field_simp [hn],
    ring_nf,
    simp [hn],
    -- false
    sorry }
end

[Alex Kontorovich]

That's great, thanks! I think I can take it from here...

[Alex Kontorovich]

Here's something very bizarre that's happening. This works perfectly well:

import data.fin.basic
import data.matrix.basic
import data.real.basic

@[simp] lemma finset.filter_eq_univ {α : Type*} [fintype α] [decidable_eq α] (i : α) :
  finset.univ.filter (eq i) = {i} :=
by ext; simp [eq_comm]

@[simp] lemma finset.filter_eq_univ' {α : Type*} [fintype α] [decidable_eq α] (i : α) :
  finset.univ.filter (λ x, x = i) = {i} :=
by ext; simp [eq_comm]

@[simp] lemma finset.filter_ne_univ {α : Type*} [fintype α] [decidable_eq α] (i : α) :
  finset.univ.filter (λ x, x ≠ i) = finset.univ \ {i} :=
by ext; simp [eq_comm]

@[simp] lemma finset.filter_ne_univ' {α : Type*} [fintype α] [decidable_eq α] (i : α) :
  finset.univ.filter (λ x, i ≠ x) = finset.univ \ {i} :=
by ext; simp [eq_comm]

@[simp] lemma finset.filter_eq_singleton {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (eq i) {k} = if i = k then {i} else ∅ :=
by split_ifs; ext; simp [h] {contextual := tt}

@[simp] lemma finset.filter_eq_singleton' {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (λ x, x = i) {k} = if i = k then {i} else ∅ :=
by split_ifs; ext; simp [h, eq_comm] {contextual := tt}

@[simp] lemma finset.filter_ne_singleton {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (ne i) {k} = if i = k then ∅ else {k} :=
by split_ifs; ext; simp [h] {contextual := tt}

@[simp] lemma finset.filter_ne_singleton' {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (λ x, x ≠ i) {k} = if i = k then ∅ else {k} :=
by split_ifs; ext; simp [h, eq_comm] {contextual := tt}

@[simp] lemma finset.filter_eq_univ_erase {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (eq i) (finset.univ \ {k}) = if i = k then ∅ else {i} :=
by { split_ifs; ext; simp [h, eq_comm] {contextual := tt}; rintro rfl; simp [ne.symm h]  }

@[simp] lemma finset.filter_eq_univ_erase' {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (λ x, x =  i) (finset.univ \ {k}) = if i = k then ∅ else {i} :=
by split_ifs; ext; simp [h, eq_comm] {contextual := tt}

@[simp] lemma finset.filter_ne_univ_erase {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (ne i) (finset.univ \ {k}) = finset.univ \ {i, k} :=
by ext; simp [not_or_distrib, eq_comm, and.comm]

@[simp] lemma finset.filter_ne_univ_erase' {α : Type*} [fintype α] [decidable_eq α] (i k : α) :
  finset.filter (λ x, x ≠ i) (finset.univ \ {k}) = finset.univ \ {i, k} :=
by ext; simp [not_or_distrib, eq_comm, and.comm]

example {n : ℕ} (hn : (n:ℝ) ≥ 3) (G : matrix (fin n) (fin n) ℝ)
  (hG : ∀ i j, G i j = ite (i=j) (-1) 1) (Ginv : matrix (fin n) (fin n) ℝ)
  (hGinv : ∀ i j, Ginv i j = if (i = j) then (-((-3 + n) / (2 * (-2 + n))))
    else (1 / (2 * (-2 + n)))) :
  Ginv.mul G = 1 :=
begin
  have : -2 + (n : ℝ) ≠ 0 := by linarith,
  ext i j,
  by_cases hh : i = j,
  { simp only [matrix.mul, matrix.dot_product, hGinv, hG, hh, one_div, mul_inv_rev, mul_ite,
      mul_neg, mul_one, matrix.one_apply_eq],
    simp [finset.sum_ite], -- squeeze_simp doesn't give full simp!!!???
    field_simp,
    ring, },
  { simp only [matrix.mul, matrix.dot_product, hGinv, hG, hh, one_div, mul_inv_rev, mul_ite,
      mul_neg, mul_one, matrix.one_apply_ne, ne.def, not_false_iff],
    simp [finset.sum_ite, hh],
    field_simp,
    ring, },
end

But if I replace those inner simps with squeeze_simps, then it tells me to use this long line:

 simp only [finset.sum_ite, finset.filter_eq_univ', finset.sum_neg_distrib, finset.filter_eq_singleton, eq_self_iff_true, if_true,
  finset.sum_const, finset.card_singleton, nsmul_eq_mul, nat.cast_one, one_mul, finset.filter_ne_singleton,
  finset.card_empty, zero_smul, add_zero, neg_neg, finset.filter_ne_univ, finset.filter_eq_univ_erase, neg_zero',
  finset.filter_ne_univ_erase, finset.pair_self_eq, finset.sum_sdiff_eq_sub, finset.subset_univ, finset.card_fin,
  zero_add],

which does not simplify to the same amount?!?!?

Thanks again for all your help!!

[Yakov Pechersky]

In cases like this, I extract the goal after the simp that works, clean it of coercion arrows, and make that goal into a suffices

[Eric Wieser]

Note we already have docs#finset.filter_ne and docs#finset.filter_eq