Проблемы с фреймворком Accelerate в Swift

Я использую алгоритм dgeev от LAPACK в структуре Accelerate для вычисления собственных значений и собственных векторов матрицы. Вот мой код:

    var matrix:[__CLPK_doublereal] = [1,2,3,4,5,6,7,8,9]

    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
    var lwork = __CLPK_integer(-1)
    // Real parts of eigenvalues
    var wr = [Double](count: Int(N), repeatedValue: 0)
    // Imaginary parts of eigenvalues
    var wi = [Double](count: Int(N), repeatedValue: 0)
    // Left eigenvectors
    var vl = [__CLPK_doublereal](count: Int(N*N), repeatedValue: 0)
    // Right eigenvectors
    var vr = [__CLPK_doublereal](count: Int(N*N), repeatedValue: 0)

    dgeev_(UnsafeMutablePointer(("V" as NSString).UTF8String), UnsafeMutablePointer(("V" as NSString).UTF8String), &N, &matrix, &N, &wr, &wi, &vl, &N, &vr, &N, &workspace, &lwork, &error)

    println("\(wr), \(vl), \(vr)")

Это печатает только массивы, содержащие нули, что означает, что они не изменяются функцией. Что я делаю неправильно?

ОБНОВЛЕНИЕ 1

Теперь у меня есть это:

    var matrix:[__CLPK_doublereal] = [1,2,3,4,5,6,7,8,9]

    var N = __CLPK_integer(sqrt(Double(matrix.count)))
    var pivots = [__CLPK_integer](count: Int(N), repeatedValue: 0)
    var workspaceQuery = [Double](count: 1, repeatedValue: 0.0)
    var error : __CLPK_integer = 0
    var lwork = __CLPK_integer(-1)
    // Real parts of eigenvalues
    var wr = [Double](count: Int(N), repeatedValue: 0)
    // Imaginary parts of eigenvalues
    var wi = [Double](count: Int(N), repeatedValue: 0)
    // Left eigenvectors
    var vl = [__CLPK_doublereal](count: Int(N*N), repeatedValue: 0)
    // Right eigenvectors
    var vr = [__CLPK_doublereal](count: Int(N*N), repeatedValue: 0)

    dgeev_(UnsafeMutablePointer(("V" as NSString).UTF8String), UnsafeMutablePointer(("V" as NSString).UTF8String), &N, &matrix, &N, &wr, &wi, &vl, &N, &vr, &N, &workspaceQuery, &lwork, &error)
    var workspace = [Double](count: Int(workspaceQuery[0]), repeatedValue: 0.0)

    dgeev_(UnsafeMutablePointer(("V" as NSString).UTF8String), UnsafeMutablePointer(("V" as NSString).UTF8String), &N, &matrix, &N, &wr, &wi, &vl, &N, &vr, &N, &workspace, &lwork, &error)

    println("\(wr), \(vl), \(vr)")

Тем не менее он печатает нули.


swift matrix lapack eigenvector accelerate-framework
person Youssef Moawad    schedule 11.01.2015    source источник

Ответы (1)


arrow_upward
3
arrow_downward

Проблема с вашей переменной lwork. Предполагается, что это размер рабочей области, которую вы предоставляете, где -1 означает, что вы выполняете «запрос рабочей области»:

 LWORK   (input) INTEGER   
        The dimension of the array WORK.  LWORK >= max(1,3*N), and   
        if JOBVL = 'V' or JOBVR = 'V', LWORK >= 4*N.  For good   
        performance, LWORK must generally be larger.   

        If LWORK = -1, then a workspace query is assumed; the routine   
        only calculates the optimal size of the WORK array, returns   
        this value as the first entry of the WORK array, and no error   
        message related to LWORK is issued by XERBLA.

Итак, вы, вероятно, хотите что-то вроде этого:

var workspaceQuery: Double = 0.0
dgeev_(UnsafeMutablePointer(("V" as NSString).UTF8String), UnsafeMutablePointer(("V" as NSString).UTF8String), &N, &matrix, &N, &wr, &wi, &vl, &N, &vr, &N, &workspaceQuery, &lwork, &error)

// prints "102.0"
println("\(workspaceQuery)")

// size workspace per the results of the query:
var workspace = [Double](count: Int(workspaceQuery), repeatedValue: 0.0)
lwork = __CLPK_integer(workspaceQuery)

dgeev_(UnsafeMutablePointer(("V" as NSString).UTF8String), UnsafeMutablePointer(("V" as NSString).UTF8String), &N, &matrix, &N, &wr, &wi, &vl, &N, &vr, &N, &workspace, &lwork, &error)

// this now prints non-zero values
println("\(wr), \(vl), \(vr)")
person Airspeed Velocity    schedule 11.01.2015
comment
Хорошо, а каким должно быть значение lwork при втором вызове? - person Youssef Moawad; 11.01.2015
comment
@YoussefSami: я думаю, это должен быть реальный размер рабочего пространства: lwork = __CLPK_integer(workspace.count). Я нашел этот пример кода C, который также демонстрирует использование: software.intel.com/sites/products/documentation/doclib/mkl_sa/ - person Martin R; 11.01.2015
comment
Да, это ошибка в моем примере кода, я должен установить ее, а также изменить размер рабочей области. Будет редактировать. - person Airspeed Velocity; 11.01.2015
comment
Огромное спасибо вам обоим! - person Youssef Moawad; 11.01.2015
comment
Небольшое примечание: для запроса рабочей области вам не нужен массив, будет достаточно var workspaceQuery : Double = 0.0. - person Martin R; 11.01.2015
comment
Хорошо, у меня так было, потому что я изначально повторно использовал workspace, затем передумал, но не принял это рассуждение полностью... - person Airspeed Velocity; 11.01.2015