13

Well, probably a hard question, but I think it's better to cry out loud :).

I noticed that MatrixExp has a Method option when writing this answer, which is undocumented. Since this tutorial claims that MatrixExp uses variable-order Padé approximation, evaluating rational matrix functions using Paterson–Stockmeyer methods or Krylov subspace approximations, I guess there might exist non-Automatic Method(s), but I can't find it out. Sadly, techniques in this post seem not to help. Does anyone know? Or is Automatic really the only available option?

Community
  • 1
xzczd
  • 65,995
  • 9
  • 163
  • 468

1 Answers1

11

Experimentation based on the documentation you quoted led to two valid Method options:

MatrixExp[{{1.2, 5.6}, {3, 4}}, Method -> "Pade"]

{{346.557, 661.735}, {354.501, 677.425}}

MatrixExp[{{1.2, 5.6}, {3, 4}}, {1, 2}, Method -> "Krylov"]

{1670.03, 1709.35}

If "Krylov" is used for the single parameter syntax it complains:

MatrixExp[{{1.2, 5.6}, {3, 4}}, Method -> "Krylov"]

MatrixExp::novec: The method Krylov requires the specification of a vector. >>

"Pade" fails on exact or symbolic input:

MatrixExp[{{1, 2}, {3, 4}}, Method -> "Pade"]

MatrixExp::invmtd: Invalid method Pade for input with precision ∞. >>

Thanks to ilian we now have the last(?) piece of the puzzle:

MatrixExp[{{2, 0, 0}, {0, 1, -1}, {0, 1, 1}}, Method -> "BlockDecomposition"]

{{E^2, 0, 0}, {0, E Cos[1], -E Sin[1]}, {0, E Sin[1], E Cos[1]}}
Community
  • 1
Mr.Wizard
  • 271,378
  • 34
  • 587
  • 1,371
  • 3
    From the OP's quote, that would seem to be the only two settings. In particular, the "Krylov" setting is valuable for exponentiating sparse matrices, and is quite useful for the action of the exponential of a sparse matrix on a vector (that is, MatrixExp[A, v, Method -> "Krylov"] for computing $\exp(\mathbf A)\mathbf v$). – J. M.'s missing motivation May 18 '15 at 09:34
  • @J. M. What method do you suppose is used for MatrixExp[{{1, 2}, {3, 4}}] as neither "Pade" nor "Krylov" is accepted? (I mean is there something that correlates to Automatic in this case.) – Mr.Wizard May 18 '15 at 09:36
  • Certainly, one does not use a Padé approximant for computing exact results. A look at the docs tells me that either Jordan decomposition or Putzer's algorithm is used. – J. M.'s missing motivation May 18 '15 at 09:58
  • @J. M. Neither "Jordan" nor "Putzer" work either. I am just curious if there is a third valid setting or if whatever is used is only accessible via Automatic. (Oh, and thank you.) – Mr.Wizard May 18 '15 at 10:00
  • 1
    I'm guessing MatrixExp[] is smart enough not to try to compute the eigenvalues of a symbolic matrix (imagine how many Root[] objects will that yield!), and will just use Putzer exclusively in that case. For matrices with exact entries, prolly a toss-up, and the precise details are known only to WRI devs. – J. M.'s missing motivation May 18 '15 at 10:06
  • 1
    On that note, I'm surprised that there is no option to use the Schur decomposition for computing MatrixExp[] for inexact matrices, but a little bird told me MatrixFunction[] does use Schur. I guess one can use that if one suspects that Padé is returning something weird. – J. M.'s missing motivation May 18 '15 at 10:28
  • 4
    The third valid setting can be accessed with Method -> "BlockDecomposition". – ilian May 18 '15 at 12:52
  • @ilian, interesting. That's the Schur-based method, I take it? – J. M.'s missing motivation May 18 '15 at 13:15
  • 1
    More along the lines of splitting the matrix into blocks and applying the Putzer or Jordan decomposition to each. – ilian May 18 '15 at 13:39
  • @ilian Thanks! Added to my answer, with credit of course. – Mr.Wizard May 18 '15 at 13:41
  • @ilian, thanks for the detail! If it isn't too revealing, can you tell how MatrixExp[] picks between Jordan and Putzer when given a matrix with exact entries? – J. M.'s missing motivation May 18 '15 at 14:13
  • 1
    Can't really go into more detail on the heuristic, but your guess above is pretty accurate in that all 'non-simple' cases tend to go to Putzer. – ilian May 18 '15 at 14:49
  • 1
    Useful discussions and reponses. I have an inquiry. I am dealing with a very large sparse matrix and its matrix exponential on a vector. I am computing the solution using MatrixExp[mat, vector, Method -> "Krylov"]; // AbsoluteTiming, but it is still a bit time-consuming due to a very large size of my matrix. So, I want to know that is there a way to speed up this computation by imposing a tolerance or compiling so as to gain some speed? A result of lower accuracy is OK for me – Faz Dec 29 '17 at 18:57
  • 1
    @Fazlollah that is probably worth posting as a separate Question if you have not done so already. I at least don't have an answer for you. – Mr.Wizard Dec 30 '17 at 04:57