Система релаксации LAPACK (прямоугольный полный пакет (RFP))

Question

Я использую формат RFP для хранения матрицы и пытаюсь решить систему уравнений. Но результаты ошибочны. = (Может кто-то пожалуйста, помогите.

я Ь = {5,5, 10, 8,5}

Но должно быть Ь = {2,875, 4,75, 3,5}

Я не понимаю, ., где я сделал ошибку Просто используйте стандартные функции: разложение на множители, а затем решение факторизованных матриц

Спасибо

#include "mkl.h" 
#include "mkl_lapacke.h" 
#include <math.h> 
#include <iostream> 
using namespace std; 
#define NI 3 
#define NJ 1 

int main(int argc, char* argv[]) 
{ 

    double a[NI][NI]; 
    double b[NI][NJ]; 

    a[0][0] = 2; a[0][1] = -1; a[0][2] = 0; 
    a[1][0] = 0; a[1][1] = 2; a[1][2] = -1; 
    a[2][0] = 0; a[2][1] = 0; a[2][2] = 2; 

    b[0][0] = 1; 
    b[0][1] = 6; 
    b[0][2] = 7; 

    cout << "A1 = \n"; 
    for(int i = 0; i < NI; i++) { 
    for(int j = 0; j < NI; j++) { 
     cout << a[i][j] << "\t"; 
    } 
    cout << "\n"; 
    } 
    cout << "\n"; 



    cout << "B1 = \n"; 
    for(int i = 0; i < NI; i++) { 
    for(int j = 0; j < NJ; j++) { 
     cout << b[i][j] << "\t"; 
    } 
    cout << "\n"; 
    } 
    cout << "\n"; 

    char transr = 'N'; 
    char uplo = 'U'; 
    lapack_int n = NI; 
    lapack_int lda = NI; //LDA is used to define the distance in memory between elements of two consecutive columns which have the same row index. 
    double * arf = new double[ n * (n + 1)/2 ]; 
    lapack_int info = -1; 

Conver.. т общая матрица для RFP-формате

info = LAPACKE_dtrttf(LAPACK_ROW_MAJOR, transr, uplo, n, *a, lda, arf); 
    //lapack_int LAPACKE_<?>trttf(int matrix_order, char transr, char uplo, lapack_int n, const <datatype>* a, lapack_int lda, <datatype>* arf); 
    cout << "LAPACKE_dtrttf = " << info << "\n"; 

    cout << "Rectangular full packed: \n"; 
    //cout.setf(std::ios::scientific); 
    for(int i = 0; i < NI; i++) { 
    for(int j = 0; j < (NI+1)/2; j++) { 
     cout << arf[i * (NI+1)/2 + j] << "\t"; 
     //cout << arf[i] << "\t"; 
    } 
    cout << "\n"; 
    } 
    cout << "\n"; 

факторизовать матрица

int matrix_order = LAPACK_ROW_MAJOR; 
    transr = 'N'; 
    uplo = 'U'; 
    n = NI; 

    info = LAPACKE_dpftrf(matrix_order, transr, uplo, n, arf); 
    //lapack_int LAPACKE_<?>pftrf(int matrix_order, char transr, char uplo, lapack_int n, <datatype>* a); 
    cout << "INFO LAPACKE_dpftrf = " << info << "\n"; 
    cout << "Factorized matrix: " << endl; 

    for(int i = 0; i < NI; i++) { 
    for(int j = 0; j < (NI+1)/2; j++) { 
     cout << arf[i*(NI+1)/2+j] << "\t"; 
    } 
    cout << "\n"; 
    } 
    cout << "\n"; 

Решить систему

lapack_int nrhs = NJ; 
    lapack_int ldb = NJ; 
    info = LAPACKE_dpftrs(matrix_order, transr, uplo, n, nrhs, arf, &b[0][0], ldb); 
    //lapack_int LAPACKE_<?>pftrs(int matrix_order, char transr, char uplo, lapack_int n, lapack_int nrhs, const <datatype>* a, <datatype>* b, lapack_int ldb); 
    cout << "INFO LAPACKE_dpftrs = " << info << "\n"; 

Результаты

cout << "Solved = \n"; 
    for(int i = 0; i < NI; i++) { 
    for(int j = 0; j < NJ; j++) { 
     cout << b[i][j] << "\t"; 
    } 
    cout << "\n"; 
    } 
    cout << "\n"; 
    delete [] arf; 

    char ch; 
    cin.get(ch); 

    return 0; 
} 

Answer 1

Вы уверены, что значение NJ является правильным для следующее:

b[NI][NJ]; 
b[0][0] = 1; 
b[0][1] = 6; 
b[0][2] = 7; 

где NI = 3 и NJ = 1 ...

Я думаю, что вы потеряли значения строк с колонкой. Если это была проблема, тогда вам, возможно, захочется узнать, почему она не дала никаких ошибок, это опять вопрос, который уже здесь: Accessing an array out of bounds gives no error, why?