3

I have two functions $f(x,y),g(x,y)$ and a matrix

$M=\begin{pmatrix}\frac{\partial f}{\partial x} & \frac{\partial f}{\partial y}\\\frac{\partial g}{\partial x} & \frac{\partial g}{\partial y}\end{pmatrix}$.

Because the matrix elements are functions of $x$ and $y$, I expect the eigenvalues of $M$ and the components of the eigenvectors of $M$ to also be functions of $x$ and $y$.

Here's a little bit of code:

f=...;
g=...;
M={{D[f,x],D[f,y]},{D[g,x],D[g,y]}};

(note that I haven't given the explicit expressions for $f$ and $g$ because they're unimportant for the purposes of this question; suffice it to say that they're complicated enough for their partial derivatives to be horrendous).

I can ask Mathematica for the eigenvalues of $M$, and I get a pair of pretty decent (if horrible-looking) answers which are indeed functions of $x$ and $y$. If I attempt to extract the eigenvectors of $M$, however, Mathematica runs for a while and ultimately gives up:

Eigenvectors[M]

Eigenvectors::eivec0: Unable to find all eigenvectors. >>

{{0,0},{0,0}}

I've therefore resorted to graphing the (absolute values of the) components of the eigenvectors for different values of $x$ and $y$, which takes forever but works. Here's my code:

Quiet[
 DensityPlot[
  Abs[Eigenvectors[M][[1]][[1]]],
  {x,x1,x2},{y,y1,y2},
  Exclusions->None,
  ColorFunction->"Rainbow",
  PlotRangePadding->0,
  PlotLabel->Row[{Subscript[v,Subscript[1,1]]}],
  PlotLegends->Automatic
 ]
]

and similarly for the other three components. Here's an example of one of the components:

enter image description here

Now, the Mathematica documentation for Eigenvalues says

If they are numeric, eigenvalues are sorted in order of decreasing absolute value.

and the documentation for Eigenvectors says

Eigenvectors with numeric eigenvalues are sorted in order of decreasing absolute value of their eigenvalues.

and

For exact or symbolic matrices $m$, the eigenvectors are not normalized.

I therefore have four questions:

  1. If one eigenvalue has a higher absolute value than the other for certain values of $x$ and $y$ but a lower absolute value than the other for other values of $x$ and $y$, does the order change? In other words, if I ask Mathematica to graph the first eigenvector (i.e. that with the highest absolute value, as per the documentation), do I get one eigenvector for certain values of $x$ and $y$ and the other one for the other values of $x$ and $y$?

  2. If the eigenvalues are not normalised, which criterion does Mathematica use to determine which multiples of the normalised eigenvectors to display? (Any linear combination of eigenvectors of a matrix is also an eigenvector of that matrix, so whatever Mathematica returns is $k_1$ times the normalised first eigenvector and $k_2$ times the normalised second eigenvector; how does Mathematica choose the $k_j$? I seem to recall always getting eigenvectors where one of the components equals 1; is this the case? If so, which component is equal to 1?)

  3. Here I will use the word "normalise" in a general sense; I know that the eigenvectors are not normalised to 1 (as per the documentation), but they might be normalised to something else (which may or may not be different for each eigenvector; that is relevant for question 2 but not for this question). If the eigenvectors are indeed normalised, does Mathematica renormalise them for each combination of values of $x$ and $y$? (If so, my graphs of the components of the eigenvectors are meaningless except as a gauge of how large the quotient between the two components is.)

  4. Regardless of whether the eigenvalues of $M$ are functions of $x$ and $y$ or not, if I call them $e_1$ and $e_2$ (in the order in which Mathematica gives them), is there any way of knowing whether the diagonalised expression for $M$ is $\begin{pmatrix}e_1 & 0\\0 & e_2\end{pmatrix}$ or $\begin{pmatrix}e_2 & 0\\0 & e_1\end{pmatrix}$?

I am aware of this question, but that concerns only the relative ordering of eigenvectors to the ordering of eigenvalues, which is already answered in the documentation (and included in the quotes above). I am also aware of this other question, but that appears to be about intentionally reordering the eigenvalues rather than finding out whether Mathematica does it by default or not.

Thanks for any help.

J. M.'s missing motivation
  • 124,525
  • 11
  • 401
  • 574
Rain
  • 636
  • 3
  • 12
  • 5
    BTW: D[{f, g}, {{x, y}}] is a simpler way to compute your Jacobian. – J. M.'s missing motivation Sep 06 '17 at 21:45
  • Thanks; I didn't know that! – Rain Sep 06 '17 at 21:47
  • 1
    Closely related: How get eigenvectors without phase jump? and Should eigenvalues be ordered?. As for your fourth question: it's possible to do what you're asking if the matrix M has symmetries according to which you can label the eigenvectors. Then the similarity transformation diagonalizing M would be a symmetry operation parametrized by x and y and could be written down independently. It's impossible to say anything else in the absence of more information about M. – Jens Sep 06 '17 at 22:35
  • @Jens Thanks. As far as I was able to understand (some of the code is arcane to me), those reordering methods seem to require that Mathematica be able to obtain explicit values or expressions for the eigenvalues/eigenvectors, which is not true in this case. – Rain Sep 06 '17 at 23:08
  • @Jens Concerning the symmetries, I'm not sure, but I wouldn't bet on my matrix having tge necessary symmetries. It's horrible enough for me not to really want to wade through all that algebra to find potential symmetries. I understand you can't really much more without any informatiok about the matrix, and I apologise for keeping you in the dark, but I think I can use physical reasoning to determine which is the correct order of the eigenvalues and eigenvectors in this particular case. – Rain Sep 06 '17 at 23:11
  • When you say "physical reasoning," that to me implies that there is some other label you're attaching to the eigenvalues. Where does that come from if not from symmetries? – Jens Sep 06 '17 at 23:13
  • @Jens I can graph $f$ and $g$ as functions of $x$ and $y$ and see in which direction(s) they change most. The eigenvalues and eigenvectors are related to the partial derivatives of $f$ and $g$ with respect to $x$ and $y$. If, for example, $f$ changes equally strongly in the $x$ and $y$ directions but $g$ changes appreciably only in the $x$ direction and the eigenvectors are $(0.995,0.100)$ and $(0.707,0.707)$, then I know that the second eigenvector corresponds to the change of $f$ and the first eigenvector corresponds to the change of $g$. Hope that made sense. – Rain Sep 06 '17 at 23:26
  • 1
    It doesn't look like a foolproof prescription, especially near degeneracies. So I'm afraid there isn't a formalizable Mathematica answer I can think of. – Jens Sep 06 '17 at 23:47

1 Answers1

1

I've run some more tests and have found the answer to questions 2 and 3. Mathematica does seem to be renormalising with each value of $x$ and $y$. Here's a graph of the norm of one of the eigenvectors:

Quiet[
 Plot3D[
  Sqrt[Abs[Eigenvectors[m][[2]][[1]]]^2+Abs[Eigenvectors[m][[2]][[2]]]^2]-1,
  {x,x1,x2},{y,y1,y2},
  Exclusions->None,
  PlotLabel->Row[{"|",Subscript[v,2],"|"}],
  PlotRangePadding->0
 ]
]

enter image description here

I'll ascribe the variation in the norm (of the order of $10^{-16}$) to machine precision (i.e. the norm of this eigenvector is $1$ throughout the domain). This seems to imply that Mathematica is not only renormalising the eigenvectors at each point, but it is actually giving me the eigenvectors normalised to 1 (contrary to what the documentation says).

Of course, this assumes that the order of the eigenvalues and eigenvectors is preserved regardless of which eigenvalue has a larger absolute value at each point. If it isn't, my whole interpretation of the above graph automatically goes out the window.

Edit:

I've found the answer to question 1 by testing a little further. I built a matrix with eigenvalues whose ordering should change in a very obvious way:

M=DiagonalMatrix[{x-2,2-x}];

Graphing the eigenvalues explicitly (i.e. $x-2$ and $2-x$) yields the following:

Plot[{2-x,x-2},{x,0,4}]

enter image description here

However, graphing the eigenvalues as eigenvalues has a different result:

Plot[{Eigenvalues[M][[1]],Eigenvalues[M][[2]]},{x,0,4}]

enter image description here

It is as I feared: Eigenvalues are reordered at each point. Clearly, then, eigenvectors are also reordered at each point, since the eigenvalue and eigenvector ordering must be consistent at all times (at least according to Mathematica). (This also renders my analysis on the answers to questions 2 and 3 wrong, although the ultimate answers (i.e. Mathematica renormalises to $1$ at every point) is still right because of what J.M. said in the comments.)

Other than having explicit expressions for the eigenvector components (which I do not have, as explained at the beginning of my original post), is there a simple way around this reordering?

Another edit:

Finally, as mentioned in the comment section, I've found a way around the problem posed in question 4, but I imagine it's highly dependent on the type of problem one is trying to solve.

In my case, I can graph $f(x,y)$ and $g(x,y)$ and see in which direction(s) they change most strongly. The eigenvalues and eigenvectors of $M$ are related to the partial derivatives of $f$ and $g$ with respect to $x$ and $y$, so I can use the (qualitative) visual information obtained from the graphs of $f$ and $g$ to "guess" which order the eigenvectors and the corresponding eigenvalues should be in (at least in regions where there is a clear difference between the directions of maximum variation of $f$ and $g$).

This is not a very satisfactory answer to question 4, but it is the only one I have so far.

Rain
  • 636
  • 3
  • 12
  • 1
    "t is actually giving me the eigenvectors normalised to 1" - since plotting functions are feeding inexact numbers to your function, inexact matrices are being generated, and indeed the internal routines return normalized eigenvectors for symmetric inexact matrices. – J. M.'s missing motivation Sep 06 '17 at 21:47
  • @J.M. So when the documentation says "For approximate numerical matrices m, the eigenvectors are normalized." it refers to what you said? – Rain Sep 06 '17 at 21:48
  • Yep. Regarding the order of eigenvalues: the eigensystem is sorted by absolute value by default, as demonstrated by evaluating Eigenvalues[DiagonalMatrix[{-1, 5}]] and Eigenvalues[DiagonalMatrix[{1, -5}]]. – J. M.'s missing motivation Sep 06 '17 at 21:57
  • @J.M. I understand that they're ordered by absolute value by default. Does that mean that they are reordered as $x$ and $y$ change if instead of being $-5$ and $1$ they are functions of $x$ and $y$ (question 1)? And how can I know whether my matrix is DiagonalMatrix[{1,-5}] or DiagonalMatrix[{-5,1}] if I look at the eigenvalues ($-5$ and $1$ in this case) Mathematica has given me (question 4)? – Rain Sep 06 '17 at 22:03
  • 1
    Looks like you answered your question yourself. :) – J. M.'s missing motivation Sep 07 '17 at 02:24