Я пытаюсь сделать следующее в CUSP:
A=[
1,1,0,0;
2,2,2,0;
0,3,3,3;
0,0,4,4];
B=[1,1,1,1]';
disp(mldivide(A,B));
который
X=[0.9167,0.0833,-0.5000,0.7500]
С другой стороны, я получаю странный ответ от CUSP
#include <cusp/dia_matrix.h>
#include <cusp/krylov/cg.h>
#include <cusp/print.h>
int main()
{
cusp::dia_matrix<int,float,cusp::host_memory> A(4,4,10,3);
A.diagonal_offsets[0] = -1;
A.diagonal_offsets[1] = 0;
A.diagonal_offsets[2] = 1;
for (int i = 0;i <3;i++)
{
for (int q = 0 ;q < A.num_cols;q++)
{
A.values(q,i)=q+1;
}
}
//copy
cusp::dia_matrix<int,float,cusp::device_memory> AA = A;
cusp::array1d<float,cusp::device_memory> BB(A.num_rows,1);
cusp::array1d<float,cusp::device_memory> XX(A.num_rows,0);
cusp::print(AA);
cusp::print(XX);
cusp::print(BB);
cusp::krylov::cg(AA,XX,BB);\
cusp::print(XX);
return 0;
}
Результат выглядит как
sparse matrix <4, 4> with 10 entries
0 0 1
0 1 1
1 0 2
1 1 2
1 2 2
2 1 3
2 2 3
2 3 3
3 2 4
3 3 4
array1d <4>
0
0
0
0
array1d <4>
1
1
1
1
array1d <4>
-39.9938
-53.436
87.9025
-30.1429
Последний выглядит не совсем правильно. Кто-нибудь знает, что я делаю не так? Я неправильно использую код или у нас должно быть действительно удачное угаданное решение + использовать предварительный кондиционер?
метод сопряженных градиентов допустимо только для использования в симметричных положительно определенных матрицах. Ваша матрица не симметрична. Вот почему он не (и не может) создать правильное решение. Либо используйте подходящую, хорошо кондиционированную SPD-матрицу, либо используйте другой численный метод.