negative (semi-)definite matrices

Stream: mathlib4

[Monica Omar]

Should we add the definitions for negative (semi-)definite matrices? It's defined in wikipedia.
I'm hesitant since it's literally just the negative of a positive (semi-)definite matrix (as in, x is negative (semi-)definite iff -x is positive (semi-)definite).

[Eric Wieser]

I think in the past @Jireh Loreaux and I discussed replacing docs#PosSemidef with 0 \le A via a StarOrderedRing instance

[Monica Omar]

that would be less general than docs#Matrix.PosSemidef right now though

[Monica Omar]

what we can do is add a partial order and StarOrderedRing instance now via docs#Matrix.PosSemidef (for real and complex)

[Eric Wieser]

Monica Omar said:

that would be less general than docs#Matrix.PosSemidef right now though

How so?

[Eric Wieser]

I can't find the old thread, but my naive proposal is

instance : LE (Matrix n n R) where
  le M N := (M - N).IsHermitian ∧ ∀ x : n → R, star x ⬝ᵥ (M *ᵥ x) ≤ star x ⬝ᵥ (N *ᵥ x)

though obviously the M - N is problematic for semirings

[Monica Omar]

It would define a LE but not a partial order for general *-rings

[Eric Wieser]

But for StarOrderedRing R would it result in StarOrderedRing (Matrix n n R)?

[Eric Wieser]

Eric Wieser said:

I can't find the old thread

#mathlib4 > Conflicting `≤` and `<` @ 💬

[Monica Omar]

Eric Wieser said:

But for StarOrderedRing R would it result in StarOrderedRing (Matrix n n R)?

I think we'd need a C*-algebra or some sort of norm structure to turn Matrix n n R into a StarOrderedRing (as far as I know, correct me if i'm wrong)

[Monica Omar]

^ in which case we have docs#CStarMatrix.instStarOrderedRing

[Jireh Loreaux]

1. No, we should not define negative semidefinite matrices too.
2. docs#CStarMatrix.instStarOrderedRing is just defined using the spectral order, but that's suboptimal in the long run.
3. We should provide the instance Eric mentions, at least for RCLike entries.
4. You can provide a StarOrderedRing instance on that, even without invoking the C⋆-norm. There's two ways to do this: via the pre-existing square root machinery, or using the continuous functional calculus for hermitian matrices (which we have thanks to @Jon Bannon).
5. At some point, I did this in a branch, but maybe only locally. Let me look for it later today.

[Jireh Loreaux]

Indeed, I claimed to have built that instance before: #PR reviews > #8809 positive semidef matrices @ 💬

[Monica Omar]

do you mean adding an instance StarOrderedRing (Matrix n n 𝕜)?

[Jireh Loreaux]

oh right, that's docs#CStarMatrix

[Jireh Loreaux]

Yes, I meant one for Matrix.

[Monica Omar]

I think using the pre-existing square root machinery would be fairly straightforward:

import Mathlib.LinearAlgebra.Matrix.PosDef

namespace Matrix
variable {n 𝕜} [RCLike 𝕜] [Fintype n] [DecidableEq n]

open ComplexOrder

instance : PartialOrder (Matrix n n 𝕜) where
  le A B := (B - A).PosSemidef
  le_refl A := sub_self A ▸ .zero
  le_trans _ _ _ h₁ h₂ := by simpa using h₁.add h₂
  le_antisymm _ _ h₁ h₂ := by
    rw [← sub_eq_zero, ← toLin'.injective.eq_iff, LinearMap.ext_iff]
    intro a
    simp only [toLin'_apply, zero_mulVec, ← h₂.dotProduct_mulVec_zero_iff]
    have := h₁.2 a
    rw [← neg_sub, Matrix.neg_mulVec, dotProduct_neg, le_neg, neg_zero] at this
    exact le_antisymm this (h₂.2 a)

theorem le_def {x y : Matrix n n 𝕜} : x ≤ y ↔ (y - x).PosSemidef := Iff.rfl

theorem nonneg_iff_posSemidef {x : Matrix n n 𝕜} : 0 ≤ x ↔ x.PosSemidef := by simp [le_def]

instance : StarOrderedRing (Matrix n n 𝕜) := .of_le_iff <| fun a b => by
  refine ⟨fun hab => ?_, fun ⟨s, hs⟩ => ?_⟩
  · simp_rw [← sub_eq_iff_eq_add']
    exact posSemidef_iff_eq_conjTranspose_mul_self.mp hab
  · simp_rw [le_def, posSemidef_iff_eq_conjTranspose_mul_self, sub_eq_iff_eq_add']
    exact ⟨_, hs⟩

end Matrix

[Jireh Loreaux]

Yes, but my preference would be to obliterate docs#Matrix.PosSemidef.sqrt in favor of docs#CFC.sqrt.

[Monica Omar]

should be fairly straightforward to deprecate docs#Matrix.PosSemidef.sqrt and use docs#Matrix.IsHermitian.cfc instead

[Jireh Loreaux]

The only minor trick is that you have to work with the ℝ functional calculus to get the StarOrderedRing instance, instead of the ℝ≥0 one (so you don't have access to docs#CFC.sqrt directly. Nevertheless, you can still use cfc Real.sqrt A and show that it is positive semidefinite.

[Jireh Loreaux]

Aha, that's a footgun, albeit a necessary one. You never want to touch Matrix.IsHermitian.cfc (okay, fine, almost never).

[Monica Omar]

lol okay! I just read the documentation saying to use the generic API instead

[Monica Omar]

so should I work on this? deprecate docs#Matrix.PosSemidef.sqrt and use the continuous functional calculus instead?

[Jireh Loreaux]

Let's make sure there's buy-in from others first, so I'm not just dictating things.

[Monica Omar]

how heavy would you say the imports would be on Matrix/PosDef if we import Matrix/HermitianFunctionalCalculus?

[Monica Omar]

there's only 3 files importing Matrix/PosDef so shouldn't be a big deal actually

[Frédéric Dupuis]

Jireh Loreaux said:

Let's make sure there's buy-in from others first, so I'm not just dictating things.

I agree with the plan, for what it's worth!

[Eric Wieser]

I think we could PR the above without first deciding what to do about sqrt?

[Monica Omar]

I don't see why not, since I don't think the proof would change after changing sqrt (right?)

[Eric Wieser]

Is positive semidefinite-ness really not a partial order on matrices of integers?

[Monica Omar]

should be fine, just coerce it as a matrix of real numbers?

[Monica Omar]

Monica Omar said:

should be fine, just coerce it as a matrix of real numbers?

I don't think it would make sense if you don't coerce it? For example, [[1,1],[1,2]] is positive as a real matrix. But it has non-integer eigenvalues. So I don't know. I think I'm overthinking it now.

[Jireh Loreaux]

Monica, it's not that simple since over ℤ your semi-definiteness condition is weaker than positive semidefiniteness of the corresponding matrix over ℝ (since you only consider vectors with integer entries). What is true is that if the ℝ matrix is positive semidefinite, then so is the ℤ matrix.

I think it would be possible to prove that positive semidefiniteness still induces a partial order, but I think you would need to go through QuadraticMap.polar (or some inline recreation of it) to show that for a positive semidefinite matrix A over ℤ with A ≤ 0, the off-diagonal entries are zero.

In any case, I claim that it's irrelevant because this order on Matrix n n ℤ is not "natural" because it doesn't satisfy StarOrderedRing (because [[1,1],[1,1]] is positive semidefinite but cannot be written as a sum of squares of matrices).

[Jireh Loreaux]

So we should only put the order on RCLike-valued matrices.

[Jireh Loreaux]

(possibly scoped)

[Monica Omar]

I see, I see

[Monica Omar]

should we PR the above star ordered ring instance (scoped to MatrixOrder or something?) for now before deciding whether to move forward with changing the matrix sqrt to use cfc?

[Eric Wieser]

Jireh Loreaux said:

In any case, I claim that it's irrelevant because this order on Matrix n n ℤ is not "natural"

Could it not be said that this means PosSemidef is also not natural?

[Jireh Loreaux]

My point with Matrix.PosSemidef.sqrt is not that we should replace its definition with CFC.sqrt, but that we should remove it entirely. For that reason, I think Eric's comment:
Eric Wieser said:

I think we could PR the above without first deciding what to do about sqrt?

doesn't make any sense because we would be adding more things which depend on something that will be ripped out anyway.

[Jireh Loreaux]

Eric Wieser said:

Could it not be said that this means PosSemidef is also not natural?

I think PosSemidef on ℤ is probably pointless, but :shrug:, but it's not doing so much harm to have a definition lying around unused in certain contexts. It's much worse to have an instance which could be (surprisingly!) available that is useless.

[Monica Omar]

Monica Omar said:

so should I work on this? deprecate docs#Matrix.PosSemidef.sqrt and use the continuous functional calculus instead?

I'll give it a few days to see if there's any objections with removing docs#Matrix.PosSemidef.sqrt, before I work on it then

[Eric Wieser]

To me, defining LE on integer matrices even though they don't form a star-ordered ring doesn't sound any worse than defining subtraction on nat even though it doesn't form an additive group

[Eric Wieser]

I think as long as the LE is defined such that A <= B \iff (A.map Int.cast) <= (B.map Int.cast) then the instance is not surprising

[Monica Omar]

yeah, that's fair

[Monica Omar]

Eric Wieser said:

I think as long as the LE is defined such that A <= B \iff (A.map Int.cast) <= (B.map Int.cast) then the instance is not surprising

so something like this?

variable {n 𝕜} [RCLike 𝕜] [Fintype n] [DecidableEq n]

instance : LE (Matrix n n 𝕜) :=
  ⟨fun x y => (y - x).PosSemidef⟩

def instLE' {R} (f : R → 𝕜) : LE (Matrix n n R) :=
  ⟨fun x y => x.map f ≤ y.map f⟩

instance : LE (Matrix n n ℤ) := instLE' <| @Int.cast ℝ _

example (x y : Matrix n n ℤ) :
    x ≤ y ↔ x.map (@Int.cast ℝ _) ≤ y.map Int.cast :=
  Iff.rfl

[Aaron Liu]

I don't think that instance works since you can't infer 𝕜

[Aaron Liu]

oh what it works that was surprising

[Aaron Liu]

but you get diamonds when it infers the wrong 𝕜

[Monica Omar]

change it to the reals or complex then perhaps?

[Alex Meiburg]

There is probably more code to take from https://github.com/Timeroot/Lean-QuantumInfo/blob/main/QuantumInfo/ForMathlib/Matrix.lean if you would like

[Monica Omar]

I think most of that is already in Mathlib or are PR-ed (like the kronecker product ones)
from skimming the code, I think we only don't have the StarOrderedRing stuff (and related stuff like the partial order etc)

[Alex Meiburg]

I guess I don't know about other PRs, but I think almost none of that is currently in Mathlib
(e.g. no hits on loogle for _ • _, Matrix.PosSemidef, Matrix.trace, Matrix.PosSemidef, Matrix.PosDef, Matrix.submatrix .. there's a lot of lexical gaps!)
but yes I mostly meant the section partialOrder as relevant to this discussion :slight_smile:

[Jireh Loreaux]

Eric Wieser said:

To me, defining LE on integer matrices even though they don't form a star-ordered ring doesn't sound any worse than defining subtraction on nat even though it doesn't form an additive group

I disagree here. There's at least a reason to define - on ℕ (which is that ℕ is ubiquitous and having all the basic operations is important). In contrast, I don't see any real reason to define a partial order on Matrix n n ℤ. Why do we need it?

[Eric Wieser]

The main benefit would be that we can drop PosSemidef entirely and use 0 ≤ M notation

[Jireh Loreaux]

Eric, when do you ever want to talk about positive semidefiniteness of a matrix with integer entries (regardless of which notation you use)?

[Kevin Buzzard]

The Fourier expansion of a degree g Siegel modular form is a sum over positive semidefinite integral binary quadratic forms in g variables, which is almost symmetric positive semidefinite gxg matrices with integer entries, except that half-integers are allowed off the diagonal.

[Jireh Loreaux]

So they are terms of A : Matrix g g ℤ[1/2] (I don't know what the Mathlib notation is for that) satisfying PosSemidef A according to the definition above? Would you like to be able to phrase the latter as 0 ≤ A?

[Kevin Buzzard]

Usually people just sum over all binary quadratic forms and then occasionally remark that the coefficient is zero unless the binary quadratic form is positive semidefinite. This is really a theory about quadratic forms not matrices, it's not Z[1/2] it's Z/2 for example

[Jireh Loreaux]

My question is mainly: despite this being an interesting example of positive semidefinite matrices over things other than ℝ or ℂ occuring in practice, you still never really want to discuss it in the language of positive semidefinite matrices, correct?

[Kevin Buzzard]

I think that formally I would rather discuss it in the language of positive semidefinite binary quadratic forms, because it seems to me to be a bit of a mess talking about positive semidefinite matrices whose diagonal entries are integers and whose off-diagonal entries are half-integers.

[Monica Omar]

Jireh Loreaux said:

My point with Matrix.PosSemidef.sqrt is not that we should replace its definition with CFC.sqrt, but that we should remove it entirely. For that reason, I think Eric's comment:
Eric Wieser said:

I think we could PR the above without first deciding what to do about sqrt?

doesn't make any sense because we would be adding more things which depend on something that will be ripped out anyway.

I'm a bit confused. Are we just literally wanting to remove the definition docs#Matrix.PosSemidef.sqrt or are we wanting to do more? Because, in order to use docs#CFC.sqrt_sq for example, we'd need PartialOrder and StarOrderedRing.
So, I'm not sure how we're supposed to remove results like docs#Matrix.PosSemidef.sqrt_sq if they're used to show PartialOrder and StarOrderedRing on matrices.
All I've done so far is remove the definition PosSemidef.sqrt and use cfc Real.sqrt _ instead in the results.

[Jireh Loreaux]

I can detail this later today or tomorrow.

[Jireh Loreaux]

It can surely be cleaned up, but adding this gives you access to the continuous functional calculus for nonnegative elements (and hence docs#CFC.sqrt).

import Mathlib.Analysis.CStarAlgebra.ContinuousFunctionalCalculus.Instances
import Mathlib.LinearAlgebra.Matrix.HermitianFunctionalCalculus
import Mathlib.LinearAlgebra.Matrix.PosDef

variable {m n 𝕜 : Type*}
variable [Fintype m] [Fintype n]
variable [RCLike 𝕜]

namespace Matrix

open scoped ComplexOrder

instance : PartialOrder (Matrix n n 𝕜) where
  le A B := (B - A).PosSemidef
  le_refl A := (sub_self A) ▸ PosSemidef.zero
  le_trans A B C h₁ h₂ := sub_add_sub_cancel C B A ▸ h₂.add h₁
  le_antisymm A B h₁ h₂ := by
    have foo := neg_sub A B ▸ h₁.trace_nonneg
    have bar := h₂.trace_nonneg
    rw [trace_neg, neg_nonneg] at foo
    have : (A - B).trace = 0 := le_antisymm foo bar
    classical
    rw [← sub_eq_zero, ← h₂.isHermitian.eigenvalues_eq_zero_iff]
    ext i
    rw [h₂.isHermitian.trace_eq_sum_eigenvalues, ← RCLike.ofReal_sum] at this
    norm_cast at this
    rw [← (Finset.univ (α := n)).sum_const_zero, eq_comm,
      Finset.sum_eq_sum_iff_of_le (by simpa using h₂.eigenvalues_nonneg)] at this
    exact this i (by simp) |>.symm

lemma le_iff {A B : Matrix n n 𝕜} : A ≤ B ↔ (B - A).PosSemidef := Iff.rfl

lemma nonneg_iff {A : Matrix n n 𝕜} : 0 ≤ A ↔ A.PosSemidef := by rw [le_iff, sub_zero]

instance : IsOrderedAddMonoid (Matrix n n 𝕜) where
  add_le_add_left _ _ _ _ := by rwa [le_iff, add_sub_add_left_eq_sub]

instance : NonnegSpectrumClass ℝ (Matrix n n 𝕜) where
  quasispectrum_nonneg_of_nonneg A hA := by
    classical
    rw [quasispectrum_eq_spectrum_union_zero ℝ A]
    simp only [Set.union_singleton, Set.mem_insert_iff, forall_eq_or_imp, le_refl, true_and]
    rw [nonneg_iff] at hA
    rw [hA.isHermitian.spectrum_real_eq_range_eigenvalues]
    rintro - ⟨x, rfl⟩
    exact hA.eigenvalues_nonneg x

instance : StarOrderedRing (Matrix n n 𝕜) :=
  .of_nonneg_iff' add_le_add_left fun A ↦ by
    constructor
    · intro hA
      have := QuasispectrumRestricts.nnreal_of_nonneg hA
      rw [nonneg_iff] at hA
      have hA' : IsSelfAdjoint A := hA.isHermitian
      classical
      obtain ⟨X, hX, -, rfl⟩ :=
        CFC.exists_sqrt_of_isSelfAdjoint_of_quasispectrumRestricts hA' this
      exact ⟨X, by rw [hX.star_eq]⟩
    · rintro ⟨A, rfl⟩
      rw [nonneg_iff, star_eq_conjTranspose]
      exact posSemidef_conjTranspose_mul_self A

[Monica Omar]

#29896
I move a lot of results around, so there's a lot of noise