我想弄清楚为什么我的稀疏(CRS)高斯赛德尔求解器这么慢。我试图在网上找到稀疏格式的 Gauss Seidel 方法的实现,但只能找到使用密集矩阵的实现。这是我的代码:
void gs(int r[], int c[], double v[], double x[], int n, double tol)
{
//x is initially b in Ax=b
double *b = new double[n];
for(int i=0;i<n;i++){b[i] = x[i];}
int ii = 0, jj = 0;
double err = 1.0;
while(err>tol && ii<MAX_ITER){
//Gauss-Seidel iteration
double sigma;
double ajj;
for(int j=0;j<n;j++){
sigma = 0.0;
ajj = 0.0; //diagonal entry a_jj
for(int k=r[j];k<r[j+1];k++){
if(c[k]!=j){
sigma = sigma + v[k]*x[c[k]];
}
else{
ajj = v[k];
}
}
x[j] = (b[j] - sigma)/ajj;
}
if(jj==4){
//err = error(ar,ac,av,x,b,n);
err = fast_error(r,c,v,x,b,n,tol);
jj = 0;
#if(DEBUG)
std::cout<<"error: "<<err<<std::endl;
#endif
}
ii++;
jj++;
}
delete[] b;
}
注意:我使用预处理器标志来选择如何计算错误,但我已经分析过,这不是问题。
在上面的代码片段中,数组r
、c
和v
表示压缩行存储中的稀疏矩阵A
,b 是 的右侧Ax=b
,n 是矩阵的维数A
,tol
是我们要解决的容错度。我的问题是:
- 这就是你将如何实现稀疏 CRS Gauss Seidel 求解器的方法吗?
- 有什么明显的东西可以改变来加快这段代码的速度吗?
为什么我认为我的实施很慢?我使用这个求解器作为我的 AMG 实现的平滑器,最高级别的 Gauss Seidel 平滑花费大部分时间(比计算 Galerkin 三重矩阵乘积要多得多),即使我可能只使用 72 次 Gauss Seidel 迭代最高级别。任何帮助是极大的赞赏。