Generalized Inverse function in GSL

Question

Is there any function to calculate Generalized inverse of a matrix using GSL? Like in R, we have ginv(X, tol = sqrt(.Machine$double.eps)).

Answer 1

No. It seems as if there is no routine to directly calculate the pseudo-inverse of a matrix (although here you can find a discussion on how one could get it).

However, the explicit pseudo-inverse itself is seldom required. Instead, gsl provides the routine

int gsl_linalg_SV_solve (const gsl_matrix * U, const gsl_matrix * V, const gsl_vector * S
                       , const gsl_vector * b, gsl_vector * x)

see the documentation here.

It solves the linear system A x = b, which is equivalent to applying the pseudo-inverse A^+ onto b and yields x = A^+ b.

Before application, the SVD must be found via the routine gsl_linalg_SV_decomp. The tolerance factor tol you mentioned can be incorporated by looping over the singular values S and setting those smaller than tol to zero.

(Further, here is a personal suggestion: drop the gsl and switch to Eigen, armadillo, or comparable modern libraries)