% #33, p. 326
A=[-6 4 0 9; 
    -3 0 1 6;
    -1 -2 1 0;
    -4 4 0 7]

% the eigenvalues of A
ev=eig(A)
rref(A-ev(1)*eye(4))
rref(A-ev(2)*eye(4))
rref(A-ev(3)*eye(4))

% here are bases for the null spaces:
null(A-ev(1)*eye(4))
null(A-ev(2)*eye(4))
null(A-ev(3)*eye(4))

% You can get the eigenvectors all at once, as follows:
[V,d]=eig(A)

% Now we'll check to see if they behave as we suspect they should:
(A-ev(3)*eye(4))*null(A-ev(3)*eye(4))
(A-ev(3)*eye(4))*V