0

我使用 CUSP 共轭梯度法来求解我的对称稀疏矩阵。而且我不知道为什么它不收敛。我使用的矩阵尺寸并不大（1K 到 100K）。相同的线性系统很容易用 MKL 求解，所以矩阵不是病态的。但是我尝试添加预处理器，但没有结果：

对角线预处理器和 AINV（不完全 Cholesky）给出了残差的无限增长（只要 cg 和 bicgstab）

这是我的代码：

cusp::csr_matrix <int, float, cusp::device_memory> A (n, n, nnz);

for (i = 0; i < n + 1; i++)
    A.row_offsets[i] = csrRowPtr[i] - 1;
for (i = 0; i < nnz; i++)
    A.values[i] = csrVal[i];
for (i = 0; i < nnz; i++)
    A.column_indices[i] = csrColInd[i] - 1;

cusp::array1d <float, cusp::device_memory> x (A.num_rows, 0);
cusp::array1d <float, cusp::device_memory> b (A.num_rows, 1);

for (i = 0; i < n; i++)
    b[i] = b_host[i];

cusp::verbose_monitor<float> monitor(b, 100, 1e-3);
cusp::identity_operator<float, MemorySpace> M(A.num_rows, A.num_rows);
    /*
    cusp::precond::diagonal<float, MemorySpace> M(A);
    cusp::precond::scaled_bridson_ainv<float, MemorySpace> M(A, .1);
    */
cusp::krylov::cg(A, x, b, monitor, M);

for (i = 0; i < n; i++)
    x_host[i] = x[i];

为什么它不能正常工作？

PS 据我了解，CUSP 假设索引从零开始，这就是我减少 csrRowPtr 和 csrColInd 的原因。当我使用 nVidia cuSparse 库时，有一个选项可以设置其他参数，例如矩阵类型和填充模式。我如何确定它们在 CUSP 中设置正确？

4

1 回答 1

1

只有上三角形的元素以 MKL 的 CSR 格式存储，但所有元素都必须以 CUSP 的 CSR 格式存储，即使您正在求解对称线性系统。

我也觉得

for (i = 0; i < n; i++)
    x_host[i] = x[i];

不是一个好主意；首先将其传输回 host_memory

cusp::array1d<float, cusp::host_memory> _x = x;

然后将其复制回 x_host 或任何您的结果数组

for (i = 0; i < n; i++)
    x_host[i] = _x[i];
于 2014-04-02T04:28:05.777 回答