The algorithm for finding the eigenvalues and eigenvectors doesn't, in fact, use the determinant. In fact, it doesn't use Gaussian Elimination like the person above says either.
It uses the QR decomposition. There is a series of Matlab codes here showing some examples. One is eigvecs and geometric multiplicity.
function [S, D] = eigvec(A)
% eigvec Eigenvectors and their geometric multiplicity.
% S = eigvec(A) returns the largest possible set of linearly
% independent eigenvectors of A.
% [S, D] = eigvec(A) also returns the corresponding eigenvalues
% in the diagonal matrix D.
% Each eigenvalue in D is repeated according to the number of its
% linearly independent eigenvectors. This is its geometric multiplicity.
% Always A*S = S*D. If S is square then A is diagonalizable and
% inv(S)*A*S = D = LAMBDA.
[m, n] = size(A);
I = eye(n);
[evalues, repeats] = eigval(A);
S = ; d = ;
for k = 1 : length(evalues);
s = nulbasis(A - evalues(k)*I);
[ms, ns] = size(s);
S = [S s];
temp = ones(ns, 1) * evalues(k);
d = [d; temp];
D = diag(d);
function [evalues, repeats] = eigval(A)
% eigval Eigenvalues and their algebraic multiplicity.
% evalues = eigval(A) returns the distinct eigenvalues of A,
% with duplicates removed and sorted in decreasing order.
% [evalues, repeats] = eigval(A) also returns the row vector
% repeats that gives the multiplicity of each eigenvalue.
% The sum of the multiplicities is n.
% Examples: Let A = eye(n) and B = diag([3 4]).
% For A, evalues is 1 and repeats is n.
% For B, evalues is [4; 3] and repeats is [1 1].
tol = sqrt(eps);
lambda = sort(eig(A));
lambda = round(lambda/tol) * tol;
% lambda gives all n eigenvalues (repetitions included).
evalues = unique(lambda);
evalues = flipud(evalues);
n = length(lambda);
d = length(evalues);
A = ones(n, 1) * evalues';
B = lambda * ones(1, d);
MATCH = abs(A-B) <= tol;
% MATCH is an n by d zero matrix except
% MATCH(i,j) = 1 when lambda(i) = evalues(j).
% Summing the columns gives the row vector repeats.
repeats = sum(MATCH);
function N = nulbasis(A)
% nulbasis Basis for nullspace.
% N = nulbasis(A) returns a basis for the nullspace of A
% in the columns of N. The basis contains the n-r special
% solutions to Ax=0. freecol is the list of free columns.
% >> A = [1 2 0 3;
% [0 0 1 4];
% >> N = nulbasis(A)
% N = [-2 -3]
% [ 1 0]
% [ 0 -4]
% [ 0 1]
% See also fourbase.
[R, pivcol] = rref(A, sqrt(eps));
[m, n] = size(A);
r = length(pivcol);
freecol = 1:n;
freecol(pivcol) = ;
N = zeros(n, n-r);
N(freecol, : ) = eye(n-r);
N(pivcol, : ) = -R(1:r, freecol);
A = [2,0,0;0,2,0;0,1,1];
[evalues, repeats] = eigval(A)
[S, D] = eigvec(A)
1 0 0
0 1 0
0 1 1
2 0 0
0 2 0
0 0 1
comparing to matlab
>> [V,D] = eigs(A)
0 1.0000 0
0 0 0.7071
1.0000 0 0.7071
1 0 0
0 2 0
0 0 2
If you're asking to which eigenvector does each eigen value belong. It goes by order I believe. However the specifics are in Netlib.