Associativity of kronecker product

Stream: new members

[Yunong Shi]

Hi, everyone. @Runzhou Tao and I are new to Lean. We want to use Lean in quantum computing. As a practice, we want to prove the associativity of tensor products of (finite dimensional) matrices in different dimensions using the following code. It seems like defining kronecker product as on line 14, it depends on the dimension of the matrices, thus we cannot express the associativity theorem easily because the 3rd matrix will have a different dimension (a failed trial is on line 26). I am wondering is there a way of defining kronecker product (and matrices) so that the associativity can be expressed and proved? Thanks!

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

def comp {m n} : fin m → fin n → fin (m * n)
| ⟨a, ha⟩ ⟨b, hb⟩ := ⟨a * n + b, by nlinarith⟩

variables (p q : ℕ)
variables (x : fin p) (y : fin q)


@[simp] def my_matrix (m n : ℕ) := matrix (fin m) (fin n) ℝ

def kron {m1 n1 m2 n2} (a : my_matrix m1 n1) (b : my_matrix m2 n2) : my_matrix (m1*m2) (n1*n2) :=
begin
    intros i j,
    cases i with i hi,
    cases j with j hj,
    let i1 := i % m1,
    let j1 := j % n1,
    let i2 := i / m1,
    let j2 := j / n1,
    exact 0, -- m1 n1 m2 n2 should not be
end

theorem kron_assoc {m1 n1 m2 n2 m3 n3} (a : my_matrix m1 n1) (b : my_matrix m2 n2) (c : my_matrix m3 n3) :
kron (kron a b) c = kron a (kron b c) := sorry

[Kevin Buzzard]

Yup, this is a standard problem with dependent type theory and there are a few ways around it, although I still haven't really internalised them. The problem is that x*(y*z) and (x*y)*z are equal, but not definitionally equal.

[Kevin Buzzard]

This is the ugly way:

def my_matrix_equiv {m₁ n₁ m₂ n₂ : ℕ} (hm : m₁ = m₂) (hn : n₁ = n₂)
  (a : my_matrix m₁ n₁) (b : my_matrix m₂ n₂) : Prop :=
∀ (i : fin m₁) (j : fin n₁), a i j = b ⟨i.1, hm ▸ i.2⟩ ⟨j.1, hn ▸ j.2⟩

theorem kron_assoc {m1 n1 m2 n2 m3 n3} (a : my_matrix m1 n1) (b : my_matrix m2 n2) (c : my_matrix m3 n3) :
my_matrix_equiv (mul_assoc _ _ _) (mul_assoc _ _ _) (kron (kron a b) c) (kron a (kron b c)) := sorry

but I think there are better ways nowadays.

[Eric Wieser]

Note that comp a b exists already as docs#fin_prod_fin_equiv, fin_prod_fin_equiv (a, b)

[Eric Wieser]

You can also use heq, which isn't always as bad as its made out to be, kron (kron a b) c == kron a (kron b c)

[Yakov Pechersky]

import data.matrix.notation

section repr

variables {m' n' : ℕ} {X : Type*} [has_repr X]

def vec_repr : Π {n' : ℕ}, (fin n' → X) → string :=
λ _ v, string.intercalate ", " ((vector.of_fn v).to_list.map repr)

instance vec_repr_instance : has_repr (fin n' → X) := ⟨vec_repr⟩

def matrix_repr : Π {m' n'}, matrix (fin m') (fin n') X → string :=
λ _ _ M, string.intercalate ";\n" ((vector.of_fn M).to_list.map repr)

instance matrix_repr_instance :
  has_repr (matrix (fin n') (fin m') X) := ⟨matrix_repr⟩

end repr

section kron

variables {K : Type*} [semiring K]

def kron {m1 n1 m2 n2 : ℕ}
  (a : matrix (fin (m1 + 1)) (fin (n1 + 1)) K) (b : matrix (fin (m2 + 1)) (fin (n2 + 1)) K) :
  matrix (fin ((m1 + 1) * (m2 + 1))) (fin ((n1 + 1) * (n2 + 1))) K :=
λ i j, ((a ⟨i / (m2 + 1), (nat.div_lt_iff_lt_mul' (nat.succ_pos')).mpr i.is_lt⟩
           ⟨j / (n2 + 1), (nat.div_lt_iff_lt_mul' (nat.succ_pos')).mpr j.is_lt⟩)
           • b)
           ⟨i % (m2 + 1), nat.mod_lt _ nat.succ_pos'⟩
           ⟨j % (n2 + 1), nat.mod_lt _ nat.succ_pos'⟩


#eval kron ![![1, 2], ![3, 4]] ![![0, 5], ![6, 7]]
/-
0, 5, 0, 10;
6, 7, 12, 14;
0, 15, 0, 20;
18, 21, 24, 28
-/

theorem kron_assoc
{m1 n1 m2 n2 m3 n3 : ℕ} {K : Type*} [semiring K]
  (a : matrix (fin (m1 + 1)) (fin (n1 + 1)) K)
  (b : matrix (fin (m2 + 1)) (fin (n2 + 1)) K)
  (c : matrix (fin (m3 + 1)) (fin (n3 + 1)) K)
  :
kron (@kron _ _ _ _ _ _ a b) c == (@kron _ _ _ _ _ _ a (@kron _ _ _ _ _ _ b c)) := sorry

[Yakov Pechersky]

Although proving associativity might be a tiny bit simpler if the mul definition is used instead of smul

[Frédéric Dupuis]

Hi and welcome! I also work in quantum information, and I've thought a bit about the best way to represent states on multiple quantum systems in Lean. I think the way to go will be to write a pi-type version of the code in linear_algebra/tensor_product.lean. What we have there is the tensor product of two spaces M and N, and it's defined as a free_add_monoid M × N quotiented by an appropriate equivalence relation, so basically arbitrary sums of vectors of the form x ⊗ y and then selecting the right equivalence classes. I think we should do the same thing with free_add_monoid (Π i : ι, f i) with f : ι → Type* in order to have more than two systems. This would have the advantage of avoiding these associativity problems, by not even having an order on the systems in the first place.

[Frédéric Dupuis]

It also lets us work at a more abstract level, which I've (somewhat counterintuitively) found to be a lot easier in Lean. I also tried to go this route when I first learned Lean (i.e. do everything with concrete matrices) and it was constant pain, with indices everywhere.

[Yunong Shi]

Thanks everyone. It will take some time for me to fully understand your replies. Let me get back to some of you later.

[Yunong Shi]

@Frédéric Dupuis Thanks for the great suggestions and glad to know someone here is in quantum computation as well! As I understand your comment, there are two suggestions: 1. we can work with abstract tensor product that defined in linear_algebra/tensor_product.lean. 2. Using Pi-type to get rid of the dimension dependence. For 1, will the abstract definition make some proof harder, for example, prove some concrete results are indeed tensor products? For 2, I think it's a great idea! Will it involve some changes in the matrix definition (matrices with arbitrary dimensions)? Maybe I should look into the implementation more.

[Alexander Bentkamp]

@Frédéric Dupuis @Yunong Shi Not sure if it's useful for your application, but I've formalized an alternative definition of tensors in here:
https://github.com/leanprover-community/mathlib/blob/1f309c52c472a6cd37e41c580a470033e63d4720/src/data/holor.lean

We called this alternative view on tensors "holors" (I think we found that word in the literature somewhere). The difference is that this definition views the tensors as multidimensional arrays, not as functions. It's like viewing linear maps as matrices. The holor file also contains a proof of associativity of the tensor product.

[Eric Wieser]

@Alexander Bentkamp, is there a reason you use \N and not fin n like matrices (usually) use?

[Alexander Bentkamp]

I guess it makes the whole definition a bit simpler, but I am not sure. In some sense, the definition of holor_index is like fin n, but on lists.

[Alexander Bentkamp]

How would the definition of holor_index look like with fin n? I guess you would need a Pi-type or something?

[Alexander Bentkamp]

For matrices, it's simpler. It's just fin n ⨉ fin m.

[Eric Wieser]

def holor {n : ℕ} (α : Type*) (ds : fin n → Type*) [∀ i, fintype (ds i)] := (Π i, ds i) → α

[Eric Wieser]

Eg holor ℝ ![fin 2, fin 3, fin 4] for a 2x3x4 tensor

[Alexander Bentkamp]

wow, what is the ![...] notation? It constructs fin n?

[Eric Wieser]

https://leanprover-community.github.io/mathlib_docs/data/matrix/notation.html

[Alexander Bentkamp]

Well, then the reason for my definition was that the fin n/matrix library wasn't as developed back then.

[Eric Wieser]

It constructs fin n → T

[Eric Wieser]

Ugh, that definition doesn't work: (https://leanprover.zulipchat.com/#narrow/stream/113488-general/topic/fintype.20instances.20for.20matrix.20notation/near/219314820)

[Eric Wieser]

On second thoughts, my definition probably doesn't work because it means 2-tensors and 3-tensors have different types

[Eric Wieser]

Maybe I'll explore trying to overhaul holor to use it and see where I get stuck

[Alexander Bentkamp]

Sure, go for it. It would be nice to have it more consistent with the matrix library.

[Eric Wieser]

Having explored a little more, I think all the holor_index stuff is still quite hard to write with fin. #4406 would make it easier.

[Frédéric Dupuis]

@Alexander Bentkamp Thanks for pointing this out! This is very close to what I want, but I think the definition I outlined above would still be worth doing, for a few reasons. I would really like to avoid having the tensor factors indexed by a specific predefined type like ℕ or fin n: for example, in a typical application in quantum information, we would have a quantum state $|φ⟩_{A^n B^n E}$ with $A^n$ and $B^n$ being two $n$-qubit states, and $E$ being some other quantum system. So this object would be a vector in a space that is the tensor product of $2n+1$ spaces with those labels, and I'd like to be able to keep those system labels and never have to artificially number them. Also, one would want to be able to support infinite-dimensional spaces, which wouldn't be possible with holors.

[Frédéric Dupuis]

@Yunong Shi I don't know what your plans are, but here's my take on doing quantum information in Lean at the moment: we're not there yet! :-) To give you an idea of the current situation, unitaries are not there, adjoints are not there (so forget Hermitian operators), there is no spectral decomposition, and we only got complex inner products and eigenvalues about two months ago. If you have a specific application that doesn't require too many of those features it might be possible to get started now, but personally my plan is to try to get those into mathlib by making PRs whenever I have a bit of time before trying to do actual quantum information.

[Frédéric Dupuis]

Depending on what you're doing, one option might be to switch to Coq and use this library: https://github.com/inQWIRE/QWIRE

[Yunong Shi]

Frédéric Dupuis said:

Yunong Shi I don't know what your plans are, but here's my take on doing quantum information in Lean at the moment: we're not there yet! :-) To give you an idea of the current situation, unitaries are not there, adjoints are not there (so forget Hermitian operators), there is no spectral decomposition, and we only got complex inner products and eigenvalues about two months ago. If you have a specific application that doesn't require too many of those features it might be possible to get started now, but personally my plan is to try to get those into mathlib by making PRs whenever I have a bit of time before trying to do actual quantum information.

Thanks for the info. I'd love to help with the dev of libraries for quantum information. I will probably look at the code you mentioned and play around with it first.

[Frédéric Dupuis]

I just managed to motivate myself to have a crack at this idea for the tensor product of an indexed family of spaces: #5311

[Eric Wieser]

What should the tensor product of an empty index be? For the sake of tensor algebras, it feels like it should be isomorphic to R, but I think your PR is isomorphic to the trivial ring.

[Kevin Buzzard]

$M^{\otimes 0}=R$ but I'm not sure if that's what you're asking

[Riccardo Brasca]

If you mean $\bigotimes_{i \in \emptyset} M_i$, this should definitely be $R$, as in Kevin's answer.

[Riccardo Brasca]

Having it the zero ring would break the usual formulas for the union of set indexes.

[Eric Wieser]

That's mostly what I'm asking, but for an indexed tensor product instead of a homogenous power - edit: yes

[Frédéric Dupuis]

My math answer would also be R, but right now the code assumes that the index type is nonempty.

[Frédéric Dupuis]

It's a bit awkward to change, because the has_scalar instance is defined by choosing one of the indices and doing smul on that coordinate.

[Frédéric Dupuis]

Also right now we can only do the tensor product over the whole type, and not a subset. It might be useful to put that in actually.

[Eric Wieser]

Storing an R and Pi i, M i term separately under your quotient would fix that base case

[Frédéric Dupuis]

Yeah, I was hoping for a more clever way but it doesn't look possible.