Matrix square root
I need a function equivalent to Matlab's sqrtm, i.e., a square root for matrices. I constructed the following definition: def sqrtm(M): (U,S,VT) = LinearAlgebra.singular_value_decomposition(M) D = MLab.diag(sqrt(S)) return matrixmultiply(matrixmultiply(U,D),VT) but this technique only works for Hermitian, positive definite matrices. The matrices I operate on don't necessarily satisfy that condition. Does anybody know of a more general technique, preferably with an accompanying Python implementation? Thanks. Andrew.
On Thursday 04 September 2003 16:42, Andrew Nesbit wrote:
I need a function equivalent to Matlab's sqrtm, i.e., a square root for matrices.
I constructed the following definition:
def sqrtm(M): (U,S,VT) = LinearAlgebra.singular_value_decomposition(M) D = MLab.diag(sqrt(S)) return matrixmultiply(matrixmultiply(U,D),VT)
but this technique only works for Hermitian, positive definite matrices. The matrices I operate on don't necessarily satisfy that condition.
I'd use an eigenvalue decomposition, then take the square root of the eigenvalues, and then apply the diagonlization matrix in reverse. If you convert to eigenvalues to complex before taking the square root, this will work for nonpositivedefinite matrices, yielding a complex result. Konrad.   Konrad Hinsen  EMail: hinsen@cnrsorleans.fr Centre de Biophysique Moleculaire (CNRS)  Tel.: +332.38.25.56.24 Rue Charles Sadron  Fax: +332.38.63.15.17 45071 Orleans Cedex 2  Deutsch/Esperanto/English/ France  Nederlands/Francais 
Konrad Hinsen
On Thursday 04 September 2003 16:42, Andrew Nesbit wrote:
I need a function equivalent to Matlab's sqrtm, i.e., a square root for matrices.
[snip]
I'd use an eigenvalue decomposition, then take the square root of the eigenvalues, and then apply the diagonlization matrix in reverse. If you convert to eigenvalues to complex before taking the square root, this will work for nonpositivedefinite matrices, yielding a complex result.
Thankyou for the advice. This gives me a good starting point. Andrew.
Andrew Nesbit wrote:
Konrad Hinsen
writes: On Thursday 04 September 2003 16:42, Andrew Nesbit wrote:
I need a function equivalent to Matlab's sqrtm, i.e., a square root for matrices
[snip]
I'd use an eigenvalue decomposition, then take the square root of the eigenvalues, and then apply the diagonlization matrix in reverse. If you convert to eigenvalues to complex before taking the square root, this will work for nonpositivedefinite matrices, yielding a complex result
Thankyou for the advice. This gives me a good starting point.
Andrew.
The only problem with this approach is that you have to find all the eigenvalues and eigenvectors to get the spectral decomposition. For large, poorly conditioned matrices this can be a source of significant errors, and a lot of work. An alternative, and often faster and more accurate approach is to use a series expansion for the square root. Evaluating the series expansion only requires multiplication. Bernie
Bernard Frankpitt
The only problem with this approach is that you have to find all the eigenvalues and eigenvectors to get the spectral decomposition. For large, poorly conditioned matrices this can be a source of significant errors, and a lot of work. An alternative, and often faster and more accurate approach is to use a series expansion for the square root. Evaluating the series expansion only requires multiplication.
This technique sounds like it's definitely worth checking out. Thankyou for the tip. Andrew.
participants (4)

Andrew Nesbit

Andrew Nesbit

Bernard Frankpitt

Konrad Hinsen