Add inverse routines to lapack64

lapack.go

@@ -28,6 +28,7 @@ type Float64 interface {
 	Dgelqf(m, n int, a []float64, lda int, tau, work []float64, lwork int)
 	Dgeqrf(m, n int, a []float64, lda int, tau, work []float64, lwork int)
 	Dgetrf(m, n int, a []float64, lda int, ipiv []int) (ok bool)
+	Dgetri(n int, a []float64, lda int, ipiv []int, work []float64, lwork int) (ok bool)
 	Dgetrs(trans blas.Transpose, n, nrhs int, a []float64, lda int, ipiv []int, b []float64, ldb int)
 	Dlantr(norm MatrixNorm, uplo blas.Uplo, diag blas.Diag, m, n int, a []float64, lda int, work []float64) float64
 	Dlange(norm MatrixNorm, m, n int, a []float64, lda int, work []float64) float64
@@ -37,6 +38,7 @@ type Float64 interface {
 	Dpocon(uplo blas.Uplo, n int, a []float64, lda int, anorm float64, work []float64, iwork []int) float64
 	Dpotrf(ul blas.Uplo, n int, a []float64, lda int) (ok bool)
 	Dtrcon(norm MatrixNorm, uplo blas.Uplo, diag blas.Diag, n int, a []float64, lda int, work []float64, iwork []int) float64
+	Dtrtri(uplo blas.Uplo, diag blas.Diag, n int, a []float64, lda int) (ok bool)
 	Dtrtrs(uplo blas.Uplo, trans blas.Transpose, diag blas.Diag, n, nrhs int, a []float64, lda int, b []float64, ldb int) (ok bool)
 }
 

lapack64/lapack64.go

@@ -163,6 +163,22 @@ func Getrf(a blas64.General, ipiv []int) bool {
 	return lapack64.Dgetrf(a.Rows, a.Cols, a.Data, a.Stride, ipiv)
 }
 
+// Getri computes the inverse of the matrix A using the LU factorization computed
+// by Getrf. On entry, a contains the PLU decomposition of A as computed by
+// Getrf and on exit contains the reciprocal of the original matrix.
+//
+// Getri will not perform the inversion if the matrix is singular, and returns
+// a boolean indicating whether the inversion was successful.
+//
+// Work is temporary storage, and lwork specifies the usable memory length.
+// At minimum, lwork >= n and this function will panic otherwise.
+// Dgetri is a blocked inversion, but the block size is limited
+// by the temporary space available. If lwork == -1, instead of performing Getri,
+// the optimal work length will be stored into work[0].
+func Getri(a blas64.General, ipiv []int, work []float64, lwork int) (ok bool) {
+	return lapack64.Dgetri(a.Cols, a.Data, a.Stride, ipiv, work, lwork)
+}
+
 // Dgetrs solves a system of equations using an LU factorization.
 // The system of equations solved is
 //  A * X = B if trans == blas.Trans
@@ -272,6 +288,15 @@ func Trcon(norm lapack.MatrixNorm, a blas64.Triangular, work []float64, iwork []
 	return lapack64.Dtrcon(norm, a.Uplo, a.Diag, a.N, a.Data, a.Stride, work, iwork)
 }
 
+// Trtri computes the inverse of a triangular matrix, storing the result in place
+// into a.
+//
+// Trtri will not perform the inversion if the matrix is singular, and returns
+// a boolean indicating whether the inversion was successful.
+func Trtri(a blas64.Triangular) (ok bool) {
+	return lapack64.Dtrtri(a.Uplo, a.Diag, a.N, a.Data, a.Stride)
+}
+
 // Trtrs solves a triangular system of the form A * X = B or A^T * X = B. Trtrs
 // returns whether the solve completed successfully. If A is singular, no solve is performed.
 func Trtrs(trans blas.Transpose, a blas64.Triangular, b blas64.General) (ok bool) {