The matrices you are looking to invert are not "valid" covariances matrices because they are not positive definite; numerically they even have some eigenvalues that are negative (but close to zero). This most probably due to machine zeros, for example the last eigenvalue of your "matlab_covariance" matrix is -0.000000016313723. To correct to positive definite you can do two things :

1. Just add some very small constant along the diagonal; a form of Tikhonov correction really ($Χ+\lambda I$ with $\lambda \rightarrow 0$).
2. (Based on what whuber proposed) Use SVD, set the "problematic" eigenvalues to a fixed small value (not zero), reconstruct your covariance matrix and then invert that. Clearly if you set some of those eigenvalues to zero you will end up with a non-negative (or semi-positive) matrix, that will not be invertible still.

A non-negative matrix does not have a inverse but it does have a pseudo inverse (all matrices with real or complex entries have a pseudo-inverse), nevertheless the Moore–Penrose pseudo-inverse is more computationally expensive than a true inverse and if the inverse exists it is equal to the pseudo-inverse. So just go for the inverse :)

Both methods practically try to handle the eigenvalues that evaluate to zero (or below zero). The first method is bit hand-wavy but probably much faster to implement. For something slightly more stable you might want to compute the SVD and then set the $\lambda$ equal to absolute of the smallest eigenvalue (so you get non-negative) plus something very small (so you get positive). Just be careful not to enforce positivity to a matrix that is obviously negative (or already positive). Both methods will alter the conditioning number of your matrix.

In statistical terms what you do by adding that $\lambda$ across the diagonal of your covariance matrix you add noise to your measurements. (Because the diagonal of the covariance matrix is the variance of each point and by adding something to those values you just say "the variance at the points I have readings for is actually a bit bigger than what I thought originally".)

A fast test for the positive-definiteness of a matrix is the existence (or not) of the Cholesky decomposition of it.

Also as a computational note:

1. Don't use floats, use doubles to store your covariance matrix.
2. Use numerical linear algebra libraries in C++ (like Eigen or Armadillo) to get inverses of matrices, matrix products, etc. It's faster, safer and more concise.

The matrices you are looking to invert are not "valid" covariances matrices because they are not positive definite; numerically they even have some eigenvalues that are negative (but close to zero). This most probably due to machine zeros, for example the last eigenvalue of your "matlab_covariance" matrix is -0.000000016313723. To correct to positive definite you can do two things :

1. Just add some very small constant along the diagonal; a form of Tikhonov correction really ($Χ+\lambda I$ with $\lambda \rightarrow 0$).
2. (Based on what whuber proposed) Use SVD, set the "problematic" eigenvalues to a fixed small value (not zero), reconstruct your covariance matrix and then invert that. Clearly if you set some of those eigenvalues to zero you will end up with a non-negative (or semi-positive) matrix, that will not be invertible still.

A non-negative matrix does not have a inverse but it does have a pseudo inverse (all matrices with real or complex entries have a pseudo-inverse), nevertheless the Moore–Penrose pseudo-inverse is more computationally expensive than a true inverse and if the inverse exists it is equal to the pseudo-inverse. So just go for the inverse :)

Both methods practically try to handle the eigenvalues that evaluate to zero (or below zero). The first method is bit hand-wavy but probably much faster to implement. For something slightly more stable you might want to compute the SVD and then set the $\lambda$ equal to absolute of the smallest eigenvalue (so you get non-negative) plus something very small (so you get positive). Just be careful not to enforce positivity to a matrix that is obviously negative (or already positive). Both methods will alter the conditioning number of your matrix.

In statistical terms what you do by adding that $\lambda$ across the diagonal of your covariance matrix you add noise to your measurements. (Because the diagonal of the covariance matrix is the variance of each point and by adding something to those values you just say "the variance at the points I have readings for is actually a bit bigger than what I thought originally".)

A fast test for the positive-definiteness of a matrix is the existence (or not) of the Cholesky decomposition of it.

Also as a computational note:

1. Don't use floats, use doubles to store your covariance matrix.
2. Use numerical linear algebra libraries in C++ (like Eigen or Armadillo) to get inverses of matrices, matrix products, etc. It's faster, safer and more concise.

EDIT: Given you have a Cholesky decomposition of your matrix $K$ such that $LL^T$ (you have to do that to check you are having a Pos.Def. matrix) you should be able to immediately solve the system $Kx=b$. You just solve Ly = b for y by forward substitution, and then L^Tx = y for x by back substitution. (In eigen just use the .solve(x) method of your Cholesky object) Thanks to bnaul and Zen for pointing out that I focused so much on getting the $K$ be Pos.Def. that I forgot why we cared about that in the first place :)