Diagonalising a Skew-Symmetric Matrix
16 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
Jonah Berean-Dutcher
am 22 Sep. 2020
Kommentiert: Christine Tobler
am 15 Apr. 2021
I have a skew-symmetric matrix B, and when I run:
[V,D] = eig(B)
the V matrix returned is not unitary, as I desire it to be. I think this is becuase I should not be using 'eig' for this purpose. What is the correct function, or algorithm to be using for diagonalizing a skew-symmetric matrix, in a complex vector space?
0 Kommentare
Akzeptierte Antwort
Christine Tobler
am 28 Sep. 2020
When EIG is called with an exactly symmetric/hermitian matrix, MATLAB falls back to a specialized algorithm that guarantees that U is orthogonal/unitary and that D is real. There is no such special algorithm choice for skew-symmetric matrices, so there is no guarantee here, even though if the problem is nicely conditioned, the result will be close to that:
>> rng default; A = randn(10); A = A - A';
>> [U, D] = eig(A);
>> max(abs(real(diag(D))))
ans =
2.1034e-16
>> norm(U'*U - eye(10))
ans =
4.7239e-15
However, if matrix B is (exactly) skew-symmetric, it implies that matrix A = 1i*B is hermitian, and passing this matrix to EIG will result in unitary eigenvectors and all-real eigenvalues, which you can then transform back:
[U, D] = eig(1i*A);
D = -1i*D;
>> max(abs(real(diag(D))))
ans =
0
>> norm(U'*U - eye(10))
ans =
1.5001e-15
1 Kommentar
Christine Tobler
am 15 Apr. 2021
As of R2021a, eig will now detect if a matrix is skew-symmetric or skew-hermitian, and will return purely imaginary eigenvalues and orthogonal eigenvectors in that case. See the Release Notes. This does require the input to be exactly skew-hermitian, not just up to round-off (you can use ishermitian(A, 'skew') to check this).
Thank you for pointing out this case as a useful special case treatment to add to EIG, it's helped us improve the functionality.
Weitere Antworten (3)
sushanth govinahallisathyanarayana
am 22 Sep. 2020
I tested this out in Matlab R2018 A
m=[0,-1;1,0];
[V,D]=eig(m)
V'*V seems to be identity for me, m is skew symmetric here. I may be missing something about your original matrix, or you could check to see if it really is skew symmetric.
Paul
am 22 Sep. 2020
Bearbeitet: Paul
am 23 Sep. 2020
For your matrix B, you can diagonalize it and get the associated trasnsformation matrix as follows:
[T,J]=jordan(B);
any(any(J-diag(diag(J)))) % prove J is diagonal
ans =
logical
0
You can inspect J and see that the diagonal elements are equal to eig(B).
max(abs([cplxpair(diag(J))-cplxpair(eig(B))]))
ans =
4.4409e-16
2 Kommentare
Paul
am 23 Sep. 2020
Bearbeitet: Paul
am 23 Sep. 2020
It certainly appears that a) the columns of T are linearly independent and b) T contains eigenvectors of B that, as we've seen, diagonalizes B:
>> rank(T)
ans =
12
>> max(abs(B*T-T*J),[],'all')
ans =
0
Furthermore, now that I think about it some more, if B is defective then it should have at least one entry of 1 on the superdiagonal of J, But it doesn't. Which again indicates that B is not defective. Let's try the symbolic approach:
[Vs,Ds]=eig(sym(B));
Vs\B*Vs
ans =
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, -2i, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 2i, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, -1i, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, -1i, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, -1i, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 1i, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1i, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1i]
>> [Vs-T]
ans =
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
So it certainly seems like B should be diagonalizable, and it is either via the transformation returned from Jordan, or via the eigenvectors returned from the Symbolic toolbox.
Try a different option for eig, which seems to work better for your case.
[V1,D1]=eig(B,'nobalance');
rank(V1)
ans =
12
max(max(abs(V1\B*V1-D1)))
ans =
1.7114e-08
John D'Errico
am 23 Sep. 2020
Bearbeitet: John D'Errico
am 23 Sep. 2020
Matrices that are defective will not have a complete set of eigenvectors.
" In particular, an n × n matrix is defective if and only if it does not have n linearly independent eigenvectors."
The classic example that I know of is
>> A = triu(ones(2));
>> [V,D] = eig(A);
V =
1 -1
0 2.22044604925031e-16
D =
1 0
0 1
As you can see, there is a duplicate eigenvalue. But you can also see the two eigenvectors (columns of V) are not orthogonal.
In this case, the output from eig still satisfies the relation that A*V == V*D, at least within floating point trash.
>> norm(A*V - V*D)
ans =
2.22044604925031e-16
How about the case of your matrix B?
>> [V,D] = eig(B);
>> norm(B*V - V*D)
ans =
8.74077768365071e-16
So again, at best we can see this norm is zero. But V is not a complete set of eigenvectors. Why not? Because B is defective.
What can I say? If you read the help for eig, all it can promise is
"[V,D] = eig(A) produces a diagonal matrix D of eigenvalues and
a full matrix V whose columns are the corresponding eigenvectors
so that A*V = V*D."
When the matrix is defective, it can do no better than that, since the set of eigenvectors you are asking it to produce apparently don't exist.
0 Kommentare
Siehe auch
Kategorien
Mehr zu Matrix Indexing finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!