# Eigenvectors are divided by zero depending on evaluation, 4x4 matrix

I calculated the Eigenvectors of the $$4\times 4$$ matrix $$m$$ with parameters $$(a,b,c,d)\in \mathbb{R_{\ge 0}}$$. If I set e.g. $$b\rightarrow 0$$ after calculation, then some of the Eigenvectors are affected by division by zero. If I set $$b\rightarrow 0$$ before I calculate the Eigenvectors I get a solution without division by zero. The situation $$b\rightarrow 0$$ is only one among of other examples when the denominator is zero. I am not looking for a solution that works only for this special case. How can I avoid division by zero if I want to replace the values $$a,b,c,d$$ of the Eigenvectors?

A more complicated case is discussed in this post.

$Assumptions = a>=0 && b>=0 && c>=0 && d>=0; m = {{0,0,a c,-a b}, {0,0,0,b d}, {a c,0,0,0}, {-a b,b d,0,0}}; (* b->0 after calculation of Eigensystem *) {l,v} = Limit[Eigensystem[m],b->0,Direction->"FromAbove"]; l (* {0,0,-a c,a c} *) v (* {{0,-1,0,1},{0,1,0,1},{Infinity,0,-Infinity,1},{-Infinity,0,-Infinity,1}} *) (* b->0 before calculation of Eigensystem *) {l,v} = Eigensystem[m /. b->0]; l (* {-a c,a c,0,0} *) v (* {{-1,0,1,0},{1,0,1,0},{0,0,0,1},{0,1,0,0}} *)  $Version
(* 13.3.1 for Microsoft Windows (64-bit) (July 24, 2023) *)


Comment to the answers, Nov. 13th, 2023

The answers use tricks to avoid division by zero. But we learn nothing from this and don't see which steps take place. Actually I thought the expression for the Eigenvectors can be reformulated so that any division by zero is avoided at all. Considering the answers it would be easier to set the parameters to zero before we calculate the Eigenvectors. Then we also learn nothing but it is faster and without tricks.

You have to use FullSimplify:

Clear["Global*"];

m = {{0, 0, a c, -a b}, {0, 0, 0, b d}, {a c, 0, 0, 0}, {-a b, b d, 0, 0}};

v = Eigenvectors[m] // FullSimplify;

v /. b -> 0


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

Limit[v, b -> 0]


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

Eigenvectors[m /. b -> 0]


{{-1, 0, 1, 0}, {1, 0, 1, 0}, {0, 0, 0, 1}, {0, 1, 0, 0}}

$Version  13.3.0 for Mac OS X ARM (64-bit) (June 3, 2023) • This cannot be correct. If we set$b\rightarrow0$before or after calculation of Eigenvalues we get always the 4 Eigenvalues$(−ac,ac,0,0)$. This means that at least two linear independent Eigenvectors must exist. However in this answer all Eigenvectors are the same. They are only linear independent if we set$b\rightarrow0$before we calculated the Eigenvectors. However the task is to set$b\rightarrow0$after calculation of Eigenvectors. Commented Nov 12, 2023 at 16:21 As you set the limit of b->0, some components will approach infinity, causing inconvenience. Fortunately, as we are dealing with eigenvectors, we can scale them arbitrarily. So why not divide each eigenvector by the element with largest absolute value in them, and then all elements of the vector will always be finite (and <= 1)! So for each eigenvector, we first use Fold to find out the maximum component as b->0 and then divide the whole eigenvector by that component to get an eigenvector with elements all less or equal to 1 as b->0. The code would be: m = {{0, 0, a c, -a b}, {0, 0, 0, b d}, {a c, 0, 0, 0}, {-a b, b d, 0, 0}};$Assumptions = a >= 0 && b >= 0 && c >= 0 && d >= 0;
(*b->0 after calculation of Eigensystem*)
{l, v} = Eigensystem[m];
Limit[l, b -> 0, Direction -> "FromAbove"]
(*{0,0,-a c,a c}*)
Limit[(#/
Fold[If[Abs[Limit[#1/#2, b -> 0, Direction -> "FromAbove"]] >=
1, #1, #2] &, #]), b -> 0, Direction -> "FromAbove"] & /@ v
(*{{0,1,0,-1},{0,1,0,1},{1,0,-1,0},{1,0,1,0}}*)


This code is undoubtfully slow as it's computing symbolic limits again and again, and it could fail if the eigenvectors have symbolic components in them. To deal with this issue, a hacky approach could be to change the judging part in the previous code to:

If[Abs[Limit[#2/#1, b -> 0, Direction -> "FromAbove"]/.Thread[{a,b,c,d}->RandomReal[{1,2},4]]] < Infinity, #1, #2]&


It might cause problems in edge cases, but should work for most applications.

• Is this a clean mathematical approach? If we get $(0,1,-\infty,\infty)$ then you suggest to divide by $\infty$ to get $(0,0,-1,1)$. Commented Nov 12, 2023 at 23:26
• −∞/∞ is not always -1, it really depends on the expresion. i.e. limit of (1/b)/(1/b^2) as b->0+ is 0 while both diverges to infinity. I think the mathematical cleaniness depends on how you define the 'limit of eigenvector'.
– Wjx
Commented Nov 13, 2023 at 0:27
• It is not a question of definition of limit as the eigenvectors exist and are unique (except of linear factor), Commented Nov 13, 2023 at 0:34
• Well my code is just dividing everything with a 'linear factor', and the linear factor is the fastest diverging element. So it should be mathematically sound except for some marginal cases.
– Wjx
Commented Nov 13, 2023 at 8:15

Do

\$Assumptions = a >= 0 && b >= 0 && c >= 0 && d >= 0;
m = {{0, 0, a c, -a b}, {0, 0, 0, b d}, {a c, 0, 0, 0}, {-a b, b d, 0, 0}};
es = Eigensystem[m];


and inspect the eigenvectors es[[2]]. You will notice that the fourth component is always 1. That is because MMA applies a Gauss algorithm to m-eigenvalue IdentityMatrix[4], then fixes the fourth eigenvector component to 1 and calculates the other components w.r.t. the fourth component. So what happens when, for some parameter combination, this component vanishes? The other components blow up! How to fix this? Just expand the eigenvectors by the product of the denominators of all components like so:

 eigenvectors = FullSimplify[Table[es[[2, i]] Times @@ Map[Denominator, es[[2, i, ;;]]], {i, 1, 4}]]


You can check with

Table[(m - es[[1, i]] IdentityMatrix[4]) . eigenvectors[[i]] //FullSimplify, {i, 1, 4}]


that eigenvectors indeed contain the eigenvectors. Also, the eigenvectors behave well when inserting a critical parameter,

eigenvectors /. b -> 0


gives

{{0, -Sqrt[a^2 c^2 - Sqrt[a^4 c^4]] (a^2 c^2 + Sqrt[a^4 c^4])^3, 0,   0},
{0, Sqrt[a^2 c^2 - Sqrt[a^4 c^4]] (a^2 c^2 + Sqrt[a^4 c^4])^3, 0, 0},
{0, -(-a^2 c^2 + Sqrt[a^4 c^4])^3 Sqrt[a^2 c^2 + Sqrt[a^4 c^4]], 0,  0},
{0, (-a^2 c^2 + Sqrt[a^4 c^4])^3 Sqrt[a^2 c^2 + Sqrt[a^4 c^4]],   0, 0}}
`
• After solving the roots all Eigenvectors become zero vectors. Commented Nov 24, 2023 at 9:03