Fast Auxiliary Space Preconditioning  1.8.4 Feb/15/2016
precond_bsr.c
Go to the documentation of this file.
1 
6 #ifdef _OPENMP
7 #include <omp.h>
8 #endif
9 
10 #include "fasp.h"
11 #include "fasp_functs.h"
12 
13 #include "mg_util.inl"
14 
15 /*---------------------------------*/
16 /*-- Public Functions --*/
17 /*---------------------------------*/
18 
38  REAL *z,
39  void *data)
40 {
41  precond_diagbsr *diag = (precond_diagbsr *)data;
42  const INT nb = diag->nb;
43 
44  switch (nb) {
45 
46  case 2:
47  fasp_precond_dbsr_diag_nc2( r, z, diag);
48  break;
49  case 3:
50  fasp_precond_dbsr_diag_nc3( r, z, diag);
51  break;
52 
53  case 5:
54  fasp_precond_dbsr_diag_nc5( r, z, diag);
55  break;
56 
57  case 7:
58  fasp_precond_dbsr_diag_nc7( r, z, diag);
59  break;
60 
61  default:
62  {
63  REAL *diagptr = diag->diag.val;
64  const INT nb2 = nb*nb;
65  const INT m = diag->diag.row/nb2;
66  unsigned INT i;
67 
68 #ifdef _OPENMP
69  if (m > OPENMP_HOLDS) {
70  INT myid, mybegin, myend;
71  INT nthreads = FASP_GET_NUM_THREADS();
72 #pragma omp parallel for private(myid, mybegin, myend, i)
73  for (myid = 0; myid < nthreads; myid++) {
74  FASP_GET_START_END(myid, nthreads, m, &mybegin, &myend);
75  for (i=mybegin; i<myend; ++i) {
76  fasp_blas_smat_mxv(&(diagptr[i*nb2]),&(r[i*nb]),&(z[i*nb]),nb);
77  }
78  }
79  }
80  else {
81 #endif
82  for (i = 0; i < m; ++i) {
83  fasp_blas_smat_mxv(&(diagptr[i*nb2]),&(r[i*nb]),&(z[i*nb]),nb);
84  }
85 #ifdef _OPENMP
86  }
87 #endif
88  break;
89  }
90  }
91 }
92 
112  REAL *z,
113  void *data)
114 {
115  precond_diagbsr *diag = (precond_diagbsr *)data;
116  REAL *diagptr = diag->diag.val;
117 
118  unsigned INT i;
119  const INT m = diag->diag.row/4;
120 
121 #ifdef _OPENMP
122  if (m > OPENMP_HOLDS) {
123  INT myid, mybegin, myend;
124  INT nthreads = FASP_GET_NUM_THREADS();
125 #pragma omp parallel for private(myid, mybegin, myend, i)
126  for (myid = 0; myid < nthreads; myid++) {
127  FASP_GET_START_END(myid, nthreads, m, &mybegin, &myend);
128  for (i = mybegin; i < myend; ++i) {
129  fasp_blas_smat_mxv_nc2(&(diagptr[i*4]),&(r[i*2]),&(z[i*2]));
130  }
131  }
132  }
133  else {
134 #endif
135  for (i = 0; i < m; ++i) {
136  fasp_blas_smat_mxv_nc2(&(diagptr[i*4]),&(r[i*2]),&(z[i*2]));
137  }
138 #ifdef _OPENMP
139  }
140 #endif
141 }
142 
162  REAL *z,
163  void *data)
164 {
165  precond_diagbsr *diag = (precond_diagbsr *)data;
166  REAL *diagptr = diag->diag.val;
167 
168  const INT m = diag->diag.row/9;
169  unsigned INT i;
170 
171 #ifdef _OPENMP
172  if (m > OPENMP_HOLDS) {
173  INT myid, mybegin, myend;
174  INT nthreads = FASP_GET_NUM_THREADS();
175 #pragma omp parallel for private(myid, mybegin, myend, i)
176  for (myid = 0; myid < nthreads; myid++) {
177  FASP_GET_START_END(myid, nthreads, m, &mybegin, &myend);
178  for (i = mybegin; i < myend; ++i) {
179  fasp_blas_smat_mxv_nc3(&(diagptr[i*9]),&(r[i*3]),&(z[i*3]));
180  }
181  }
182  }
183  else {
184 #endif
185  for (i = 0; i < m; ++i) {
186  fasp_blas_smat_mxv_nc3(&(diagptr[i*9]),&(r[i*3]),&(z[i*3]));
187  }
188 #ifdef _OPENMP
189  }
190 #endif
191 }
192 
212  REAL *z,
213  void *data)
214 {
215  precond_diagbsr *diag = (precond_diagbsr *)data;
216  REAL *diagptr = diag->diag.val;
217 
218  unsigned INT i;
219  const INT m = diag->diag.row/25;
220 
221 #ifdef _OPENMP
222  if (m > OPENMP_HOLDS) {
223  INT myid, mybegin, myend;
224  INT nthreads = FASP_GET_NUM_THREADS();
225 #pragma omp parallel for private(myid, mybegin, myend, i)
226  for (myid = 0; myid < nthreads; myid++) {
227  FASP_GET_START_END(myid, nthreads, m, &mybegin, &myend);
228  for (i = mybegin; i < myend; ++i) {
229  fasp_blas_smat_mxv_nc5(&(diagptr[i*25]),&(r[i*5]),&(z[i*5]));
230  }
231  }
232  }
233  else {
234 #endif
235  for (i = 0; i < m; ++i) {
236  fasp_blas_smat_mxv_nc5(&(diagptr[i*25]),&(r[i*5]),&(z[i*5]));
237  }
238 #ifdef _OPENMP
239  }
240 #endif
241 }
242 
261  REAL *z,
262  void *data)
263 {
264  precond_diagbsr *diag = (precond_diagbsr *)data;
265  REAL *diagptr = diag->diag.val;
266 
267  unsigned INT i;
268  const INT m = diag->diag.row/49;
269 
270 #ifdef _OPENMP
271  if (m > OPENMP_HOLDS) {
272  INT myid, mybegin, myend;
273  INT nthreads = FASP_GET_NUM_THREADS();
274 #pragma omp parallel for private(myid, mybegin, myend, i)
275  for (myid = 0; myid < nthreads; myid++) {
276  FASP_GET_START_END(myid, nthreads, m, &mybegin, &myend);
277  for (i = mybegin; i < myend; ++i) {
278  fasp_blas_smat_mxv_nc7(&(diagptr[i*49]),&(r[i*7]),&(z[i*7]));
279  }
280  }
281  }
282  else {
283 #endif
284  for (i = 0; i < m; ++i) {
285  fasp_blas_smat_mxv_nc7(&(diagptr[i*49]),&(r[i*7]),&(z[i*7]));
286  }
287 #ifdef _OPENMP
288  }
289 #endif
290 }
291 
307  REAL *z,
308  void *data)
309 {
310  const ILU_data *iludata=(ILU_data *)data;
311  const INT m=iludata->row, mm1=m-1, mm2=m-2, memneed=2*m;
312  const INT nb=iludata->nb, nb2=nb*nb, size=m*nb;
313 
314  INT *ijlu=iludata->ijlu;
315  REAL *lu=iludata->luval;
316 
317  INT ib, ibstart,ibstart1;
318  INT i, j, jj, begin_row, end_row;
319  REAL *zz, *zr, *mult;
320 
321  if (iludata->nwork<memneed) {
322  printf("### ERROR: Need %d memory, only %d available!\n", memneed, iludata->nwork);
323  exit(ERROR_ALLOC_MEM);
324  }
325 
326  zz = iludata->work;
327  zr = zz + size;
328  mult = zr + size;
329 
330  memcpy(zr, r, size*sizeof(REAL));
331 
332  switch (nb) {
333 
334  case 1:
335 
336  // forward sweep: solve unit lower matrix equation L*zz=zr
337  zz[0]=zr[0];
338  for (i=1;i<=mm1;++i) {
339  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
340  for (j=begin_row;j<=end_row;++j) {
341  jj=ijlu[j];
342  if (jj<i) zr[i]-=lu[j]*zz[jj];
343  else break;
344  }
345  zz[i]=zr[i];
346  }
347 
348  // backward sweep: solve upper matrix equation U*z=zz
349  z[mm1]=zz[mm1]*lu[mm1];
350  for (i=mm2;i>=0;i--) {
351  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
352  for (j=end_row;j>=begin_row;j--) {
353  jj=ijlu[j];
354  if (jj>i) zz[i]-=lu[j]*z[jj];
355  else break;
356  }
357  z[i]=zz[i]*lu[i];
358  }
359 
360  break; //end (if nb==1)
361 
362  case 3:
363 
364  // forward sweep: solve unit lower matrix equation L*zz=zr
365  zz[0] = zr[0];
366  zz[1] = zr[1];
367  zz[2] = zr[2];
368 
369  for (i=1;i<=mm1;++i) {
370  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
371  ibstart=i*nb;
372  for (j=begin_row;j<=end_row;++j) {
373  jj=ijlu[j];
374  if (jj<i)
375  {
376  fasp_blas_smat_mxv_nc3(&(lu[j*nb2]),&(zz[jj*nb]),mult);
377  for (ib=0;ib<nb;++ib) zr[ibstart+ib]-=mult[ib];
378  }
379  else break;
380  }
381 
382  zz[ibstart] = zr[ibstart];
383  zz[ibstart+1] = zr[ibstart+1];
384  zz[ibstart+2] = zr[ibstart+2];
385  }
386 
387  // backward sweep: solve upper matrix equation U*z=zz
388  ibstart=mm1*nb2;
389  ibstart1=mm1*nb;
390  fasp_blas_smat_mxv_nc3(&(lu[ibstart]),&(zz[ibstart1]),&(z[ibstart1]));
391 
392  for (i=mm2;i>=0;i--) {
393  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
394  ibstart=i*nb2;
395  ibstart1=i*nb;
396  for (j=end_row;j>=begin_row;j--) {
397  jj=ijlu[j];
398  if (jj>i) {
399  fasp_blas_smat_mxv_nc3(&(lu[j*nb2]),&(z[jj*nb]),mult);
400  for (ib=0;ib<nb;++ib) zz[ibstart1+ib]-=mult[ib];
401  }
402 
403  else break;
404  }
405 
406  fasp_blas_smat_mxv_nc3(&(lu[ibstart]),&(zz[ibstart1]),&(z[ibstart1]));
407 
408  }
409 
410  break; // end (if nb=3)
411 
412  case 5:
413 
414  // forward sweep: solve unit lower matrix equation L*zz=zr
415  fasp_array_cp(nb,&(zr[0]),&(zz[0]));
416 
417  for (i=1;i<=mm1;++i) {
418  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
419  ibstart=i*nb;
420  for (j=begin_row;j<=end_row;++j) {
421  jj=ijlu[j];
422  if (jj<i) {
423  fasp_blas_smat_mxv_nc5(&(lu[j*nb2]),&(zz[jj*nb]),mult);
424  for (ib=0;ib<nb;++ib) zr[ibstart+ib]-=mult[ib];
425  }
426  else break;
427  }
428 
429  fasp_array_cp(nb,&(zr[ibstart]),&(zz[ibstart]));
430  }
431 
432  // backward sweep: solve upper matrix equation U*z=zz
433  ibstart=mm1*nb2;
434  ibstart1=mm1*nb;
435  fasp_blas_smat_mxv_nc5(&(lu[ibstart]),&(zz[ibstart1]),&(z[ibstart1]));
436 
437  for (i=mm2;i>=0;i--) {
438  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
439  ibstart=i*nb2;
440  ibstart1=i*nb;
441  for (j=end_row;j>=begin_row;j--) {
442  jj=ijlu[j];
443  if (jj>i) {
444  fasp_blas_smat_mxv_nc5(&(lu[j*nb2]),&(z[jj*nb]),mult);
445  for (ib=0;ib<nb;++ib) zz[ibstart1+ib]-=mult[ib];
446  }
447 
448  else break;
449  }
450 
451  fasp_blas_smat_mxv_nc5(&(lu[ibstart]),&(zz[ibstart1]),&(z[ibstart1]));
452 
453  }
454 
455  break; //end (if nb==5)
456 
457  case 7:
458 
459  // forward sweep: solve unit lower matrix equation L*zz=zr
460  fasp_array_cp(nb,&(zr[0]),&(zz[0]));
461 
462  for (i=1;i<=mm1;++i) {
463  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
464  ibstart=i*nb;
465  for (j=begin_row;j<=end_row;++j) {
466  jj=ijlu[j];
467  if (jj<i) {
468  fasp_blas_smat_mxv_nc7(&(lu[j*nb2]),&(zz[jj*nb]),mult);
469  for (ib=0;ib<nb;++ib) zr[ibstart+ib]-=mult[ib];
470  }
471  else break;
472  }
473 
474  fasp_array_cp(nb,&(zr[ibstart]),&(zz[ibstart]));
475  }
476 
477  // backward sweep: solve upper matrix equation U*z=zz
478  ibstart=mm1*nb2;
479  ibstart1=mm1*nb;
480  fasp_blas_smat_mxv_nc7(&(lu[ibstart]),&(zz[ibstart1]),&(z[ibstart1]));
481 
482  for (i=mm2;i>=0;i--) {
483  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
484  ibstart=i*nb2;
485  ibstart1=i*nb;
486  for (j=end_row;j>=begin_row;j--) {
487  jj=ijlu[j];
488  if (jj>i) {
489  fasp_blas_smat_mxv_nc7(&(lu[j*nb2]),&(z[jj*nb]),mult);
490  for (ib=0;ib<nb;++ib) zz[ibstart1+ib]-=mult[ib];
491  }
492 
493  else break;
494  }
495 
496  fasp_blas_smat_mxv_nc7(&(lu[ibstart]),&(zz[ibstart1]),&(z[ibstart1]));
497 
498  }
499 
500  break; //end (if nb==7)
501 
502  default:
503 
504  // forward sweep: solve unit lower matrix equation L*zz=zr
505  fasp_array_cp(nb,&(zr[0]),&(zz[0]));
506 
507  for (i=1;i<=mm1;++i) {
508  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
509  ibstart=i*nb;
510  for (j=begin_row;j<=end_row;++j) {
511  jj=ijlu[j];
512  if (jj<i) {
513  fasp_blas_smat_mxv(&(lu[j*nb2]),&(zz[jj*nb]),mult,nb);
514  for (ib=0;ib<nb;++ib) zr[ibstart+ib]-=mult[ib];
515  }
516  else break;
517  }
518 
519  fasp_array_cp(nb,&(zr[ibstart]),&(zz[ibstart]));
520  }
521 
522  // backward sweep: solve upper matrix equation U*z=zz
523  ibstart=mm1*nb2;
524  ibstart1=mm1*nb;
525  fasp_blas_smat_mxv(&(lu[ibstart]),&(zz[ibstart1]),&(z[ibstart1]),nb);
526 
527  for (i=mm2;i>=0;i--) {
528  begin_row=ijlu[i]; end_row=ijlu[i+1]-1;
529  ibstart=i*nb2;
530  ibstart1=i*nb;
531  for (j=end_row;j>=begin_row;j--) {
532  jj=ijlu[j];
533  if (jj>i) {
534  fasp_blas_smat_mxv(&(lu[j*nb2]),&(z[jj*nb]),mult,nb);
535  for (ib=0;ib<nb;++ib) zz[ibstart1+ib]-=mult[ib];
536  }
537 
538  else break;
539  }
540 
541  fasp_blas_smat_mxv(&(lu[ibstart]),&(zz[ibstart1]),&(z[ibstart1]),nb);
542 
543  }
544 
545  break; // end everything else
546  }
547 
548  return;
549 }
550 
564  REAL *z,
565  void *data)
566 {
567  precond_data_bsr *predata=(precond_data_bsr *)data;
568  const INT row=predata->mgl_data[0].A.ROW;
569  const INT nb = predata->mgl_data[0].A.nb;
570  const INT maxit=predata->maxit;
571  unsigned INT i;
572  const INT m = row*nb;
573 
574  AMG_param amgparam; fasp_param_amg_init(&amgparam);
575  amgparam.cycle_type = predata->cycle_type;
576  amgparam.smoother = predata->smoother;
577  amgparam.smooth_order = predata->smooth_order;
578  amgparam.presmooth_iter = predata->presmooth_iter;
579  amgparam.postsmooth_iter = predata->postsmooth_iter;
580  amgparam.relaxation = predata->relaxation;
581  amgparam.coarse_scaling = predata->coarse_scaling;
582  amgparam.tentative_smooth = predata->tentative_smooth;
583  amgparam.ILU_levels = predata->mgl_data->ILU_levels;
584 
585  AMG_data_bsr *mgl = predata->mgl_data;
586  mgl->b.row=m; fasp_array_cp(m,r,mgl->b.val); // residual is an input
587  mgl->x.row=m; fasp_dvec_set(m,&mgl->x,0.0);
588 
589  for (i=0;i<maxit;++i) fasp_solver_mgcycle_bsr(mgl,&amgparam);
590  //fasp_solver_mgrecurmgl,&amgparam,0); //
591 
592  fasp_array_cp(m,mgl->x.val,z);
593 }
594 
608  REAL *z,
609  void *data)
610 {
611  precond_data_bsr *pcdata=(precond_data_bsr *)data;
612  const INT row=pcdata->mgl_data[0].A.ROW;
613  const INT nb=pcdata->mgl_data[0].A.nb;
614  const INT maxit=pcdata->maxit;
615  const SHORT num_levels=pcdata->max_levels;
616  const INT m=row*nb;
617  unsigned INT i;
618 
619  AMG_param amgparam; fasp_param_amg_init(&amgparam);
620  fasp_param_prec_to_amg_bsr(&amgparam,pcdata);
621 
622  AMG_data_bsr *mgl = pcdata->mgl_data;
623  mgl->b.row=m; fasp_array_cp(m,r,mgl->b.val); // residual is an input
624  mgl->x.row=m; fasp_dvec_set(m,&mgl->x,0.0);
625 
626  for (i=0;i<maxit;++i) fasp_solver_nl_amli_bsr(mgl,&amgparam,0, num_levels);
627 
628  fasp_array_cp(m,mgl->x.val,z);
629 }
630 
644  REAL *z,
645  void *data)
646 {
647  precond_data_bsr *predata=(precond_data_bsr *)data;
648  const INT row=predata->mgl_data[0].A.ROW;
649  const INT nb = predata->mgl_data[0].A.nb;
650  const INT maxit=predata->maxit;
651  unsigned INT i;
652  const INT m = row*nb;
653 
654  dCSRmat *A_nk = predata->A_nk;
655  dCSRmat *P_nk = predata->P_nk;
656  dCSRmat *R_nk = predata->R_nk;
657 
658  fasp_array_set(m, z, 0.0);
659 
660  // local variables
661  dvector r_nk, z_nk;
662  fasp_dvec_alloc(A_nk->row, &r_nk);
663  fasp_dvec_alloc(A_nk->row, &z_nk);
664 
665  //----------------------
666  // extra kernel solve
667  //----------------------
668  // r_nk = R_nk*r
669  fasp_blas_dcsr_mxv(R_nk, r, r_nk.val);
670 
671  // z_nk = A_nk^{-1}*r_nk
672 #if WITH_UMFPACK // use UMFPACK directly
673  fasp_solver_umfpack(A_nk, &r_nk, &z_nk, 0);
674 #else
675  fasp_coarse_itsolver(A_nk, &r_nk, &z_nk, 1e-12, 0);
676 #endif
677 
678  // z = z + P_nk*z_nk;
679  fasp_blas_dcsr_aAxpy(1.0, P_nk, z_nk.val, z);
680 
681  //----------------------
682  // AMG solve
683  //----------------------
684  AMG_param amgparam; fasp_param_amg_init(&amgparam);
685  amgparam.cycle_type = predata->cycle_type;
686  amgparam.smoother = predata->smoother;
687  amgparam.smooth_order = predata->smooth_order;
688  amgparam.presmooth_iter = predata->presmooth_iter;
689  amgparam.postsmooth_iter = predata->postsmooth_iter;
690  amgparam.relaxation = predata->relaxation;
691  amgparam.coarse_scaling = predata->coarse_scaling;
692  amgparam.tentative_smooth = predata->tentative_smooth;
693  amgparam.ILU_levels = predata->mgl_data->ILU_levels;
694 
695  AMG_data_bsr *mgl = predata->mgl_data;
696  mgl->b.row=m; fasp_array_cp(m,r,mgl->b.val); // residual is an input
697  mgl->x.row=m; //fasp_dvec_set(m,&mgl->x,0.0);
698  fasp_array_cp(m, z, mgl->x.val);
699 
700  for (i=0;i<maxit;++i) fasp_solver_mgcycle_bsr(mgl,&amgparam);
701  //fasp_solver_mgrecurmgl,&amgparam,0); //
702 
703  fasp_array_cp(m,mgl->x.val,z);
704 
705  //----------------------
706  // extra kernel solve
707  //----------------------
708  // r = r - A*z
709  fasp_blas_dbsr_aAxpy(-1.0, &(predata->mgl_data[0].A), z, mgl->b.val);
710 
711  // r_nk = R_nk*r
712  fasp_blas_dcsr_mxv(R_nk, mgl->b.val, r_nk.val);
713 
714  // z_nk = A_nk^{-1}*r_nk
715 #if WITH_UMFPACK // use UMFPACK directly
716  fasp_solver_umfpack(A_nk, &r_nk, &z_nk, 0);
717 #else
718  fasp_coarse_itsolver(A_nk, &r_nk, &z_nk, 1e-12, 0);
719 #endif
720 
721  // z = z + P_nk*z_nk;
722  fasp_blas_dcsr_aAxpy(1.0, P_nk, z_nk.val, z);
723 }
724 
725 /*---------------------------------*/
726 /*-- End of File --*/
727 /*---------------------------------*/
void fasp_precond_dbsr_amg_nk(REAL *r, REAL *z, void *data)
AMG with extra near kernel solve preconditioner.
Definition: precond_bsr.c:643
dCSRmat * A_nk
Matrix data for near kernal.
Definition: fasp_block.h:382
INT * ijlu
integer array of row pointers and column indexes, the size is nzlu
Definition: fasp.h:412
dvector x
pointer to the iterative solution at level level_num
Definition: fasp_block.h:219
Parameters for AMG solver.
Definition: fasp.h:583
SHORT postsmooth_iter
number of postsmoothers
Definition: fasp.h:619
REAL tentative_smooth
relaxation parameter for smoothing the tentative prolongation
Definition: fasp.h:676
SHORT smoother
smoother type
Definition: fasp.h:610
SHORT smooth_order
AMG smoother ordering.
Definition: fasp_block.h:335
SHORT smooth_order
smoother order
Definition: fasp.h:613
#define REAL
Definition: fasp.h:67
void fasp_precond_dbsr_diag_nc2(REAL *r, REAL *z, void *data)
Diagonal preconditioner z=inv(D)*r.
Definition: precond_bsr.c:111
SHORT presmooth_iter
number of presmoothers
Definition: fasp.h:616
dBSRmat A
pointer to the matrix at level level_num
Definition: fasp_block.h:207
void fasp_precond_dbsr_diag_nc7(REAL *r, REAL *z, void *data)
Diagonal preconditioner z=inv(D)*r.
Definition: precond_bsr.c:260
void fasp_blas_smat_mxv_nc3(REAL *a, REAL *b, REAL *c)
Compute the product of a 3*3 matrix a and a array b, stored in c.
Definition: blas_smat.c:105
INT nb
dimension of each sub-block
Definition: fasp_block.h:56
INT fasp_solver_umfpack(dCSRmat *ptrA, dvector *b, dvector *u, const SHORT prtlvl)
Solve Au=b by UMFpack.
REAL * val
actual vector entries
Definition: fasp.h:348
dCSRmat * P_nk
Prolongation for near kernal.
Definition: fasp_block.h:385
Vector with n entries of REAL type.
Definition: fasp.h:342
SHORT cycle_type
AMG cycle type.
Definition: fasp_block.h:329
INT max_levels
max number of AMG levels
Definition: fasp_block.h:323
#define INT
Definition: fasp.h:64
Data for ILU setup.
Definition: fasp.h:400
void fasp_solver_nl_amli_bsr(AMG_data_bsr *mgl, AMG_param *param, INT level, INT num_levels)
Solve Ax=b with recursive nonlinear AMLI-cycle.
Definition: amlirecur.c:508
void fasp_solver_mgcycle_bsr(AMG_data_bsr *mgl, AMG_param *param)
Solve Ax=b with non-recursive multigrid cycle.
Definition: mgcycle.c:265
INT nwork
work space size
Definition: fasp.h:421
void FASP_GET_START_END(INT procid, INT nprocs, INT n, INT *start, INT *end)
Assign Load to each thread.
Definition: threads.c:83
Data for multigrid levels. (BSR format)
Definition: fasp_block.h:198
void fasp_array_set(const INT n, REAL *x, const REAL val)
Set initial value for an array to be x=val.
Definition: array.c:48
void fasp_precond_dbsr_nl_amli(REAL *r, REAL *z, void *data)
Nonlinear AMLI-cycle AMG preconditioner.
Definition: precond_bsr.c:607
Data passed to the preconditioners.
Definition: fasp_block.h:311
void fasp_precond_dbsr_diag_nc3(REAL *r, REAL *z, void *data)
Diagonal preconditioner z=inv(D)*r.
Definition: precond_bsr.c:161
void fasp_blas_smat_mxv_nc7(REAL *a, REAL *b, REAL *c)
Compute the product of a 7*7 matrix a and a array b, stored in c.
Definition: blas_smat.c:154
#define OPENMP_HOLDS
Definition: fasp_const.h:248
dvector b
pointer to the right-hand side at level level_num
Definition: fasp_block.h:216
REAL * luval
nonzero entries of LU
Definition: fasp.h:415
SHORT ILU_levels
number of levels use ILU smoother
Definition: fasp.h:682
Main header file for FASP.
SHORT coarse_scaling
switch of scaling of the coarse grid correction
Definition: fasp_block.h:353
#define ERROR_ALLOC_MEM
Definition: fasp_const.h:37
INT row
row number of matrix LU, m
Definition: fasp.h:403
SHORT smoother
AMG smoother type.
Definition: fasp_block.h:332
Data passed to diagnal preconditioner for dBSRmat matrices.
Definition: fasp_block.h:293
REAL relaxation
relaxation parameter for SOR smoother
Definition: fasp.h:622
INT row
row number of matrix A, m
Definition: fasp.h:151
void fasp_precond_dbsr_ilu(REAL *r, REAL *z, void *data)
ILU preconditioner.
Definition: precond_bsr.c:306
Sparse matrix of REAL type in CSR format.
Definition: fasp.h:148
AMG_data_bsr * mgl_data
AMG preconditioner data.
Definition: fasp_block.h:368
INT row
number of rows
Definition: fasp.h:345
void fasp_param_amg_init(AMG_param *amgparam)
Initialize AMG parameters.
Definition: parameters.c:390
SHORT cycle_type
type of AMG cycle
Definition: fasp.h:604
void fasp_precond_dbsr_diag_nc5(REAL *r, REAL *z, void *data)
Diagonal preconditioner z=inv(D)*r.
Definition: precond_bsr.c:211
#define SHORT
FASP integer and floating point numbers.
Definition: fasp.h:63
INT ILU_levels
number of levels use ILU smoother
Definition: fasp_block.h:255
dvector diag
diagnal elements
Definition: fasp_block.h:299
void fasp_precond_dbsr_amg(REAL *r, REAL *z, void *data)
AMG preconditioner.
Definition: precond_bsr.c:563
void fasp_param_prec_to_amg_bsr(AMG_param *amgparam, precond_data_bsr *pcdata)
Set AMG_param with precond_data.
Definition: parameters.c:767
void fasp_blas_dcsr_aAxpy(const REAL alpha, dCSRmat *A, REAL *x, REAL *y)
Matrix-vector multiplication y = alpha*A*x + y.
Definition: blas_csr.c:479
REAL tentative_smooth
smooth factor for smoothing the tentative prolongation
Definition: fasp_block.h:362
void fasp_dvec_set(INT n, dvector *x, REAL val)
Initialize dvector x[i]=val for i=0:n-1.
Definition: vec.c:235
void fasp_dvec_alloc(const INT m, dvector *u)
Create dvector data space of REAL type.
Definition: vec.c:99
INT maxit
max number of iterations of AMG preconditioner
Definition: fasp_block.h:320
void fasp_blas_smat_mxv_nc2(REAL *a, REAL *b, REAL *c)
Compute the product of a 2*2 matrix a and a array b, stored in c.
Definition: blas_smat.c:83
INT ROW
number of rows of sub-blocks in matrix A, M
Definition: fasp_block.h:47
INT nb
block size for BSR type only
Definition: fasp.h:418
void fasp_array_cp(const INT n, REAL *x, REAL *y)
Copy an array to the other y=x.
Definition: array.c:165
SHORT presmooth_iter
number of presmoothing
Definition: fasp_block.h:338
REAL * work
work space
Definition: fasp.h:424
void fasp_precond_dbsr_diag(REAL *r, REAL *z, void *data)
Diagonal preconditioner z=inv(D)*r.
Definition: precond_bsr.c:37
SHORT postsmooth_iter
number of postsmoothing
Definition: fasp_block.h:341
INT nb
dimension of each sub-block
Definition: fasp_block.h:296
dCSRmat * R_nk
Resriction for near kernal.
Definition: fasp_block.h:388
void fasp_blas_smat_mxv_nc5(REAL *a, REAL *b, REAL *c)
Compute the product of a 5*5 matrix a and a array b, stored in c.
Definition: blas_smat.c:128
REAL relaxation
relaxation parameter for SOR smoother
Definition: fasp_block.h:347
void fasp_blas_smat_mxv(REAL *a, REAL *b, REAL *c, const INT n)
Compute the product of a small full matrix a and a array b, stored in c.
Definition: blas_smat.c:183
void fasp_blas_dbsr_aAxpy(const REAL alpha, dBSRmat *A, REAL *x, REAL *y)
Compute y := alpha*A*x + y.
Definition: blas_bsr.c:337
SHORT coarse_scaling
switch of scaling of the coarse grid correction
Definition: fasp.h:631
void fasp_blas_dcsr_mxv(dCSRmat *A, REAL *x, REAL *y)
Matrix-vector multiplication y = A*x.
Definition: blas_csr.c:225