Unitary Transformation of the HF equations I am working on a code to preform a HF procedure and I am a bit confused about the unitary transformation of the basis. My instructions are as follows. 


specify the basis and geometry (done)
Evaluate and store overlap, Kinetic energy, nuclear attraction, and electron repulsion integrals (done)
diagonalize $ \mathbf{S}$ and find $\mathbf{S}^{-1/2}$ which is the unitary transformation matrix. 
Use the core hamiltonian to find the initial Fock matrix:
$$ \mathbf{S^{-1/2}}^{T} \mathbf{H}~\mathbf{S^{-1/2}}\equiv\mathbf{F}_0$$
Diagonalize the Fock matrix to find the mo coefficients:
$$\mathbf{F}_0 \mathbf{C}^{'}_0=\mathbf{C}^{'}_0\epsilon_0$$
Back transform the coefficients.
$$\mathbf{C^{'}_0}~\mathbf{S^{-1/2}}=\mathbf{C}_0$$
compute the density, use the density to evaluate new fock matrix, find new coefficients repeat until convergeance.


Now step 7 kinds glazes over a lot of the work but I understand that part, so It isn't important to the question.

I am confused by the $\mathbf{S}^{-1/2}$ matrix.  Is this literally one over the square root of each element of $\mathbf{S}$? Or is there some decomposition I should preform that I just do not understand.  I can carry out the whole procedure but without orthonormalization my results are nonsense.
You are going to want to check out Szabo and Ostlund page 143, they explain all the nitty-gritty for this matrix. The math itself is on page 22 and 23. But if you don't have a copy...

The $\mathbf{S}^{-1/2}$ comes from the need to orthogonalize your basis so that you get your equation (5) instead of the $\mathbf{FC}=\epsilon\mathbf{SC}$. It's one method of several to achieve this orthogonalization, and it's called symmetric orthogonalization. 

So how do you compute this function of a Hermitian matrix? Well, functions on diagonal matrices are super easy; you just apply the function to each diagonal element. For non-diagonal matrices, this is not true. More to your question, you can't just take one over the square root of each element of $\mathbf{S}$. But you can if it's diagonal. And you can make any Hermitian matrix diagonal via a unitary transformation.

The trick, then, is to diagonalize $\mathbf{S}$ first, apply the inverse square root to each eigenvalue, then transform back using the eigenvectors. So you're going to want to


Solve $\mathbf{SU} = \mathbf{sU}$. Where $\mathbf{s}$ are the eigenvalues of $\mathbf{S}$, and $\mathbf{U}$ are the eigenvectors.
For each $s$ in $\mathbf{s}$ do $s^{-1/2}$.
Transform back to original basis by $\mathbf{S}^{-1/2} = \mathbf{Us}^{-1/2}\mathbf{U}^{\dagger}$