1

Suppose we have a Hermitian matrix A that is known to have an inverse.

I know that ZGETRF and ZGETRI subroutines in LAPACK library can compute the inverse matrix.

Is there any subroutine in LAPACK or BLAS library can calculate A^{-1/2} directly or any other way to compute A^{-1/2}?

Ian Bush
  • 6,996
  • 1
  • 21
  • 27
Kieran
  • 123
  • 1
  • 10
  • No, there are no LAPACK routines to do this directly, though you can probably use a few routines together to get what you want, assuming it exists. Does A have any other nice properties, e.g. Hermitian? And *why* do you want this? – Ian Bush Apr 13 '22 at 17:01
  • Does this answer your question? [exponential of a matrix in Fortran using Taylor Series expansion](https://stackoverflow.com/questions/70728907/exponential-of-a-matrix-in-fortran-using-taylor-series-expansion). Oh, no sorry, thats `e^matrix`, not `matrix^n`. I guess do the same thing, but with `(e_i)^{-1/}` rather than `e^{e_i}`. – veryreverie Apr 13 '22 at 17:09
  • @Ian Thank you for the reply. The complex matrix A is Hermitian. I need to compute A^{-1/2} as the intermediate data for the following calculation. – Kieran Apr 13 '22 at 17:10
  • @veryreverie Thank you for the reply. I only want compute the exponential value of the complex matrix and the exponential is -1/2. That post explains how to compute the Taylor expansion with a matrix as the exponential term. It is different from what I want to calculate. – Kieran Apr 13 '22 at 17:14
  • @Kieran I've taken the liberty of changing your title, as I believe the "exponential of a matrix" refers to `e^A`, not `A^n`. See e.g. [here](https://en.wikipedia.org/wiki/Matrix_exponential). Feel free to rollback if you disagree. – veryreverie Apr 13 '22 at 17:16
  • 4
    OK, diagonalising it with ZHEEV and find it that way is as easy as any - I don't claim it is optimal but it will get the job done. But *what* is the next step? If e.g. this is part of solving a generalised eigenvalue problem and A is also positive definite there are better ways to solve it than this. – Ian Bush Apr 13 '22 at 17:19

1 Answers1

1

You can raise a matrix to a power following a similar procedure to taking the exponential of a matrix:

  1. Diagonalise the matrix, to give the eigenvectors v_i and corresponding eigenvalues e_i.
  2. Raise the eigenvalues to the power, {e_i}^{-1/2}.
  3. Construct the matrix whose eigenalues are {e_i}^{-1/2} and whose eigenvectors are v_i.

It's worth noting that, as described here, this problem does not have a unique solution. In step 2 above, both {e_i}^{-1/2} and -{e_i}^{-1/2} will lead to valid solutions, so an N*N matrix A will have at least 2^N matrices B such that B^{-2}=A. If any of the eigenvalues are degenerate then there will be a continuous space of valid solutions.

veryreverie
  • 2,871
  • 2
  • 13
  • 26
  • Thank you for the suggested solution. Can I ask one more question? If A is a 2000 by 2000 complex matrix. I think that diagonalising it may encounters the 'out-of-RAM' issue. Do you have any other idea to compute A^{-1/2}? Thank you again. – Kieran Apr 14 '22 at 16:41
  • 1
    https://en.wikipedia.org/wiki/Square_root_of_a_matrix has some ideas. I can't give any more concrete suggestions without knowing more about your problem. I suggest asking a new Stack Overflow question with more details. – veryreverie Apr 14 '22 at 17:19