Add dgetrf_

Author: ptoffy

README.md

@@ -39,7 +39,7 @@ The package is structured as follows:
     - `LAPACK` - LAPACK functions:
         - [x] `dgesv_` - Solve a system of linear equations
         - [x] `dgesvd_` - Singular Value Decomposition
-        - [ ] `dgetrf_` - LU Decomposition
+        - [x] `dgetrf_` - LU Decomposition
         - [ ] `dgetri_` - Inverse of a matrix
         - [ ] `dgeev_` - Eigenvalues and eigenvectors
         - [ ] `dpotrf_` - Cholesky decomposition

Sources/AccelerateLinux/MatrixOps/LAPACK.swift

@@ -54,4 +54,24 @@ public func dgesvd_(
     _ __lwork: UnsafeMutablePointer<__CLPK_integer>!,
     _ __info: UnsafeMutablePointer<__CLPK_integer>!
 )
+
+/// DGETRF computes an LU factorization of a general M-by-N matrix A
+/// using partial pivoting with row interchanges.
+///
+/// The factorization has the form
+///    A = P * L * U
+/// where P is a permutation matrix, L is lower triangular with unit
+/// diagonal elements (lower trapezoidal if m > n), and U is upper
+/// triangular (upper trapezoidal if m < n).
+///
+/// This is the right-looking Level 3 BLAS version of the algorithm.
+@_silgen_name("dgetrf_")
+public func dgetrf_(
+    _ __m: UnsafeMutablePointer<__CLPK_integer>!,
+    _ __n: UnsafeMutablePointer<__CLPK_integer>!,
+    _ __a: UnsafeMutablePointer<__CLPK_doublereal>!,
+    _ __lda: UnsafeMutablePointer<__CLPK_integer>!,
+    _ __ipiv: UnsafeMutablePointer<__CLPK_integer>!,
+    _ __info: UnsafeMutablePointer<__CLPK_integer>!
+) -> Int32
 #endif  // canImport(Accelerate)

Tests/AccelerateLinuxTests/MatrixTests/LAPACKTests.swift

@@ -113,4 +113,39 @@ struct LAPACKTests {
 
         // Note: if we care about orthogonality of U and V we should check those too
     }
+
+    // https://numericalalgorithmsgroup.github.io/LAPACK_Examples/examples/doc/dgetrf_example.html
+    @Test("dgetrf_")
+    func test_dgetrf_() {
+        let M = 4
+        let N = 4
+        let LDA = M
+
+        var m = __CLPK_integer(M)
+        var n = __CLPK_integer(N)
+        var lda = __CLPK_integer(LDA)
+        var info = __CLPK_integer(0)
+
+        var ipiv = [__CLPK_integer](repeating: 0, count: Int(min(m, n)))
+
+        var a: [__CLPK_doublereal] = [
+            1.80, 5.25, 1.58, -1.11,
+            2.88, -2.95, -2.69, -0.66,
+            2.05, -0.95, -2.90, -0.59,
+            -0.89, -3.80, -1.04, 0.80,
+        ]
+
+        _ = dgetrf_(&m, &n, &a, &lda, &ipiv, &info)
+
+        #expect(
+            a.map {
+                ($0 * pow(10, 2)).rounded() / pow(10, 2)
+            } == [
+                5.25, 0.34, 0.30, -0.21,
+                -2.95, 3.89, -0.46, -0.33,
+                -0.95, 2.38, -1.51, 0.00,
+                -3.8, 0.41, 0.29, 0.13,
+            ]
+        )
+    }
 }