小编Alp*_*jet的帖子

我用LU预处理编写了一个共轭梯度求解器(用于线性方程组),我使用了Maxim Naumov博士关于nvidia研究界的论文作为指导,残差更新步骤,需要求解下三角矩阵系统然后求解上三角矩阵系统分为两个阶段:

1. 分析阶段(利用稀疏模式并决定并行化级别).
2. 解决方案阶段本身.

根据与该主题相关的所有帖子,以及Naumov的论文本身,分析阶段明显慢于解决阶段,但它执行了一次,所以考虑到整个执行时间不应该是一个问题,但是,在我的程序中,分析阶段需要整个解决时间的约35-45%(执行所有迭代所需的时间!),这非常烦人,另一件事是我很确定矩阵的稀疏模式不会允许大量并行化,因为几乎所有元素都需要知道先前的元素(因为在CFD中每个节点至少需要相邻的6个节点(六面体体积),并且每个节点都重复),因此分析阶段反正不会很有用!

此代码中的 matrixLU包含上三角矩阵和下三角矩阵,上三角矩阵使用原始矩阵对角线,下三角矩阵具有统一对角线(LU分解).

// z = inv(matrixLU)*r
cusparseMatDescr_t descrL = 0 ;
cusparseMatDescr_t descrU = 0 ;
cusparseStatus = cusparseCreateMatDescr(&descrL) ;
cusparseStatus = cusparseCreateMatDescr(&descrU) ;

cusparseSetMatType(descrL,CUSPARSE_MATRIX_TYPE_GENERAL) ;
cusparseSetMatIndexBase(descrL,CUSPARSE_INDEX_BASE_ONE) ;
cusparseSetMatDiagType(descrL,CUSPARSE_DIAG_TYPE_UNIT) ;
cusparseSetMatFillMode(descrL,CUSPARSE_FILL_MODE_LOWER) ;

cusparseSetMatType(descrU,CUSPARSE_MATRIX_TYPE_GENERAL) ;
cusparseSetMatIndexBase(descrU,CUSPARSE_INDEX_BASE_ONE) ;
cusparseSetMatDiagType(descrU,CUSPARSE_DIAG_TYPE_NON_UNIT) ;
cusparseSetMatFillMode(descrU,CUSPARSE_FILL_MODE_UPPER) ;

cusparseSolveAnalysisInfo_t inforL = 0 ;
cusparseSolveAnalysisInfo_t inforU = 0 ;
cusparseStatus = cusparseCreateSolveAnalysisInfo(&inforL) ;
cusparseStatus = cusparseCreateSolveAnalysisInfo(&inforU) ;

double startFA = omp_get_wtime() ;
cusparseStatus = cusparseDcsrsv_analysis(cusparseHandle, CUSPARSE_OPERATION_NON_TRANSPOSE, N, NZ, descrL, matrixLU, iRow, …