2015-05-13 7 views
1

Я пытаюсь решить линейную систему, используя следующий код:Использования LAPACKE_zgetrs с LAPACK_ROW_MAJOR вызывает незаконный доступ к памяти

#include <stdio.h> 
#include <lapacke.h> 

int main() { 
    lapack_complex_double mat[4]; 
    lapack_complex_double vec[2]; 
    lapack_int p[2]; 

    mat[0] = lapack_make_complex_double(1,0); 
    mat[1] = lapack_make_complex_double(1,0); 
    mat[2] = lapack_make_complex_double(1,0); 
    mat[3] = lapack_make_complex_double(-1,0); 

    vec[0] = lapack_make_complex_double(1,0); 
    vec[1] = lapack_make_complex_double(1,0); 

    LAPACKE_zgetrf(LAPACK_ROW_MAJOR, 2, 2, mat, 2, p); 
    LAPACKE_zgetrs(LAPACK_ROW_MAJOR, 'N', 2, 1, mat, 2, p, vec, 2); 

    printf("%g %g\n", lapack_complex_double_real(vec[0]), 
     lapack_complex_double_imag(vec[0])); 
    return 0; 
} 

По некоторым причинам, это вызывает незаконный доступ к памяти в LAPACKE_zgetrs (в обнаруженной valgrind и мой большая ошибка программы в zgetrs из-за «glibc обнаружил коррупцию или двойной свободный»). Я не включал это в свой SSCCE для краткости, но все процедуры, которые возвращают, возвращают 0.

Тот же код с работает безотказно.

Мой лапак, лапак и т. Д. Является самомонтирующимся для Ubuntu 12.04. Я использовал следующие настройки в файле LaPack CMake:

BUILD_COMPLEX  ON 
BUILD_COMPLEX16  ON 
BUILD_DOUBLE  ON 
BUILD_SHARED_LIBS ON 
BUILD_SINGLE  ON 
BUILD_STATIC_LIBS ON 
BUILD_TESTING  ON 
CMAKE_BUILD_TYPE Release 
LAPACKE    ON 
LAPACKE_WITH_TMG ON 

, а остальные (оптимизированный Блас/LAPACK и xblas) выкл. Во время сборки ошибок не было, и все тесты были успешными.

Где я испортил?

Редактировать: Я просто попробовал это с Fedora21 и упакованным лапакцем. Он сделал не воспроизвести ошибку.

Edit 2: Хотя он не воспроизводит память терпит неудачу, он производит неправильное решение, а именно (1 + 0I, 1 + 0I) для выше входа (должен быть (1,0))

+0

Я нашел [эту тему] (https://software.intel.com/en-us/forums/topic/499788), которая указывает на ошибку в LAPACKE с большой упаковкой строк при вычислении холесной факторизации. Я не уверен, что 'zgetrf'performs один, и если он вызывает аналогичную ошибку, но это кажется возможным, если вы хотите иметь разложение LU. – martin

+0

@martin холесная факторизация - это разложение LL^T. – ztik

+0

@ctheo yes, который в некоторых случаях является особым случаем LU-разложения. Но вы правы, матрица не является положительно определенной. – martin

ответ

0

После еще несколько исследований и overthinking вещи, я нашел виновник:

Использование LAPACK_ROW_MAJOR переключает значение основных параметров измерения ld*. В то время как ведущим измерением нормального массива Fortran являются номера строк, переключение на ROW_MAJOR переключает его значение на число столбцов. Поэтому правильные звонки (давать правильные результаты) будет:

LAPACKE_zgetrs(LAPACK_ROW_MAJOR, 'N', 2, 1, mat, 2, p, vec, 1); 

где второй 2 это число столбцов (не строки!) Из mat, а последний параметр должен равняться числу правых частей nrhs (а не число переменных!). Я выделил этот самый вызов, потому что все другие вызовы в моем проекте касались квадратных матриц, поэтому «неправильные» вызовы не имеют никакого отрицательного эффекта из-за симметрии.

Как обычно, если вы пропускаете столбцы в конце, ведущие размеры соответственно увеличиваются, так как они будут пропускать строки в обычной настройке.

Очевидно, что это не упоминается в документах Fortran. К сожалению, я не видел такого замечания в документации Лапакца, которая бы спасла мне пару часов моей жизни. :)