Let $\hat{e_{\mu }}$ and $\hat{e^{\mu }}$ be the co- and contravariant basis vectors, respectively, for an arbitrary coordinate system Is it true that sometimes, $\hat{e_\mu } \cdot \hat{e^\nu } \neq \delta _{\mu} ^{\nu}$? I seem to have found a counterexample:
Try choosing covariant basis $\hat{e_0 }$=(1, 1) $\hat{e_1 }$=(1, 3).
Then $g_{\mu \nu} = \hat{e_\mu } \cdot \hat{e_\nu }=\pmatrix{2 & 4 \\ 4& 10}$
The contravariant basis is the solution to set of linear equations $g_{\mu \nu} \hat{e^\nu }=\hat{e_\mu }$. Wolfram gives
$\hat{e^0 }=$(-1/2, 3/2)
$\hat{e^1 }=$(-5/12, 13/2)
These do not satisfy the identity $\hat{e_\mu } \cdot \hat{e^\nu } = \delta _{\mu} ^{\nu}$.
Your calculation of $e^0$ and $e^1$ is wrong. The easiest way to solve is to find the matrix of the inverse metric tensor:
$$(g^{ij}) = \pmatrix{2 & 4 \\ 4 & 10}^{-1} = \frac{1}{\det(g_{ij})}\pmatrix{10 & -4 \\ -4 & 2} = \pmatrix{\frac 52 & -1 \\ -1 & \frac 12}$$
Then
$$(e^0)^T = g^{0j}e_j = \frac 52\pmatrix{1 \\ 1} -1\pmatrix{1 \\ 3} = \pmatrix{\frac 32 \\ -\frac12} \\ (e^1)^T = g^{1j}e_j = -1\pmatrix{1 \\ 1} +\frac 12\pmatrix{1 \\ 3} = \pmatrix{-\frac 12 \\ \frac12}$$
Now with the correct dual vectors we can test the biorthonormality condition:
$$e^0\cdot e_0 = \pmatrix{\frac 32 & -\frac 12}\pmatrix{1 \\ 1} = 1 \\ e^0 \cdot e_1 = \pmatrix{\frac 32 & -\frac 12}\pmatrix{1 \\ 3} = 0 \\ e^1\cdot e_0 = \pmatrix{-\frac 12 & \frac 12}\pmatrix{1 \\ 1} = 0 \\ e^1 \cdot e_1 = \pmatrix{-\frac 12 & \frac 12}\pmatrix{1 \\ 3} = 1$$
So we get $e^i\cdot e_j = {\delta^i}_j$ as expected.