Linear combination of unitaries

Stream: Quantum information

[Maris Ozols]

I'm playing around with a small toy example in Lean and I'm running into an annoying issue that Matrix.unitaryGroup does not support scalar multiplication and addition. Of course, these operations do not preserve unitarity, however it would still be very natural for them to just return an object of type Matrix n n α. I know that I can always write U.1 to convert my unitary matrix U to a regular matrix, but this could just be automated away.

Would implementing smul and + for unitaryGroup but having output in Matrix go against Mathlib / Lean philosophy? Or has it not been implemented simply because nobody has run into this issue before?

[Alex Meiburg]

I think that + would definitely be against Mathlib's philosophy. ("Lean's philosophy" is a little bit more flexible and debatable, but I think it would still typically be frowned upon.)

You could define an HAdd (Matrix.unitaryGroup n α) (Matrix.unitaryGroup n α) (Matrix n n α), sure. This would mean you could then write U + V and get out a matrix. But there would be a few drawbacks:

• Now if you want to add a unitary U to an arbitrary matrix A, you'll need a different type for that. And again for A + U. These would be three different definitions. And then you would, presumably, need a bunch of lemmas to move between these.
• This addition is of heterogeneous types, so you can't say things like "this is the addition in a group". That is, you can't give an AddGroup instance for this operation, because addition gives a whole different type than the inputs. Similarly, there's no Ring or Algebra structure here either.
• Since you can't give those structures, a whole bunch of existing theorems won't apply. For instance, there's a lemma add_sub_cancel_left that would simplify A + B - A to B. But this can't apply here, since it's not an AddCommGroup. So any kind of rewriting or simplification with this notion of addition will cause additional pain.
• Similarly, tactics like abel and ring can't work with it, because they expect groups/rings respectively.

You could cirumvent some of the issues by marking your addition operation as reducible. This would mean that mean that many tactics can "see through" your definition of addition and that it's actually just defined as addition on Matrix n n α, which is a group and so better behaved, and then maybe take advantage of that. But this will still lead to some of the issues, and marking too many things as reducible can make things get very messy, because tactics are always messing with what's "inside" instead of encapsulating anything.

[Alex Meiburg]

Writing U.1 or U.val is probably the way to go. If that feels very unpleasant, another option is give type ascription: (U : Matrix n n α) will tell it that U should be coerced to a Matrix, and automatically add the ↑U that is actually just U.val.

But it can also propagate this information downwards in other ways, so if it knows that a Matrix is expected, you don't need a hint.

Here are some examples

[Alex Meiburg]

In general, Subtype.val is a bit ugly to write explicitly, sure -- but then most simplification you want can happen automatically. Like if turns out U = ⟨A, hA⟩, then ⟨A, hA⟩.val will turn into A. If you have U = V * W, then (V * W).val will turn into V.val * W.val. And so on.

[Alex Meiburg]

I _do_ think that something like this could be reasonable (as something that would type check I mean, of course not a theorem):

example (U V : Matrix.unitaryGroup (Fin 3) ℂ) (φ : unitary ℂ) : U = φ • V := by
  sorry

and currently it's unsupported. That is, if s is a unit-norm complex number, I should be able to multiply that on to a unitary matrix, to get a new unitary with a different phase. This makes sense because it's a construction I can do entirely without leaving "the unitary world". So, I could do ⟨Complex.I,⋯⟩ • U or something too.

That's currently absent from Mathlib

[Frédéric Dupuis]

The usual way to do this is to use coercions. For example, suppose you want to divide two natural numbers x y : ℕ and get a rational out, you would write (x : ℚ) / (y : ℚ) to first convert the two naturals to rational numbers before you divide (and in fact you only need to write (X : ℚ) / y, it figures out automatically that it should also convert y as well).
The same applies here, you can write (A : Matrix n n α) + (B : Matrix n n α). Then, depending on what else is around, Lean can sometimes figure it out automatically, for instance this doesn't require any type ascriptions:

import Mathlib.LinearAlgebra.UnitaryGroup

open Matrix
variable {n α : Type*} [DecidableEq n] [Fintype n] [CommRing α] [StarRing α]

example (A B : unitaryGroup n α) (C : Matrix n n α) :
    A + B = C := by sorry

[Jireh Loreaux]

Alex Meiburg said:

That's currently absent from Mathlib

Sometimes instances like this can be a bit dangerous, but I thihk the example below show that it's probably going to be okay. I think this is a reasonable instance to have for the unitary groups where you have a StarModule structure. That is,

import Mathlib

variable (R A : Type*) [CommRing R] [Ring A] [Algebra R A] [StarRing R] [StarRing A] [StarModule R A]

instance unitary.instMulAction :
    MulAction (unitary R) (unitary A) where
  smul r a := ⟨r • a, by
    rw [Submonoid.smul_def]
    simp [unitary.mem_iff]
    simp [← Submonoid.smul_def, ← unitary.coe_star, smul_smul]⟩
  one_smul _ := Subtype.ext <| one_smul ..
  mul_smul _ _ _ := Subtype.ext <| mul_smul ..

variable {R A}

lemma unitary.coe_smul (r : unitary R) (a : unitary A) : ↑(r • a) = r • (a : A) := rfl

-- perhaps this doesn't create any diamonds and life is okay
example : (unitary.instMulAction R R).toSMul = Mul.toSMul (unitary R) := by
  with_reducible_and_instances rfl

[Jireh Loreaux]

@Maris Ozols I should mention that you may also be interested in some work I have done on the unitary group in C⋆-algebras in the LeanOA repo. You can see it here. Namely, I prove that the span of the unitary elements is the entire algebra, and I prove that the unitary group is locally path connected and I provide the standard characterization of the path component of the identity (in terms of finite products of exponential unitaries). With that, you can prove that, in the matrix case, the unitary group is path connected.

[Jireh Loreaux]

Moreover, we already (almost) have the CStarAlgebra instance on matrices available as a scoped instance. To get it, you only need to add this instance:

import Mathlib

open scoped Matrix.Norms.L2Operator
variable {n : Type*} [Fintype n] [DecidableEq n]

noncomputable def Matrix.instCStarAlgebra : CStarAlgebra (Matrix n n ℂ) where

scoped[Matrix.Norms.L2Operator] attribute [instance] Matrix.instCStarAlgebra

[Maris Ozols]

Thanks everyone for your extensive answers! They helped me to understand the issue better. I think the main problem is that I'm still getting used to how strictly Lean enforces correct types... :smile:

Here is a code snippet where the issue with scalar multiplication originally came up. I'm trying to describe the structure of the single-qutrit Pauli group: for any P in the Pauli group I was hoping to write something like P = ω ^ o • (X ^ x * Z ^ z) where X and Z are Pauli matrices and ω is the 3rd root of unity.

But first I had to define UX and UZ to promote X and Z to unitaries so that I can even define the subgroup of the unitary group they generate. But then I need to demote UX and UZ back to matrices so that the scalar multiplication makes sense (or instead demote P and use X and Z on the right-hand side). Basically, I'm just wondering what is the cleanest way to do this.

-- Complex numbers
import Mathlib.Data.Complex.Basic
-- Matrices
import Mathlib.Data.Matrix.Basic
-- Complex exp and log
import Mathlib.Analysis.SpecialFunctions.Complex.Log
-- Unitary group
import Mathlib.LinearAlgebra.UnitaryGroup

open Matrix Complex Real

-- Third root of unity
noncomputable def ω : ℂ := exp (2 * π * I / 3)

-- Single-qutrit Pauli X matrix
def X : Matrix (Fin 3) (Fin 3) ℂ :=
  !![0, 0, 1;
     1, 0, 0;
     0, 1, 0]

-- Single-qutrit Pauli Z matrix
noncomputable def Z : Matrix (Fin 3) (Fin 3) ℂ :=
  !![1, 0, 0;
     0, ω, 0;
     0, 0, ω^2]

-- Unitary group U(3)
abbrev U3 := unitaryGroup (Fin 3) ℂ
abbrev M3 := Matrix (Fin 3) (Fin 3) ℂ

-- Unitary versions of X and Z
def UX : U3 := ⟨X, sorry⟩
noncomputable def UZ : U3 := ⟨Z, sorry⟩

-- Single-qutrit Pauli group
def Pauli_group : Subgroup U3 := Subgroup.closure {UX, UZ}

-- Integers modulo 3
abbrev Z3 := ZMod 3

-- Structure of qutrit Pauli group
theorem Pauli_group_structure (P : U3) (h : P ∈ Pauli_group) :
  ∃ (o x z : Z3), P = ω ^ o.val • (UX ^ x.val * UZ ^ z.val).1 := by
  sorry

theorem Pauli_group_structure' (P : U3) (h : P ∈ Pauli_group) :
  ∃ (o x z : Z3), P.1 = ω ^ o.val • X ^ x.val * Z ^ z.val := by
  sorry

Of course, I'm also annoyed that I need to write o.val etc. but that's a separate story.

[Alex Meiburg]

I see! Then, yes, I think maybe the smul instance in Jireh's message would improve this. I would then, maybe, switch the type of omega from just Complex to unitary Complex, then it would come bundled with the relevant information.

If you don't want to add the instance, you could instead by multiplying by ⟨omega * 1, ...⟩ to construct the unitary matrix from the phase.

[Alex Meiburg]

Then it will be unitaries on either side of the equality and no need for .Val's.

[Alex Meiburg]

Huh, there's nothing connections "roots of unity" to "unitary": https://loogle.lean-lang.org/?q=rootsOfUnity%2C+unitary

[Jireh Loreaux]

I'm not surprised that connection doesn't exist. Most of the library doesn't phrase things involving the unit circle in terms of unitary \C.

Per the code snippet above: you should add my instance, define omega, X, and Z all as unitary elements from the beginning, and then this will all be relatively trivial. If I have a few minutes later tonight I'll show you what I mean.

[Jireh Loreaux]

This is the kind of thing I was thinking. Potentially, you could have an unbundled X and Z, but if you're just working in the unitary group, staying bundled is nice. Note, I think ZMod 3 is a red-herring here, because your ^ operation is still only taking ℕ arguments, so ZMod 3 is just an obscure way of saying those arguments are less than 3.

import Mathlib

variable (R A : Type*) [CommRing R] [Ring A] [Algebra R A] [StarRing R] [StarRing A] [StarModule R A]

instance unitary.instMulAction :
    MulAction (unitary R) (unitary A) where
  smul r a := ⟨r • a, by
    rw [Submonoid.smul_def]
    simp [unitary.mem_iff]
    simp [← Submonoid.smul_def, ← unitary.coe_star, smul_smul]⟩
  one_smul _ := Subtype.ext <| one_smul ..
  mul_smul _ _ _ := Subtype.ext <| mul_smul ..

variable {R A}

lemma unitary.coe_smul (r : unitary R) (a : unitary A) : ↑(r • a) = r • (a : A) := rfl

-- perhaps this doesn't create any diamonds and life is okay
example : (unitary.instMulAction R R).toSMul = Mul.toSMul (unitary R) := by
  with_reducible_and_instances rfl

open Matrix Complex Real selfAdjoint

-- There's several different ways to get this, I'm not claiming this one is necessarily the best.
-- And definitely you shouldn't use `ω` as the name, but you could make it scoped notation or something.
-- Likewise for `X` and `Z`.
noncomputable def ω : unitary ℂ := expUnitary <| selfAdjointEquiv.symm (2 * π / 3)

-- Third root of unity
lemma coe_ω : ω = exp (2 * π * I / 3) := by
  simp [ω, exp_eq_exp_ℂ]
  grind

lemma orderOf_ω : orderOf ω = 3 := by
  refine orderOf_eq_prime ?_ ?_
  · sorry
  · sorry

scoped[Quantum] notation "U(" n ")" => unitaryGroup (Fin n) ℂ

open Quantum

-- Single-qutrit Pauli X matrix
@[simps]
def X : U(3) :=
  ⟨!![0, 0, 1;
     1, 0, 0;
     0, 1, 0],
    sorry⟩

lemma orderOf_X : orderOf X = 3 := by
  refine orderOf_eq_prime ?_ ?_
  · sorry
  · sorry

lemma star_X_eq_X_sq : star X = X ^ 2 := sorry

-- Single-qutrit Pauli Z matrix
@[simps]
noncomputable def Z : U(3) :=
  ⟨!![1, 0, 0;
     0, ω, 0;
     0, 0, ω^2],
    sorry⟩

lemma orderOf_Z : orderOf Z = 3 := by
  refine orderOf_eq_prime ?_ ?_
  · sorry
  · sorry

lemma star_Z_eq_Z_sq : star Z = Z ^ 2 := sorry

lemma Z_mul_X_eq_ω_smul_X_mul_Z : Z * X = ω • (X * Z) := by
  sorry

-- Single-qutrit Pauli group
def pauliGroup : Subgroup U(3) := Subgroup.closure {X, Z}

-- Structure of qutrit Pauli group
theorem pauliGroup_structure (P : U(3)) (h : P ∈ pauliGroup) :
    ∃ o x z : ℕ, o < 3 → x < 3 → z < 3 → P = ω ^ o • (X ^ x * Z ^ z) := by
  sorry

/-
In order to prove the last result, I suggest proving (which you should PR to Mathlib)
`(Subgroup.closure S).toSubmonoid = Submonoid.closure S`
whenever `S : Set G` with `Finite G`, which itself you can prove by using
`Subgroup.closure_toSubmonoid` and showing that `S⁻¹ ⊆ Submonoid.closure S`.
-/

[Jireh Loreaux]

I got bored:

import Mathlib.GroupTheory.OrderOfElement

variable {G : Type*} [Group G]

lemma Subgroup.closure_toSubmonoid_of_isOfFinOrder {s : Set G} (hs : ∀ x ∈ s, IsOfFinOrder x) :
    (closure s).toSubmonoid = Submonoid.closure s := by
  rw [closure_toSubmonoid]
  refine le_antisymm (Submonoid.closure_le.mpr ?_) (Submonoid.closure_mono <| by simp)
  refine Set.union_subset Submonoid.subset_closure ?_
  intro x hx
  have : Submonoid.powers x⁻¹ ≤ Submonoid.closure s :=
    Submonoid.powers_le.mpr (Submonoid.subset_closure hx)
  apply this
  simp [(hs _ hx).mem_powers_iff_mem_zpowers]

lemma Subgroup.closure_toSubmonoid_of_finite [Finite G] {s : Set G} :
    (closure s).toSubmonoid = Submonoid.closure s :=
  closure_toSubmonoid_of_isOfFinOrder <| by simp [isOfFinOrder_of_finite]

[Alex Meiburg]

By the way, this code looks like it has a _tiny_ bit of overlap with this file, which does a bit of the same with Z/X unitaries for qubit.

@Maris Ozols, if you want, you could try to write the appropriate Fin n version (probably using Equiv.Perm.permMatrix for X and Matrix.diagonal for Z) that generalizes both of these. That would definitely be a very welcome PR! (Even if you just want to PR the qutrit version, I'd be happy with that too)

[Alex Meiburg]

I guess the even-more-general version would be to define the generic Rz and Rx rotations, then prove that these simplify down as given for the specific angle of π. And then there's also the generalized d-qudit Hadamard which transforms between these. That would be pretty cool to write down.

(There is a slight annoyance that "4x4 Hadamard unitary" is used in literature to both refer to the d=4 qudit Fourier transform, Matrix.of fun x y => (Complex.I)^(x*y), but also the tensor product of two copies of H_2, which is a purely real matrix, and describes acting on each of a pair of qubits with a hadamard gate. But I don't think this should cause issue for us, I think the former is the obviously correct one in this setting; the latter is a different thing that just happens to also be a Hadamard matrix.)

... anyway, don't let this be a deterrent, I'm happy to start with the simplest ones first :D

[Jireh Loreaux]

#28309: induced actions between unitary groups

[Jireh Loreaux]

#28303: {Subgroup, Submonoid}.closure equalities