Matrix inversion in Swift using Accelerate Framework

Question

Following the good instructions that I found here: https://github.com/haginile/SwiftAccelerate I verified that matrix inversion works. In fact it did for the example given. But I get a EXC_BAD_ACCESS error for any other matrix (bigger than 2x2) for example the following 2D matrix (converted as a 1D array) has been tested in matlab and python successfully and it does not work

m = [0.55481645013013, -1.15522603580724, 0.962090414322894, -0.530226035807236, 0.168545207161447, -0.38627124296868, 0.93401699437494, -0.999999999999995, 0.684016994374945, -0.23176274578121, 0.123606797749979, -0.323606797749979, 0.432893622827287, -0.323606797749979, 0.123606797749979, 0.231762745781211, -0.684016994374948, 1.0, -0.934016994374947, 0.386271242968684, 0.168545207161448, -0.530226035807237, 0.962090414322895, -1.15522603580724, 0.554816450130132]

Its inverted matrix should be

inv(AA)

ans =

  Columns 1 through 3

          -262796763616197          -656991909040516          4.90007819375216
          -162417332048282          -406043330120712          14.6405748712708
         0.718958226823704          7.87760147961979          30.4010295628018
           162417332048287           406043330120730          46.1614842543337
           262796763616208           656991909040536          55.9019809318537

  Columns 4 through 5

          -656991909040528           262796763616211
          -406043330120721           162417332048287
         -4.28281034550088        -0.718958226823794
           406043330120704          -162417332048283
           656991909040497          -262796763616196

Could you please give me another way of matrix inversion in Swift? Or explain me how to fix this? I really don't understand why it does not work.

score 10 · Accepted Answer · edited Oct 11 '17 at 18:30

10

It doesn't work because the instructions that you found are not so good. Specifically, both pivots and workspace need to be Arrays, not scalar values; it was only working for two-by-two matrices by random chance.

Here's a modified version of the invert function that allocates the workspaces correctly:

func invert(matrix : [Double]) -> [Double] {
  var inMatrix = matrix
  var N = __CLPK_integer(sqrt(Double(matrix.count)))
  var pivots = [__CLPK_integer](count: Int(N), repeatedValue: 0)
  var workspace = [Double](count: Int(N), repeatedValue: 0.0)
  var error : __CLPK_integer = 0
  dgetrf_(&N, &N, &inMatrix, &N, &pivots, &error)
  dgetri_(&N, &inMatrix, &N, &pivots, &workspace, &N, &error)
  return inMatrix
}

I should also note that your 5x5 matrix is extremely ill-conditioned, so even when you can compute the "inverse" the error of that computation will be very large, and the inverse really shouldn't be used.

A Swift 4 version:

func invert(matrix : [Double]) -> [Double] {
    var inMatrix = matrix
    var N = __CLPK_integer(sqrt(Double(matrix.count)))
    var pivots = [__CLPK_integer](repeating: 0, count: Int(N))
    var workspace = [Double](repeating: 0.0, count: Int(N))
    var error : __CLPK_integer = 0

    withUnsafeMutablePointer(to: &N) {
        dgetrf_($0, $0, &inMatrix, $0, &pivots, &error)
        dgetri_($0, &inMatrix, $0, &pivots, &workspace, $0, &error)
    }
    return inMatrix
}

edited Oct 11 '17 at 18:30

Martin R

529,903
94
1,240
1,382

answered Nov 08 '14 at 00:10

Stephen Canon

103,815
19
183
269

Genius! Thank you so much! My matrix is ill-conditioned since it is small 5x5. I'll do other tests, but I need to use at least 21x21 matrices for my computations. – Nicholas Nov 08 '14 at 00:13
Xcode 9 fails with `Simultaneous accesses to var 'N', but modification requires exclusive access; consider copying to a local variable`. Is this a Xcode 9 bug or should we write things differently – Guig Jun 19 '17 at 22:43
1

@Guig: that's a consequence of https://github.com/apple/swift-evolution/blob/master/proposals/0176-enforce-exclusive-access-to-memory.md, though it's not obvious that it's intended to apply to C-style unsafe pointer APIs like LAPACK. You can workaround it by just making a couple copies of `N`. – Stephen Canon Jun 25 '17 at 00:50
Thanks for the link. None of the `&N` is expected to be modified by `dgetrf_` or `dgetri_` right? – Guig Jun 25 '17 at 02:36
@Guig: The LAPACK library declares those parameters as `__CLPK_integer *__n` and not as `const __CLPK_integer *__n`, even though they will not be modified. That might have historical reasons, but as a consequence, they are imported to Swift as `UnsafeMutablePointer<...>`, and the Swift compiler complains about possible simultaneous access. – Martin R Oct 11 '17 at 18:23
4

@StephenCanon: I have taken the liberty to add a Swift 4 version, hope you don't mind. – Martin R Oct 11 '17 at 18:25

Scott · Answer 2 · 2016-03-01T17:26:50.217

2

I have written a library for linear algebra in Swift. I call this library swix and it includes functions to invert matrices (this function is called inv).

Example use case:

var b = ones(10)
var A = rand((10, 10))
var AI = inv(A)
var x = AI.dot(b)

Source: https://github.com/stsievert/swix

Documentation: http://scottsievert.com/swix/

edited Mar 01 '16 at 17:26

answered Feb 29 '16 at 21:03

Scott

2,568
1
27
39

where can I get this library? – Nicholas Mar 01 '16 at 15:11
Oh gosh, I added links for source and documentation – Scott Mar 01 '16 at 17:27
genius! The link for the comparison between Swift and Matlab/Python does not work on github can you repair? I'll soon try your library for my purposes – Nicholas Mar 01 '16 at 20:39
your sample project works very well, however it is not possible to include your library in a brand new project! Your instructions are incomplete! – Nicholas Mar 16 '16 at 22:34
Another user had the same issue, and I have updated the [install instructions](http://scottsievert.com/swix/installation.html) through github [issue #24](https://github.com/stsievert/swix/issues/24). These are slight modifications; this summer, I'll revamp for the Swift Package Manager – Scott Mar 19 '16 at 21:03

Matrix inversion in Swift using Accelerate Framework

2 Answers2