Реализация подпрограммы LAPACK hseqr () в Julia (создание ошибки оболочки: метод не соответствует Array {Float64,1} (:: Int64)) - PullRequest
1 голос
/ 23 апреля 2020

Я пытаюсь реализовать процедуру hseqr () от LAPACK (ссылка: http://www.netlib.org/lapack/explore-html/da/dba/group__double_o_t_h_e_rcomputational_gacb35e85b362ce8ccf9d653cc3f8fb89c.html#gacb35e85b362ce8ccf9d653cc3f8fb89c) до Джулии, но я получаю эту ошибку при вызове функции:

error: no method matching Array{Float64,1}(::Int64)

Это мой код:

for (hseqr, elty) in
    ((:dhseqr_,:Float64),
    (:shseqr_,:Float32))
    @eval begin
        """
     JOB    
          JOB is CHARACTER*1
           = 'E':  compute eigenvalues only;
           = 'S':  compute eigenvalues and the Schur form T.
    COMPZ   
          COMPZ is CHARACTER*1
           = 'N':  no Schur vectors are computed;
           = 'I':  Z is initialized to the unit matrix and the matrix Z
                   of Schur vectors of H is returned;
           = 'V':  Z must contain an orthogonal matrix Q on entry, and
                   the product Q*Z is returned.
        """
        function hseqr!(job::Char,compz::Char,ilo::Integer,ihi::Integer,H::StridedMatrix{$elty}, Z::StridedMatrix{$elty})
            N=size(H,1)
            ldh=size(H,1)
            ldz=size(Z,1)
            work = Vector{$elty}(1)
            lwork = BlasInt(-1)
            info = Ref{BlasInt}()
            for i = 1:2  # first call returns lwork as work[1]
                ccall((@blasfunc($hseqr), liblapack),Void,
                    (Ref{UInt8},Ref{UInt8}, Ref{BlasInt},Ref{BlasInt},Ref{BlasInt},
                        Ptr{$elty}, Ref{BlasInt},Ptr{$elty}, Ptr{$elty}, Ptr{$elty}, 
                        Ref{BlasInt}, Ptr{BlasInt}, Ref{BlasInt}, Ref{BlasInt}),
                        job,compz, N,ilo,ihi,H,ldh,wr,wi,Z,ldz,work, 
                        lwork,info) 
                chklapackerror(info[])
                if i == 1
                    lwork = BlasInt(real(work[1]))
                    resize!(work, lwork)
                end
            end
            return wr,wi,Z
        end
    end
end
hseqr!(job::Char,compz::Char,H::StridedMatrix{},Z::StridedMatrix{}) = hseqr!(job,compz,1, size(H, 1), H,Z)

Это мой вызов: hseqr!('E','N',H,Matrix{Float64}(I,5,5)) (H - матрица Гессенберга с размерностью 5).

Я не уверен, что я правильно понимаю, как сделать обертки, поэтому любой намек поможет.

1 Ответ

2 голосов
/ 23 апреля 2020

Это отличный пример для вопроса, который может быть полезен, если его немного разобрать, чтобы изолировать проблему.

Если я правильно прочитал ваш код, строка

work = Vector{$elty}(1)

вызывается с Float64 в качестве типа элемента, поэтому все сводится к следующему:

julia> Vector{Float64}(1)
ERROR: MethodError: no method matching Array{Float64,1}(::Int64)

Вы неправильно используете конструктор для векторов. Теперь я не уверен, что вы пытаетесь сделать здесь, но если вы хотите создать вектор из Float64 чисел длины один, вы можете искать один из них:

julia> Array{Float64}(undef, 1)
1-element Array{Float64,1}:
 6.9274314537094e-310

julia> zeros(1)
1-element Array{Float64,1}:
 0.0

Обратите внимание, что при первом вызове массив остается неинициализированным, т. Е. Он заполняется всем, что раньше было в памяти, поэтому вы обязательно должны использовать его только в том случае, если впоследствии его перезаписываете.

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...