I'm using LAPACK via Accelerate to invert a matrix. In the code below, the matrix is represented as a Matrix struct and the inversion is performed by the inverseMatrix function. A precondition is used to ensure the input matrix is square. The error variable stores the INFO output from the dgetrf and dgetri LAPACK functions.
The error variable indicates if the matrix inversion is successful or not. The error may also indicate if the matrix is ill-conditioned therefore the inverse should not be used. Right now I'm just printing out the integers that represent these errors. Is there a better way to handle these errors and convey to the user what the error actually means?
Also, the matrix inversion may fail so how should I handle that failure? Should I make the return value optional and return nil if it fails? Or is there some other way I should handle it?
Oh, one more thing, I build/run the code with the ACCELERATE_NEW_LAPACK=1 and ACCELERATE_LAPACK_ILP64=1.
import Accelerate
struct Matrix<T> {
let rows: Int
let columns: Int
var data: [T]
func copy() -> Matrix {
Matrix(rows: self.rows, columns: self.columns, data: self.data)
}
}
func inverseMatrix(_ a: Matrix<Double>) -> Matrix<Double> {
precondition(a.rows == a.columns, "Matrix must be square")
var mat = a.copy()
var n = a.rows
var pivots = [Int](repeating: 0, count: a.rows)
var workspace = [Double](repeating: 0.0, count: a.rows)
var error = 0
// Perform LU decomposition
dgetrf_(&n, &n, &mat.data, &n, &pivots, &error)
if error != 0 {
print("LU decomposition error: \(error)")
error = 0
}
// Calculate inverse based on LU decomposition
dgetri_(&n, &mat.data, &n, &pivots, &workspace, &n, &error)
if error != 0 {
print("Inversion error: \(error)")
}
return mat
}
// --- Examples
let A = Matrix<Double>(rows: 2, columns: 2, data: [1, 2, 3, 4])
let inverseA = inverseMatrix(A)
print(inverseA, "\n")
let B = Matrix<Double>(rows: 3, columns: 3, data: [1, 2, 3, 4, 5, 6, 7, 8, 9])
let inverseB = inverseMatrix(B)
print(inverseB, "\n")