Mathematica has the function Eigensystem, which diagonalizes a matrix. Beware however that the function Eigensystem does not deliver an orthonormal set of eigenvectors, nor orders the eigenenergies. This, although mathematically formally correct behavior is extreme tedious and for most physicist / chemist totally contra intuitive behavior. So I wrote the function OrthonormalEigensystem to replace the function Eigensystem.

OrthonormalEigensystem[Opp_,n_Integer]:=Module[{maxval,val,fun,valOrder},
maxval=Chop[Eigenvalues[Opp,1][[1]]];
If[Re[maxval]>0,
maxval+=Eigenvalues[Opp-SparseArray[{{i_,i_}->maxval},{Length[Opp],Length[Opp]}],1][[1]]
];
{val,fun}=Simplify[Eigensystem[Opp-SparseArray[{{i_,i_}->maxval},{Length[Opp],Length[Opp]}],-n]];
valOrder=Ordering[val];
If[$VersionNumber>= 6,
fun=Simplify[Orthogonalize[fun]],
fun=Simplify[GramSchmidt[fun,InnerProduct->(Dot[#1, #2\[Conjugate]]&)]]
];
Simplify[{val[[valOrder]]+maxval,fun[[valOrder]]}]
];
OrthonormalEigensystem[Opp_]:=OrthonormalEigensystem[Opp,Length[Opp]];

3 by 3 graphically and numerically

For a single bond there are nice ways to graph things. It becomes a bit more tricky for three orbitals binding together. (CO2, H2O, etc.)

Here I choose to place the onsite energies of the 3 orbitals/atoms before binding at the horizontal scale on the bottom of the figure. The energies of the states after binding are represented by vertical lines.

In principle one can read this graph as a density of states. On the horizontal axes the energy of the states is plotted. On the vertical axes the number of states found in a particular energy range is shown. (delta functions for a molecule)

ClearAll[t01,t02,t12,Δ1,Δ2];
Hloc=(({{0, t01, t02}, {t01, Δ1, t12}, {t02, t12, Δ2}}));
Manipulate[{valloc,funloc}=Eigensystem[(({{0, t01, t02}, {t01, Δ1, t12}, {t02, t12, Δ2}}))];Show[Graphics[{
{RGBColor[0,0,1],
Rectangle[{valloc[[1]]-0.1,1},{valloc[[1]]-0.1+0.3((funloc[[1,1]])^2),2}],Rectangle[{valloc[[2]]-0.1,1},{valloc[[2]]-0.1+0.3((funloc[[2,1]])^2),2}],Rectangle[{valloc[[3]]-0.1,1},{valloc[[3]]-0.1+0.3((funloc[[3,1]])^2),2}]},
{RGBColor[1,0,0],
Rectangle[{valloc[[1]]-0.1+0.3((funloc[[1,1]])^2),1},{valloc[[1]]-0.1+0.3(((funloc[[1,1]])^2)+((funloc[[1,2]])^2)),2}],Rectangle[{valloc[[2]]-0.1+0.3((funloc[[2,1]])^2),1},{valloc[[2]]-0.1+0.3(((funloc[[2,1]])^2)+((funloc[[2,2]])^2)),2}],Rectangle[{valloc[[3]]-0.1+0.3((funloc[[3,1]])^2),1},{valloc[[3]]-0.1+0.3(((funloc[[3,1]])^2)+((funloc[[3,2]])^2)),2}]},
{RGBColor[0,1,0],
Rectangle[{valloc[[1]]-0.1+0.3(((funloc[[1,1]])^2)+((funloc[[1,2]])^2)),1},{valloc[[1]]-0.1+0.3,2}],Rectangle[{valloc[[2]]-0.1+0.3(((funloc[[2,1]])^2)+((funloc[[2,2]])^2)),1},{valloc[[2]]-0.1+0.3,2}],Rectangle[{valloc[[3]]-0.1+0.3(((funloc[[3,1]])^2)+((funloc[[3,2]])^2)),1},{valloc[[3]]-0.1+0.3,2}]},{Thickness[0.02],{RGBColor[0,0,1],Line[{{0,0},{0,0.5}}]},{RGBColor[1,0,0],Line[{{Δ1,0},{Δ1,0.5}}]},{RGBColor[0,1,0],Line[{{Δ2,0},{Δ2,0.5}}]}},
Text[Style[MatrixForm[Hloc],Large],{-9,3},{-1,1}],
Text[Style[MatrixForm[NumberForm[(({{0, t01, t02}, {t01, Δ1, t12}, {t02, t12, Δ2}})),{3,1},NumberPadding->{"","0"},NumberSigns->{"-","+"}]],Large],{9,3},{1,1}]}],AspectRatio->1,PlotRange->{{-9,9},{0,3}}],{{t01,0},-3.,3.},{{t02,0},-3.,3.},{{t12,0},-3.,3.},{{Δ1,-3},-5.,5},{{Δ2,3},-5.,5.}]