Почему LAPACKE_dsygvd возвращает ошибку после изменения размера матрицы? - PullRequest
Новая
       52

Почему LAPACKE_dsygvd возвращает ошибку после изменения размера матрицы?

0 голосов
/ 11 марта 2020

Я пытаюсь решить обобщенную проблему собственных значений для атома водорода, используя LAPACKE_dsygvd. Для параметров функций генератора я использую интервал, который начинается с 0,01 и N шагов 0,01. Что я изменяю, так это значение N. Все хорошо для N = 14 и ниже, где я получаю ответы из аналитического решения. Однако, когда я выбираю N = 15 и выше, я получаю ошибку, и информация возвращается со значением> N. После прочтения документации из LAPACK, она говорит следующее:

N: если INFO = N + i, для 1 <= i <= N, то старший минор порядка i в B не является положительно определенным. Факторизация B не может быть завершена, и никакие собственные значения или собственные векторы не были вычислены. </p>

Но я проверил свою матрицу B, и она положительно определена. Я не знаю, что не так.

Ниже я показываю свои сценарии 

#include <cmath>
#include <cstdio>
#include <stdio.h>
#include <stdlib.h>
#include <iostream>
#include "library.h"
#include "mkl.h"

using namespace std;

double Superposition(const double ai, const double aj, const int m);

double Hamiltonian(const double ai, const double aj, const int m);

void print_matrix(double *A, int n) {
    for (int i = 0; i < n; i++) {
        for (int j = 0; j < n; j++) {
            printf("%.7f ", A[i*n + j]);
        }
        cout << "\n";
    }
}

void print_vector(double *vec, int n) {
    for (int i = 0; i < n; i++) {
        cout << vec[i] << " ";
    }
    cout << "\n";
}

double* interval(double min, double step) {
    double *result;
    result = (double *)mkl_malloc( N*sizeof( double ), 64 );
    for (int i = 0; i < N; i++) {
        result[i] = min + i*step;
    }
    return result;
}   
int main() {
    cout << Ry << "\n";
    double *S, *H, *I, *eigenvalues;
    double alpha, beta;
    int i, j, info;

    char* uplo = "U"; char* jobz = "V";

    I = interval(0.01, 0.01);

    alpha = 1.0; beta = 0.0;

    S = (double *)mkl_malloc( N*N*sizeof( double ), 64 );
    H = (double *)mkl_malloc( N*N*sizeof( double ), 64 );
    eigenvalues = (double *)mkl_malloc( N*sizeof( double ), 64 );

    for (i = 0; i < N; i++) {
        for (j = 0; j < N; j++) {
            int index = i*N + j;
            if (j < i) {
                S[index] = 0.0;
                H[index] = 0.0;
            }
            else {
                S[index] = Superposition(I[i], I[j], m);
                H[index] = Hamiltonian(I[i], I[j], m);
            }
        }
    }

    print_matrix(S, N); cout << "\n";
    print_matrix(H, N); cout << "\n" << "\n";

    info = LAPACKE_dsygv(LAPACK_ROW_MAJOR, 1, *jobz, *uplo, N,
        H, N, S, N, eigenvalues);

    //print_matrix(H, N); cout << "\n";

    //for (i = 0; i < N; i++) {
    //  eigenvalues[i] /= Ry;
    //}
    cout << info << "\n" << "\n";
    print_matrix(H, N); cout << "\n";
    print_vector(eigenvalues, N);

    mkl_free(S);
    mkl_free(H);
    mkl_free(I);
    mkl_free(eigenvalues);
}

* Редактирование: я использовал dsygvd, как включенный в MKL, и той же ошибки не происходит. Однако я получаю очень разные результаты для обеих функций, используя одни и те же входные данные.

...