Использование 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 для краткости, но все возвращаемые подпрограммы LAPACKE возвращают 0.

Тот же код с LAPACK_COL_MAJOR работает и valgrind работает безупречно.

Мой lapacke, lapack и т. д. собраны самостоятельно для 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

а остальное (оптимизированный blas/lapack и xblas) отключил. При сборке ошибок не было, все тесты прошли успешно.

Где я накосячил?

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

Редактировать 2: хотя он не воспроизводит сбой памяти, он выдает неправильное решение, а именно (1 + 0I, 1 + 0I) для приведенного выше ввода (должно быть (1,0))


c lapack lapacke
person Baum mit Augen    schedule 13.05.2015    source источник
comment
Я нашел эту тему, которая указывает на ошибку в LAPACKE с основной упаковкой строк при вычислении факторизации Холецкого. Я не уверен, что zgetrfвыполняет один, и вызывает ли он аналогичную ошибку, но это кажется возможным, учитывая, что вы хотите иметь декомпозицию LU.   -  person martin    schedule 13.05.2015
comment
@martin факторизация холецкого - это разложение LL ^ T.   -  person ztik    schedule 13.05.2015
comment
@ctheo да, что по-прежнему является частным случаем разложения LU, полезным в некоторых случаях. Но вы правы, матрица не является положительно определенной.   -  person martin    schedule 13.05.2015
comment
В любом случае, в потоке упоминается ошибка в транспонировании в случае симметричных матриц, которая все еще может быть (и я думаю, что это должна быть ошибка либо в LAPACK, либо в LAPACKE, потому что в противном случае не имеет смысла, что это работает в основной режим столбца, но не основной режим строки)   -  person martin    schedule 13.05.2015
comment
@martin На самом деле это не особый случай, поскольку один из результатов LU-разложения будет равен 1 по диагонали, в отличие от холецкого. Но да, наверное баг.   -  person Baum mit Augen    schedule 13.05.2015
comment
Вы пробовали указанную версию с включенным исправлением?   -  person martin    schedule 13.05.2015
comment
@martin Да, использование LAPACK_COL_MAJOR и решение 'T' вместо 'N' похоже работает.   -  person Baum mit Augen    schedule 13.05.2015
comment
Я проверил это с моей сборкой и отлично работает в обоих случаях. Также valgrind ничего не сообщает. Какую версию LAPACK вы используете?   -  person ztik    schedule 13.05.2015
comment
@ктео 3.5. как скачал и построил вчера. Возможно, я сделал что-то не так при сборке библиотеки, поэтому я включил детали сборки в опции.   -  person Baum mit Augen    schedule 13.05.2015
comment
Просто в качестве уточнения: работают как ROW_MAJOR с T, так и COL_MAJOR с N? А другие комбинации не работают?   -  person martin    schedule 13.05.2015
comment
Нет, мне нужен был обходной путь, чтобы добиться своей цели, поэтому я полностью отказался от ROW_MAJOR и притворился, что матрица равна COL_MAJOR, а затем исправил эту ложь, используя вместо нее 'T'. Я считаю, что это математически то же самое, в конце концов. Для моего приложения меня заботило только решение, а не разложение   -  person Baum mit Augen    schedule 13.05.2015

Ответы (1)


arrow_upward
0
arrow_downward

После еще нескольких исследований и размышлений я нашел виновника:

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

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

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

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

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

person Baum mit Augen    schedule 13.05.2015