In [2]: 
from sympy import *
from IPython.display import display
init_printing(use_latex='mathjax')
In [3]: 
J=Matrix([[-3,0,0,0,0,0,0,0],
          [0,0,0,0,0,0,0,0],
          [0,0,2,0,0,0,0,0],
          [0,0,0,4,1,0,0,0],
          [0,0,0,0,4,0,0,0],
          [0,0,0,0,0,4,1,0],
          [0,0,0,0,0,0,4,1],
          [0,0,0,0,0,0,0,4]])
H=Matrix([[1,2,3,4,5,6,7,8],
          [0,1,2,3,4,5,6,7],
          [0,0,1,2,3,4,5,6],
          [0,0,0,1,2,3,4,5],
          [0,0,0,0,1,2,3,4],
          [0,0,0,0,0,1,2,3],
          [0,0,0,0,0,0,1,2],
          [1,0,0,0,0,0,0,1]])
Hi=H.inv()
A=H*J*Hi
display(A)
$$\left[\begin{matrix}2 & -4 & 8 & 4 & -2 & -8 & 10 & -5\\4 & -8 & 8 & 4 & -1 & -6 & 8 & -4\\3 & -6 & 5 & 4 & 0 & -4 & 6 & -3\\2 & -4 & 2 & 4 & 1 & -2 & 4 & -2\\1 & -2 & 1 & 0 & 4 & 0 & 2 & -1\\0 & 0 & 0 & 0 & 0 & 4 & 1 & 0\\-1 & 2 & -1 & 0 & 0 & 0 & 4 & 1\\-7 & 14 & -7 & 0 & 0 & 0 & 0 & 4\end{matrix}\right]$$
In [4]: 
lmbda=var('lambda')
I=Matrix.eye(8)
B=(A-lmbda*I)
display(B)
d=B.det()
display(d)
$$\left[\begin{matrix}- \lambda + 2 & -4 & 8 & 4 & -2 & -8 & 10 & -5\\4 & - \lambda - 8 & 8 & 4 & -1 & -6 & 8 & -4\\3 & -6 & - \lambda + 5 & 4 & 0 & -4 & 6 & -3\\2 & -4 & 2 & - \lambda + 4 & 1 & -2 & 4 & -2\\1 & -2 & 1 & 0 & - \lambda + 4 & 0 & 2 & -1\\0 & 0 & 0 & 0 & 0 & - \lambda + 4 & 1 & 0\\-1 & 2 & -1 & 0 & 0 & 0 & - \lambda + 4 & 1\\-7 & 14 & -7 & 0 & 0 & 0 & 0 & - \lambda + 4\end{matrix}\right]$$
$$\lambda^{8} - 19 \lambda^{7} + 134 \lambda^{6} - 360 \lambda^{5} - 320 \lambda^{4} + 4096 \lambda^{3} - 8704 \lambda^{2} + 6144 \lambda$$
In [5]: 
fd=factor(d)
display(fd)
roots=solve(fd,lmbda)
display(roots)
#eigenvalues -3, 0, 2 each have algebraic multiplicity 1
#eigenvalues 4 has algebraic multiplicity 5
$$\lambda \left(\lambda - 4\right)^{5} \left(\lambda - 2\right) \left(\lambda + 3\right)$$
$$\left [ -3, \quad 0, \quad 2, \quad 4\right ]$$
In [6]: 
height = lambda M: len(M[:,0])
width = lambda M: len(M[0,:])
def kernel_basis(A):
    '''Returns matrix whose column are a basis for the null space of A.'''
    #1. n x n matrix identity matrix B.
    #2. m x n matrix R = rref of A.
    #3. Subtract each row i of R from row i of B.
    #4. For each pivot column j of R, delete column j of B.
    m=width(A)
    B=eye(m)
    tmp=Matrix(A)
    R,pivots=tmp.rref()
    cols=range(m)
    for (p,e) in zip(pivots,range(len(pivots))):
        for c in cols:
            B[p,c]=B[p,c]-R[e,c]
    for p in reversed(pivots):
        B.col_del(p)
    return(B)
In [10]: 
V={}
for r in roots:
    #find basis of kernel of A-r*I
    V[r]=(A-r*I).nullspace()
    display(V[r])

#eigenvalues -3, 0, 2 each have geometric multiplicity 1
#eigenvalues 4 has geometric multiplicity 2
$$\left [ \left[\begin{matrix}1\\0\\0\\0\\0\\0\\0\\1\end{matrix}\right]\right ]$$
$$\left [ \left[\begin{matrix}2\\1\\0\\0\\0\\0\\0\\0\end{matrix}\right]\right ]$$
$$\left [ \left[\begin{matrix}3\\2\\1\\0\\0\\0\\0\\0\end{matrix}\right]\right ]$$
$$\left [ \left[\begin{matrix}4\\3\\2\\1\\0\\0\\0\\0\end{matrix}\right], \quad \left[\begin{matrix}-6\\-4\\-2\\0\\2\\1\\0\\0\end{matrix}\right]\right ]$$
In [11]: 
(h,j)=A.jordan_form()
display([h,j])

#Since eigenvalue 4 has geometric multiplicity 2,
#the Jordan form will have 2 Jordan blocks
#for that eigenvector and 2 associated
#Jordan chains.

#A Jordan chain of length 3:
#A-4I transform the 6th column of h to the 5th column
#and transforms the 5th column to the 4th column.
#The associated Jordan block is 3x3
display((A-4*I)*h[:,6-1]==h[:,5-1])
display((A-4*I)*h[:,5-1]==h[:,4-1])
display(j[4-1:6,4-1:6])

#A Jordan chain of length 2:
#A-4I transform the 8th column of h to the 7th column
#The associated Jordan block is 2x2
display((A-4*I)*h[:,8-1]==h[:,7-1])
display(j[7-1:8,7-1:8])

#These two chain lengths sum to 5, 
#the algebraic multiplicity of eigenvalue 4.

#The eigenvalues -3, 0, 2 have trivial
#1x1 Jordan blocks, [-3], [0], [2], respectively,
#and Jordan chains of length one, namely,
#columns 1, 2, and 3 of h, respectively.
$$\left [ \left[\begin{matrix}1 & 2 & 3 & 6 & -5 & 0 & \frac{142}{5} & - \frac{1}{5}\\0 & 1 & 2 & 5 & -4 & 0 & \frac{109}{5} & 0\\0 & 0 & 1 & 4 & -3 & 0 & \frac{76}{5} & \frac{1}{5}\\0 & 0 & 0 & 3 & -2 & 0 & \frac{43}{5} & \frac{2}{5}\\0 & 0 & 0 & 2 & -1 & 0 & 2 & \frac{3}{5}\\0 & 0 & 0 & 1 & 0 & 0 & 1 & -2\\0 & 0 & 0 & 0 & 1 & 0 & 0 & 1\\1 & 0 & 0 & 0 & 0 & 1 & 0 & 0\end{matrix}\right], \quad \left[\begin{matrix}-3 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 2 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 4 & 1 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 4 & 1 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 4 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 4 & 1\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 4\end{matrix}\right]\right ]$$ 
True

True
$$\left[\begin{matrix}4 & 1 & 0\\0 & 4 & 1\\0 & 0 & 4\end{matrix}\right]$$ 
True
$$\left[\begin{matrix}4 & 1\\0 & 4\end{matrix}\right]$$
In [12]: 
Matrix([[1,5,0],[0,1,0],[0,0,1]]).jordan_form()
Out[12]:
$$\left ( \left[\begin{matrix}5 & 0 & 0\\0 & 1 & 0\\0 & 0 & 1\end{matrix}\right], \quad \left[\begin{matrix}1 & 1 & 0\\0 & 1 & 0\\0 & 0 & 1\end{matrix}\right]\right )$$
In []: