create random diagonalisable matrix
11 visualizaciones (últimos 30 días)
Mostrar comentarios más antiguos
Gary Soh
el 18 de Sept. de 2018
Comentada: David Goodmanson
el 19 de Sept. de 2018
hi.. I would like to create a random diagonalisable integer matrix. Is there any code for that? thereafter I would want to create matrix X such that each the columns represent the eigenvectors.
0 comentarios
Respuesta aceptada
David Goodmanson
el 19 de Sept. de 2018
Hi Gary,
another way:
n = 7 % A is nxn
m = 9 % random integers from 1 to m
X = randi(m,n,n)
D = round(det(X))
lam = 1:n % some vector of unique integer eigenvalues, all nonzero
lamD = lam*D % final eigenvalues
A = round(X*diag(lamD)/X)
A*X - X*diag(lamD) % check
If n is too large and m is too small, this doesn't work sometimes because X comes up as a singular matrix.
3 comentarios
Bruno Luong
el 19 de Sept. de 2018
Editada: Bruno Luong
el 19 de Sept. de 2018
Actually there is no problem of lam to have null element(s). One can also select it randomly in the above code if the spectral probability is matter.
p = 5; % eg
lam = randi(p,1,n)
Más respuestas (2)
Bruno Luong
el 18 de Sept. de 2018
Editada: Bruno Luong
el 19 de Sept. de 2018
Code for both A and X are integer.
I edit the 1st version of the code (if you happens to see t) essentially a bug correction and better generation and simplification. Second edit: fix issue with non-simple eigen-value.
% Building A random (n x n) integer matrix
% and X (n x n) integer eigen-matrix of A
% meaning A*X = diag(lambda)*X
n = 4;
m = 5;
p = 5;
d = randi(2*m+1,[1,n])-m-1;
C = diag(d);
while true
P = randi(2*p+1,[n,n])-p-1;
detP = round(det(P));
if detP ~= 0
break
end
end
Q = round(detP * inv(P));
A = P*C*Q;
g = 0;
for i=1:n*n
g = gcd(g,abs(A(i)));
end
A = A/g;
lambda = sort(d)*(detP/g);
I = eye(n);
X = zeros(n);
s = 0;
for k=1:n
Ak = A-lambda(k)*I;
r = rank(Ak);
[~,~,E] = qr(Ak);
[p,~] = find(E);
j1 = p(r+1:end);
j2 = p(1:r);
[~,~,E] = qr(Ak(:,j2)');
[p,~] = find(E);
i1 = p(r+1:end);
i2 = p(1:r);
Asub = Ak(i2,j2);
s = mod(s,length(j1))+1;
x = Ak(:,j2) \ Ak(:,j1(s));
y = zeros(n-r,1);
y(s) = -1;
x = round([x; y]*det(Asub));
g = 0;
for i=1:n
g = gcd(g,abs(x(i)));
end
X([j2;j1],k) = x/g;
end
D = diag(lambda);
A
X
% % Verification A*X = X*D
A*X
X*D
0 comentarios
Matt J
el 18 de Sept. de 2018
Editada: Matt J
el 18 de Sept. de 2018
How about this,
A=randi(m,n);
A=A+A.';
[X,~]=eig(A,'vector');
4 comentarios
Matt J
el 18 de Sept. de 2018
Editada: Matt J
el 18 de Sept. de 2018
I don't think the problem is specified well enough. Eigenvectors are always unique only up to a scale factor and, in finite precision computer math, can always be made integer if you multiply them by a large enough scaling constant.
Ver también
Categorías
Más información sobre Creating and Concatenating Matrices en Help Center y File Exchange.
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!