create random diagonalisable matrix

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.

Answer 1

[David Goodmanson]

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.

[Bruno Luong]

Great idea to generate X first

[Bruno Luong]

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)

[David Goodmanson]

spectral variation does seem like a good idea.

Answer 2

[Bruno Luong]

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

Answer 3

[Matt J]

How about this,

A=randi(m,n);

A=A+A.';

 [X,~]=eig(A,'vector');

[Gary Soh]

hmm.. I was thinking matrix X to be an integer matrix

[Matt J]

Both X and A?

[Gary Soh]

yes

[Matt J]

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.