Opened 4 years ago
Closed 3 years ago
#23424 closed enhancement (duplicate)
Compute square root of positive semidefinite matrix
Reported by: | stumpc5 | Owned by: | |
---|---|---|---|
Priority: | major | Milestone: | sage-8.1 |
Component: | linear algebra | Keywords: | matrix root, days93.51 |
Cc: | dimpase, nthiery, jmichel, tscrim | Merged in: | |
Authors: | Reviewers: | ||
Report Upstream: | N/A | Work issues: | |
Branch: | Commit: | ||
Dependencies: | Stopgaps: |
Description (last modified by )
This is more a question, I do not have checked how to actually do it.
An (Hermitian) positive semidefinite matrix M has a unique positive semidefinite square root T, this is, M = T^{2. }
I would like to have this square root implemented as I need it for some further computations. This is already possible numerically as in
sage: W = ReflectionGroup(['B',5]) sage: M = W.invariant_form() sage: M [ 1 -1 0 0 0] [-1 2 -1 0 0] [ 0 -1 2 -1 0] [ 0 0 -1 2 -1] [ 0 0 0 -1 2] sage: D, S = M.eigenmatrix_right() sage: S [ 1 1 1 1 1] [ 0.9189859472289948? 0.3097214678905701? -0.7153703234534297? -1.830830026003773? -2.682507065662362?] [ 0.7635211184333675? -0.5943511444371404? -1.203615623775565? 0.5211085581132027? 3.513337091666136?] [ 0.5462003494572026? -1.088155921225222? 0.3727855977717922? 1.397877389115792? -3.228707415119565?] [ 0.2846296765465703? -0.8308300260037728? 1.309721467890570? -1.682507065662363? 1.918985947228995?] sage: T = S*diagonal_matrix([sqrt(x) for x in D.diagonal()])*S^-1 sage: T^2 == M True sage: T.is_positive_definite() True sage: T.is_hermitian() True
I wonder whether this guarantees the matrix square root to be positive definite whenever one starts with a positive definite matrix.
In this and a few other examples it actually is, but the documentation does not say anything about it.
Change History (5)
comment:1 Changed 4 years ago by
- Cc jmichel added
comment:2 Changed 4 years ago by
- Status changed from new to needs_info
comment:3 Changed 4 years ago by
- Description modified (diff)
Oh, I indeed missed this. Updated description.
comment:4 Changed 4 years ago by
IMHO it is more of a mathematical question whether you get a p.s.d. T.
In the basis of the columns of S, one has the quadratic form, given by M in the standard basis, diagonal, given by D, with nonnegative entries (i.e. M=SDS^{-1}). The square root T is a power series in M, and thus given by T=SD^{1/2}S^{-1}, a p.s.d. matrix with eigenvalues given by D^{1/2}, thus all positive whenever D>0.
What can still go wrong is that T is not equal to T^{*}.
However we know that D=D^{*} and so M^{*}=M=SDS^{-1}=S^{*-1}DS^{*}, implying S^{*}SD=DS^{*}S. Thus any power of D commutes with S^{*}S, thus S^{*}SD^{1/2}=D^{1/2}S^{*}S implying T=T^{*}.
Anything I still miss?
comment:5 Changed 3 years ago by
- Cc tscrim added
- Keywords days93.51 added
- Resolution set to duplicate
- Status changed from needs_info to closed
duplicate of #25464
The example in the ticket description is actually an exact computation, just displayed in a slightly misleading way. Indeed, one computes in
QQbar
here:Thus, what is the question? This way it seems like an implementation is already there.