81 |
dim_t i,k,j; |
dim_t i,k,j; |
82 |
Paso_SparseMatrix * schur=NULL; |
Paso_SparseMatrix * schur=NULL; |
83 |
Paso_SparseMatrix * schur_withFillIn=NULL; |
Paso_SparseMatrix * schur_withFillIn=NULL; |
84 |
double time0=0,time1=0,time2=0,S=0; |
double S=0; |
85 |
/*Paso_Pattern* test;*/ |
/*Paso_Pattern* test;*/ |
86 |
|
|
87 |
/* identify independend set of rows/columns */ |
/* identify independend set of rows/columns */ |
112 |
/* identify independend set of rows/columns */ |
/* identify independend set of rows/columns */ |
113 |
#pragma omp parallel for private(i) schedule(static) |
#pragma omp parallel for private(i) schedule(static) |
114 |
for (i=0;i<n;++i) mis_marker[i]=-1; |
for (i=0;i<n;++i) mis_marker[i]=-1; |
115 |
Paso_Pattern_coup(A_p,mis_marker,couplingParam); |
/*Paso_Pattern_coup(A_p,mis_marker,couplingParam);*/ |
116 |
|
Paso_Pattern_RS(A_p,mis_marker,couplingParam); |
117 |
|
/*Paso_Pattern_Aggregiation(A_p,mis_marker,couplingParam);*/ |
118 |
if (Paso_noError()) { |
if (Paso_noError()) { |
119 |
#pragma omp parallel for private(i) schedule(static) |
#pragma omp parallel for private(i) schedule(static) |
120 |
for (i = 0; i < n; ++i) counter[i]=mis_marker[i]; |
for (i = 0; i < n; ++i) counter[i]=mis_marker[i]; |
154 |
if( Paso_noError()) { |
if( Paso_noError()) { |
155 |
/* if there are no nodes in the coarse level there is no more work to do */ |
/* if there are no nodes in the coarse level there is no more work to do */ |
156 |
out->n_C=n-out->n_F; |
out->n_C=n-out->n_F; |
157 |
if (level>0 && out->n_F>500) { |
if (level>0 && out->n_C>500) { |
158 |
/*if (out->n_F>500) {*/ |
/*if (out->n_F>500) {*/ |
159 |
out->rows_in_C=MEMALLOC(out->n_C,index_t); |
out->rows_in_C=MEMALLOC(out->n_C,index_t); |
160 |
out->mask_C=MEMALLOC(n,index_t); |
out->mask_C=MEMALLOC(n,index_t); |
249 |
if (Paso_noError()) { |
if (Paso_noError()) { |
250 |
if (verbose) { |
if (verbose) { |
251 |
printf("AMG: %d unknowns eliminated. %d left.\n",out->n_F,n-out->n_F); |
printf("AMG: %d unknowns eliminated. %d left.\n",out->n_F,n-out->n_F); |
252 |
if (level>0 && out->n_F>500) { |
/*if (level>0 && out->n_F>500) {*/ |
253 |
/* if (out->n_F<500) {*/ |
/* if (out->n_F<500) {*/ |
254 |
printf("timing: AMG: MIS/reordering/elemination : %e/%e/%e\n",time2,time0,time1); |
/* printf("timing: AMG: MIS/reordering/elemination : %e/%e\n",time0,time1); |
255 |
} else { |
} else { |
256 |
printf("timing: AMG: MIS: %e\n",time2); |
printf("timing: AMG: MIS: %e\n",time2); |
257 |
} |
}*/ |
258 |
} |
} |
259 |
return out; |
return out; |
260 |
} else { |
} else { |
269 |
|
|
270 |
in fact it solves |
in fact it solves |
271 |
|
|
272 |
[ I 0 ] [ A_FF 0 ] [ I invA_FF*A_FF ] [ x_F ] = [b_F] |
[ I 0 ] [ A_FF 0 ] [ I invA_FF*A_FC ] [ x_F ] = [b_F] |
273 |
[ A_CF*invA_FF I ] [ 0 S ] [ 0 I ] [ x_C ] = [b_C] |
[ A_CF*invA_FF I ] [ 0 S ] [ 0 I ] [ x_C ] = [b_C] |
274 |
|
|
275 |
in the form |
in the form |
290 |
void Paso_Solver_solveAMG(Paso_Solver_AMG * amg, double * x, double * b) { |
void Paso_Solver_solveAMG(Paso_Solver_AMG * amg, double * x, double * b) { |
291 |
dim_t i; |
dim_t i; |
292 |
double *r=MEMALLOC(amg->n,double); |
double *r=MEMALLOC(amg->n,double); |
|
/*Paso_Solver_GS* GS=NULL;*/ |
|
293 |
double *x0=MEMALLOC(amg->n,double); |
double *x0=MEMALLOC(amg->n,double); |
294 |
double time0=0; |
/*double time0=0;*/ |
295 |
|
|
296 |
if (amg->level==0 || amg->n_F<=500) { |
if (amg->level==0 || amg->n_C<=500) { |
297 |
/*if (amg->n_F<=500) {*/ |
/*if (amg->n_F<=500) {*/ |
298 |
time0=Paso_timer(); |
/*time0=Paso_timer();*/ |
299 |
|
|
300 |
Paso_Solver_solveJacobi(amg->GS,x,b); |
Paso_Solver_solveJacobi(amg->GS,x,b); |
301 |
|
|
308 |
x[i]+=x0[i]; |
x[i]+=x0[i]; |
309 |
} |
} |
310 |
*/ |
*/ |
311 |
|
|
312 |
/*Paso_UMFPACK1(amg->A,x,b,0);*/ |
/*Paso_UMFPACK1(amg->A,x,b,0);*/ |
313 |
|
|
314 |
time0=Paso_timer()-time0; |
/*time0=Paso_timer()-time0; |
315 |
/*fprintf(stderr,"timing: DIRECT SOLVER: %e/\n",time0);*/ |
fprintf(stderr,"timing: DIRECT SOLVER: %e/\n",time0);*/ |
316 |
} else { |
} else { |
317 |
|
|
318 |
/* presmoothing */ |
/* presmoothing */ |