[swift]返回到线性矩阵方程的最小二乘解的函数

标签: Python Swift
发布时间: 2017/3/6 23:01:05
注意事项: 本文中文内容可能为机器翻译,如要查看英文原文请点击上面连接.

我必须重写来自到 Python 代码 Swift 我却对功能应回归线性矩阵方程的最小二乘解。有谁的你知道写的库, Swift 具有一种等效方法 numpy.linalg.lstsq ?我会感谢你的帮助。

Python 代码︰

a = numpy.array([[p2.x-p1.x,p2.y-p1.y],[p4.x-p3.x,p4.y-p3.y],[p4.x-p2.x,p4.y-p2.y],[p3.x-p1.x,p3.y-p1.y]])
b = numpy.array([number1,number2,number3,number4])
res = numpy.linalg.lstsq(a,b) 
result = [float(res[0][0]),float(res[0][1])]
return result

Swift到目前为止的代码︰

var matrix1 = [[p2.x-p1.x, p2.y-p1.y],[p4.x-p3.x, p4.y-p3.y], [p4.x-p2.x, p4.y-p2.y], [p3.x-p1.x, p3.y-p1.y]]
var matrix2 = [number1, number2, number3, number4]

解决方法 1:

加快框架包括LAPACK线性代数软件包,有一个DGELS函数来解决下-或超定线性系统。从文档︰

DGELS 解决了方程组或欠定实线性系统涉及 M 由 N 矩阵 A 或其转置,使用 A.QR 或 LQ 分解 假定 A 有满秩。

这里有一个例子,如何可以从使用该函数 Swift 。 它是本质上是一个翻译此 C 示例代码

func solveLeastSquare(A A: [[Double]], B: [Double]) -> [Double]? {
    precondition(A.count == B.count, "Non-matching dimensions")

    var mode = Int8(bitPattern: UInt8(ascii: "N")) // "Normal" mode
    var nrows = CInt(A.count)
    var ncols = CInt(A[0].count)
    var nrhs = CInt(1)
    var ldb = max(nrows, ncols)

    // Flattened columns of matrix A
    var localA = (0 ..< nrows * ncols).map {
        A[Int($0 % nrows)][Int($0 / nrows)]
    }

    // Vector B, expanded by zeros if ncols > nrows
    var localB = B
    if ldb > nrows {
        localB.appendContentsOf([Double](count: ldb - nrows, repeatedValue: 0.0))
    }

    var wkopt = 0.0
    var lwork: CInt = -1
    var info: CInt = 0

    // First call to determine optimal workspace size
    dgels_(&mode, &nrows, &ncols, &nrhs, &localA, &nrows, &localB, &ldb, &wkopt, &lwork, &info)
    lwork = Int32(wkopt)

    // Allocate workspace and do actual calculation
    var work = [Double](count: Int(lwork), repeatedValue: 0.0)
    dgels_(&mode, &nrows, &ncols, &nrhs, &localA, &nrows, &localB, &ldb, &work, &lwork, &info)

    if info != 0 {
        print("A does not have full rank; the least squares solution could not be computed.")
        return nil
    }
    return Array(localB.prefix(Int(ncols)))
}

一些注意事项︰

  • dgels_()修改数据传递的矩阵和向量,并期望矩阵作为"平面"的数组,其中包含的列 A 。 右侧也预计作为数组长度 max(M, N) 。 为此,输入的数据复制到本地变量第一次。
  • 所有参数都必须都通过引用向 dgels_() ,这就是为什么它们被存储在 var s。
  • C 整数是 32 位整数,这使得一些转换之间 IntCInt 有必要。

示例 1:http://www.seas.ucla.edu/~vandenbe/103/lectures/ls.pdf的超定的系统。

let A = [[ 2.0, 0.0 ],
         [ -1.0, 1.0 ],
         [ 0.0, 2.0 ]]
let B = [ 1.0, 0.0, -1.0 ]
if let x = solveLeastSquare(A: A, B: B) {
    print(x) // [0.33333333333333326, -0.33333333333333343]
}

示例 2:欠定的系统、 极小范数解到 x_1 + x_2 + x_3 = 1.0

let A = [[ 1.0, 1.0, 1.0 ]]
let B = [ 1.0 ]
if let x = solveLeastSquare(A: A, B: B) {
    print(x) // [0.33333333333333337, 0.33333333333333337, 0.33333333333333337]
}
赞助商