Fast Auxiliary Space Preconditioning 2.7.7 Aug/28/2022
SolCSR.c
Go to the documentation of this file.
1
18#include <time.h>
19
20#ifdef _OPENMP
21#include <omp.h>
22#endif
23
24#include "fasp.h"
25#include "fasp_functs.h"
26
27/*---------------------------------*/
28/*-- Declare Private Functions --*/
29/*---------------------------------*/
30
31#include "KryUtil.inl"
32
33/*---------------------------------*/
34/*-- Public Functions --*/
35/*---------------------------------*/
36
57 dvector *b,
58 dvector *x,
59 precond *pc,
60 ITS_param *itparam)
61{
62 const SHORT prtlvl = itparam->print_level;
63 const SHORT itsolver_type = itparam->itsolver_type;
64 const SHORT stop_type = itparam->stop_type;
65 const SHORT restart = itparam->restart;
66 const INT MaxIt = itparam->maxit;
67 const REAL tol = itparam->tol;
68
69 /* Local Variables */
70 REAL solve_start, solve_end;
71 INT iter;
72
73#if DEBUG_MODE > 0
74 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
75 printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
76#endif
77
78 fasp_gettime(&solve_start);
79
80 /* check matrix data */
82
83 /* Safe-guard checks on parameters */
84 ITS_CHECK ( MaxIt, tol );
85
86 /* Choose a desirable Krylov iterative solver */
87 switch ( itsolver_type ) {
88 case SOLVER_CG:
89 iter = fasp_solver_dcsr_pcg(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
90 break;
91
92 case SOLVER_BiCGstab:
93 iter = fasp_solver_dcsr_pbcgs(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
94 break;
95
96 case SOLVER_MinRes:
97 iter = fasp_solver_dcsr_pminres(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
98 break;
99
100 case SOLVER_GMRES:
101 iter = fasp_solver_dcsr_pgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
102 break;
103
104 case SOLVER_VGMRES:
105 iter = fasp_solver_dcsr_pvgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
106 break;
107
108 case SOLVER_VFGMRES:
109 iter = fasp_solver_dcsr_pvfgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
110 break;
111
112 case SOLVER_GCG:
113 iter = fasp_solver_dcsr_pgcg(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
114 break;
115
116 case SOLVER_GCR:
117 iter = fasp_solver_dcsr_pgcr(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
118 break;
119
120 default:
121 printf("### ERROR: Unknown iterative solver type %d! [%s]\n",
122 itsolver_type, __FUNCTION__);
123 return ERROR_SOLVER_TYPE;
124
125 }
126
127 if ( (prtlvl >= PRINT_SOME) && (iter >= 0) ) {
128 fasp_gettime(&solve_end);
129 fasp_cputime("Iterative method", solve_end - solve_start);
130 }
131
132#if DEBUG_MODE > 0
133 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
134#endif
135
136 return iter;
137}
138
159 dvector *b,
160 dvector *x,
161 precond *pc,
162 ITS_param *itparam)
163{
164 const SHORT prtlvl = itparam->print_level;
165 const SHORT itsolver_type = itparam->itsolver_type;
166 const SHORT stop_type = itparam->stop_type;
167 const SHORT restart = itparam->restart;
168 const INT MaxIt = itparam->maxit;
169 const REAL tol = itparam->tol;
170
171 /* Local Variables */
172 REAL solve_start, solve_end;
173 INT iter;
174
175#if DEBUG_MODE > 0
176 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
177 printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
178#endif
179
180 fasp_gettime(&solve_start);
181
182 /* check matrix data */
184
185 /* Safe-guard checks on parameters */
186 ITS_CHECK ( MaxIt, tol );
187
188 /* Choose a desirable Krylov iterative solver */
189 switch ( itsolver_type ) {
190 case SOLVER_CG:
191 iter = fasp_solver_dcsr_spcg(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
192 break;
193
194 case SOLVER_BiCGstab:
195 iter = fasp_solver_dcsr_spbcgs(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
196 break;
197
198 case SOLVER_MinRes:
199 iter = fasp_solver_dcsr_spminres(A, b, x, pc, tol, MaxIt, stop_type, prtlvl);
200 break;
201
202 case SOLVER_GMRES:
203 iter = fasp_solver_dcsr_spgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
204 break;
205
206 case SOLVER_VGMRES:
207 iter = fasp_solver_dcsr_spvgmres(A, b, x, pc, tol, MaxIt, restart, stop_type, prtlvl);
208 break;
209
210 default:
211 printf("### ERROR: Unknown iterative solver type %d! [%s]\n",
212 itsolver_type, __FUNCTION__);
213 return ERROR_SOLVER_TYPE;
214
215 }
216
217 if ( (prtlvl >= PRINT_SOME) && (iter >= 0) ) {
218 fasp_gettime(&solve_end);
219 fasp_cputime("Iterative method", solve_end - solve_start);
220 }
221
222#if DEBUG_MODE > 0
223 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
224#endif
225
226 return iter;
227}
228
246 dvector *b,
247 dvector *x,
248 ITS_param *itparam)
249{
250 const SHORT prtlvl = itparam->print_level;
251
252 /* Local Variables */
253 INT status = FASP_SUCCESS;
254 REAL solve_start, solve_end;
255
256#if DEBUG_MODE > 0
257 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
258 printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
259 printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
260#endif
261
262 fasp_gettime(&solve_start);
263
264 status = fasp_solver_dcsr_itsolver(A,b,x,NULL,itparam);
265
266 if ( prtlvl >= PRINT_MIN ) {
267 fasp_gettime(&solve_end);
268 fasp_cputime("Krylov method totally", solve_end - solve_start);
269 }
270
271#if DEBUG_MODE > 0
272 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
273#endif
274
275 return status;
276}
277
295 dvector *b,
296 dvector *x,
297 ITS_param *itparam)
298{
299 const SHORT prtlvl = itparam->print_level;
300
301 /* Local Variables */
302 INT status = FASP_SUCCESS;
303 REAL solve_start, solve_end;
304
305#if DEBUG_MODE > 0
306 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
307 printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
308 printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
309#endif
310
311 fasp_gettime(&solve_start);
312
313 status = fasp_solver_dcsr_itsolver_s (A,b,x,NULL,itparam);
314
315 if ( prtlvl >= PRINT_MIN ) {
316 fasp_gettime(&solve_end);
317 fasp_cputime("Krylov method totally", solve_end - solve_start);
318 }
319
320#if DEBUG_MODE > 0
321 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
322#endif
323
324 return status;
325}
326
344 dvector *b,
345 dvector *x,
346 ITS_param *itparam)
347{
348 const SHORT prtlvl = itparam->print_level;
349
350 /* Local Variables */
351 INT status = FASP_SUCCESS;
352 REAL solve_start, solve_end;
353
354#if DEBUG_MODE > 0
355 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
356 printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
357 printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
358#endif
359
360 fasp_gettime(&solve_start);
361
362 // setup preconditioner
363 dvector diag; fasp_dcsr_getdiag(0,A,&diag);
364
365 precond pc;
366 pc.data = &diag;
368
369 // call iterative solver
370 status = fasp_solver_dcsr_itsolver(A,b,x,&pc,itparam);
371
372 if ( prtlvl >= PRINT_MIN ) {
373 fasp_gettime(&solve_end);
374 fasp_cputime("Diag_Krylov method totally", solve_end - solve_start);
375 }
376
377 fasp_dvec_free(&diag);
378
379#if DEBUG_MODE > 0
380 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
381#endif
382
383 return status;
384}
385
406 dvector *b,
407 dvector *x,
408 ITS_param *itparam,
409 SWZ_param *schparam)
410{
411 SWZ_param swzparam;
412 swzparam.SWZ_mmsize = schparam->SWZ_mmsize;
413 swzparam.SWZ_maxlvl = schparam->SWZ_maxlvl;
414 swzparam.SWZ_type = schparam->SWZ_type;
415 swzparam.SWZ_blksolver = schparam->SWZ_blksolver;
416
417 const SHORT prtlvl = itparam->print_level;
418
419 REAL setup_start, setup_end;
420 REAL solve_start, solve_end;
421 INT status = FASP_SUCCESS;
422
423#if DEBUG_MODE > 0
424 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
425 printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
426 printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
427#endif
428
429 fasp_gettime(&solve_start);
430 fasp_gettime(&setup_start);
431
432 // setup preconditioner
434
435 // symmetrize the matrix (get rid of this later)
437
438 // construct Schwarz precondtioner
440 fasp_swz_dcsr_setup (&SWZ_data, &swzparam);
441
442 fasp_gettime (&setup_end);
443 printf("SWZ_Krylov method setup %f seconds.\n", setup_end - setup_start);
444
445 precond prec;
446 prec.data = &SWZ_data;
447 prec.fct = fasp_precond_swz;
448
449 // solver part
450 status = fasp_solver_dcsr_itsolver(A,b,x,&prec,itparam);
451
452 if ( prtlvl > PRINT_NONE ) {
453 fasp_gettime(&solve_end);
454 fasp_cputime("SWZ_Krylov method totally", solve_end - solve_start);
455 }
456
457#if DEBUG_MODE > 0
458 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
459#endif
460
462
463 return status;
464}
465
484 dvector *b,
485 dvector *x,
486 ITS_param *itparam,
487 AMG_param *amgparam)
488{
489 const SHORT prtlvl = itparam->print_level;
490 const SHORT max_levels = amgparam->max_levels;
491 const INT nnz = A->nnz, m = A->row, n = A->col;
492
493 /* Local Variables */
494 INT status = FASP_SUCCESS;
495 REAL solve_start, solve_end;
496
497#if MULTI_COLOR_ORDER
498 A->color = 0;
499 A->IC = NULL;
500 A->ICMAP = NULL;
501#endif
502
503#if DEBUG_MODE > 0
504 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
505 printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
506 printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
507#endif
508
509 fasp_gettime(&solve_start);
510
511 // initialize A, b, x for mgl[0]
512 AMG_data *mgl=fasp_amg_data_create(max_levels);
513 mgl[0].A=fasp_dcsr_create(m,n,nnz); fasp_dcsr_cp(A,&mgl[0].A);
514 mgl[0].b=fasp_dvec_create(n); mgl[0].x=fasp_dvec_create(n);
515
516 // setup preconditioner
517 switch (amgparam->AMG_type) {
518
519 case SA_AMG: // Smoothed Aggregation AMG
520 status = fasp_amg_setup_sa(mgl, amgparam); break;
521
522 case UA_AMG: // Unsmoothed Aggregation AMG
523 status = fasp_amg_setup_ua(mgl, amgparam); break;
524
525 default: // Classical AMG
526 status = fasp_amg_setup_rs(mgl, amgparam);
527
528 }
529
530#if DEBUG_MODE > 1
532#endif
533
534 if (status < 0) goto FINISHED;
535
536 // setup preconditioner
537 precond_data pcdata;
538 fasp_param_amg_to_prec(&pcdata,amgparam);
539 pcdata.max_levels = mgl[0].num_levels;
540 pcdata.mgl_data = mgl;
541
542 precond pc; pc.data = &pcdata;
543
544 if (itparam->precond_type == PREC_FMG) {
545 pc.fct = fasp_precond_famg; // Full AMG
546 }
547 else {
548 switch (amgparam->cycle_type) {
549 case AMLI_CYCLE: // AMLI cycle
550 pc.fct = fasp_precond_amli; break;
551 case NL_AMLI_CYCLE: // Nonlinear AMLI
552 pc.fct = fasp_precond_namli; break;
553 default: // V,W-cycles or hybrid cycles
555 }
556 }
557
558 // call iterative solver
559 status = fasp_solver_dcsr_itsolver(A, b, x, &pc, itparam);
560
561 if ( prtlvl >= PRINT_MIN ) {
562 fasp_gettime(&solve_end);
563 fasp_cputime("AMG_Krylov method totally", solve_end - solve_start);
564 }
565
566FINISHED:
567 fasp_amg_data_free(mgl, amgparam);
568
569#if DEBUG_MODE > 0
570 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
571#endif
572
573 return status;
574}
575
594 dvector *b,
595 dvector *x,
596 ITS_param *itparam,
597 ILU_param *iluparam)
598{
599 const SHORT prtlvl = itparam->print_level;
600
601 /* Local Variables */
602 INT status = FASP_SUCCESS;
603 REAL solve_start, solve_end, solve_time;
604
605#if DEBUG_MODE > 0
606 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
607 printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
608 printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
609#endif
610
611 fasp_gettime(&solve_start);
612
613 // ILU setup for whole matrix
614 ILU_data LU;
615 if ( (status = fasp_ilu_dcsr_setup(A, &LU, iluparam)) < 0 ) goto FINISHED;
616
617 // check iludata
618 if ( (status = fasp_mem_iludata_check(&LU)) < 0 ) goto FINISHED;
619
620 // set preconditioner
621 precond pc;
622 pc.data = &LU;
624
625 // call iterative solver
626 status = fasp_solver_dcsr_itsolver(A, b, x, &pc, itparam);
627
628 if ( prtlvl >= PRINT_MIN ) {
629 fasp_gettime(&solve_end);
630 solve_time = solve_end - solve_start;
631
632 switch (iluparam->ILU_type) {
633 case ILUt:
634 fasp_cputime("ILUt_Krylov method totally", solve_time);
635 break;
636 case ILUtp:
637 fasp_cputime("ILUtp_Krylov method totally", solve_time);
638 break;
639 default: // ILUk
640 fasp_cputime("ILUk_Krylov method totally", solve_time);
641 }
642 }
643
644FINISHED:
646
647#if DEBUG_MODE > 0
648 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
649#endif
650
651 return status;
652}
653
677 dvector *b,
678 dvector *x,
679 ITS_param *itparam,
680 ILU_param *iluparam,
681 dCSRmat *M)
682{
683 const SHORT prtlvl = itparam->print_level;
684
685 /* Local Variables */
686 REAL solve_start, solve_end, solve_time;
687 INT status = FASP_SUCCESS;
688
689#if DEBUG_MODE > 0
690 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
691 printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
692 printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
693#endif
694
695 fasp_gettime(&solve_start);
696
697 // ILU setup for M
698 ILU_data LU;
699 if ( (status = fasp_ilu_dcsr_setup(M, &LU, iluparam)) < 0 ) goto FINISHED;
700
701 // check iludata
702 if ( (status = fasp_mem_iludata_check(&LU)) < 0 ) goto FINISHED;
703
704 // set precondtioner
705 precond pc;
706 pc.data = &LU;
708
709 // call iterative solver
710 status = fasp_solver_dcsr_itsolver(A, b, x, &pc, itparam);
711
712 if ( prtlvl >= PRINT_MIN ) {
713 fasp_gettime(&solve_end);
714 solve_time = solve_end - solve_start;
715
716 switch (iluparam->ILU_type) {
717 case ILUt:
718 fasp_cputime("ILUt_Krylov method", solve_time);
719 break;
720 case ILUtp:
721 fasp_cputime("ILUtp_Krylov method", solve_time);
722 break;
723 default: // ILUk
724 fasp_cputime("ILUk_Krylov method", solve_time);
725 }
726 }
727
728FINISHED:
730
731#if DEBUG_MODE > 0
732 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
733#endif
734
735 return status;
736}
737
760 dvector *b,
761 dvector *x,
762 ITS_param *itparam,
763 AMG_param *amgparam,
764 dCSRmat *A_nk,
765 dCSRmat *P_nk,
766 dCSRmat *R_nk)
767{
768 const SHORT prtlvl = itparam->print_level;
769 const SHORT max_levels = amgparam->max_levels;
770 const INT nnz = A->nnz, m = A->row, n = A->col;
771
772 /* Local Variables */
773 INT status = FASP_SUCCESS;
774 REAL solve_start, solve_end, solve_time;
775
776#if MULTI_COLOR_ORDER
777 A->color = 0;
778 A->IC = NULL;
779 A->ICMAP = NULL;
780#endif
781
782#if DEBUG_MODE > 0
783 printf("### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
784 printf("### DEBUG: matrix size: %d %d %d\n", A->row, A->col, A->nnz);
785 printf("### DEBUG: rhs/sol size: %d %d\n", b->row, x->row);
786#endif
787
788 fasp_gettime(&solve_start);
789
790 // initialize A, b, x for mgl[0]
791 AMG_data *mgl=fasp_amg_data_create(max_levels);
792 mgl[0].A=fasp_dcsr_create(m,n,nnz); fasp_dcsr_cp(A,&mgl[0].A);
793 mgl[0].b=fasp_dvec_create(n); mgl[0].x=fasp_dvec_create(n);
794
795 // setup preconditioner
796 switch (amgparam->AMG_type) {
797
798 case SA_AMG: // Smoothed Aggregation AMG
799 status = fasp_amg_setup_sa(mgl, amgparam); break;
800
801 case UA_AMG: // Unsmoothed Aggregation AMG
802 status = fasp_amg_setup_ua(mgl, amgparam); break;
803
804 default: // Classical AMG
805 status = fasp_amg_setup_rs(mgl, amgparam);
806
807 }
808
809#if DEBUG_MODE > 1
811#endif
812
813 if (status < 0) goto FINISHED;
814
815 // setup preconditioner
816 precond_data pcdata;
817 fasp_param_amg_to_prec(&pcdata,amgparam);
818 pcdata.max_levels = mgl[0].num_levels;
819 pcdata.mgl_data = mgl;
820
821 // near kernel space
822#if WITH_UMFPACK // use UMFPACK directly
823 dCSRmat A_tran;
824 A_tran = fasp_dcsr_create(A_nk->row, A_nk->col, A_nk->nnz);
825 fasp_dcsr_transz(A_nk, NULL, &A_tran);
826 // It is equivalent to do transpose and then sort
827 // fasp_dcsr_trans(A_nk, &A_tran);
828 // fasp_dcsr_sort(&A_tran);
829 pcdata.A_nk = &A_tran;
830#else
831 pcdata.A_nk = A_nk;
832#endif
833
834 pcdata.P_nk = P_nk;
835 pcdata.R_nk = R_nk;
836
837 precond pc; pc.data = &pcdata;
839
840 // call iterative solver
841 status = fasp_solver_dcsr_itsolver(A, b, x, &pc, itparam);
842
843 if ( prtlvl >= PRINT_MIN ) {
844 fasp_gettime(&solve_end);
845 solve_time = solve_end - solve_start;
846 fasp_cputime("AMG_NK_Krylov method", solve_time);
847 }
848
849FINISHED:
850 fasp_amg_data_free(mgl, amgparam);
851
852#if WITH_UMFPACK // use UMFPACK directly
853 fasp_dcsr_free(&A_tran);
854#endif
855
856#if DEBUG_MODE > 0
857 printf("### DEBUG: [--End--] %s ...\n", __FUNCTION__);
858#endif
859
860 return status;
861}
862
863/*---------------------------------*/
864/*-- End of File --*/
865/*---------------------------------*/
void fasp_mem_usage(void)
Show total allocated memory currently.
Definition: AuxMemory.c:185
SHORT fasp_mem_iludata_check(const ILU_data *iludata)
Check wether a ILU_data has enough work space.
Definition: AuxMemory.c:205
void fasp_cputime(const char *message, const REAL cputime)
Print CPU walltime.
Definition: AuxMessage.c:179
void fasp_param_amg_to_prec(precond_data *pcdata, const AMG_param *amgparam)
Set precond_data with AMG_param.
Definition: AuxParam.c:687
void fasp_gettime(REAL *time)
Get system time.
Definition: AuxTiming.c:36
void fasp_dvec_free(dvector *u)
Free vector data space of REAL type.
Definition: AuxVector.c:145
dvector fasp_dvec_create(const INT m)
Create dvector data space of REAL type.
Definition: AuxVector.c:62
SHORT fasp_ilu_dcsr_setup(dCSRmat *A, ILU_data *iludata, ILU_param *iluparam)
Get ILU decomposition of a CSR matrix A.
INT fasp_swz_dcsr_setup(SWZ_data *swzdata, SWZ_param *swzparam)
Setup phase for the Schwarz methods.
dCSRmat fasp_dcsr_create(const INT m, const INT n, const INT nnz)
Create CSR sparse matrix data memory space.
Definition: BlaSparseCSR.c:47
void fasp_dcsr_shift(dCSRmat *A, const INT offset)
Re-index a REAL matrix in CSR format to make the index starting from 0 or 1.
void fasp_dcsr_free(dCSRmat *A)
Free CSR sparse matrix data memory space.
Definition: BlaSparseCSR.c:184
void fasp_dcsr_transz(dCSRmat *A, INT *p, dCSRmat *AT)
Generalized transpose of A: (n x m) matrix given in dCSRmat format.
void fasp_dcsr_cp(const dCSRmat *A, dCSRmat *B)
copy a dCSRmat to a new one B=A
Definition: BlaSparseCSR.c:851
void fasp_dcsr_getdiag(INT n, const dCSRmat *A, dvector *diag)
Get first n diagonal entries of a CSR matrix A.
Definition: BlaSparseCSR.c:537
dCSRmat fasp_dcsr_sympart(dCSRmat *A)
Get symmetric part of a dCSRmat matrix.
void fasp_check_dCSRmat(const dCSRmat *A)
Check whether an dCSRmat matrix is supported or not.
INT fasp_solver_dcsr_pbcgs(dCSRmat *A, dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT StopType, const SHORT PrtLvl)
Preconditioned BiCGstab method for solving Au=b for CSR matrix.
Definition: KryPbcgs.c:62
INT fasp_solver_dcsr_pcg(dCSRmat *A, dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT StopType, const SHORT PrtLvl)
Preconditioned conjugate gradient method for solving Au=b.
Definition: KryPcg.c:98
INT fasp_solver_dcsr_pgcg(dCSRmat *A, dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT StopType, const SHORT PrtLvl)
Preconditioned generilzed conjugate gradient (GCG) method for solving Au=b.
Definition: KryPgcg.c:60
INT fasp_solver_dcsr_pgcr(dCSRmat *A, dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, const SHORT restart, const SHORT StopType, const SHORT PrtLvl)
A preconditioned GCR method for solving Au=b.
Definition: KryPgcr.c:55
INT fasp_solver_dcsr_pgmres(dCSRmat *A, dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, const SHORT restart, const SHORT StopType, const SHORT PrtLvl)
Right preconditioned GMRES method for solving Au=b.
Definition: KryPgmres.c:67
INT fasp_solver_dcsr_pminres(dCSRmat *A, dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT StopType, const SHORT PrtLvl)
A preconditioned minimal residual (Minres) method for solving Au=b.
Definition: KryPminres.c:62
INT fasp_solver_dcsr_pvfgmres(dCSRmat *A, dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, const SHORT restart, const SHORT StopType, const SHORT PrtLvl)
Solve "Ax=b" using PFGMRES(right preconditioned) iterative method in which the restart parameter can ...
Definition: KryPvfgmres.c:67
INT fasp_solver_dcsr_pvgmres(dCSRmat *A, dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, const SHORT restart, const SHORT StopType, const SHORT PrtLvl)
Right preconditioned GMRES method in which the restart parameter can be adaptively modified during it...
Definition: KryPvgmres.c:66
INT fasp_solver_dcsr_spbcgs(const dCSRmat *A, const dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT StopType, const SHORT PrtLvl)
Preconditioned BiCGstab method for solving Au=b with safety net.
Definition: KrySPbcgs.c:61
INT fasp_solver_dcsr_spcg(const dCSRmat *A, const dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT StopType, const SHORT PrtLvl)
Preconditioned conjugate gradient method for solving Au=b with safety net.
Definition: KrySPcg.c:60
INT fasp_solver_dcsr_spgmres(const dCSRmat *A, const dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, SHORT restart, const SHORT StopType, const SHORT PrtLvl)
Preconditioned GMRES method for solving Au=b with safe-guard.
Definition: KrySPgmres.c:66
INT fasp_solver_dcsr_spminres(const dCSRmat *A, const dvector *b, dvector *u, precond *pc, const REAL tol, const INT MaxIt, const SHORT StopType, const SHORT PrtLvl)
A preconditioned minimal residual (Minres) method for solving Au=b with safety net.
Definition: KrySPminres.c:60
INT fasp_solver_dcsr_spvgmres(const dCSRmat *A, const dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, SHORT restart, const SHORT StopType, const SHORT PrtLvl)
Solve "Ax=b" using PGMRES(right preconditioned) iterative method in which the restart parameter can b...
Definition: KrySPvgmres.c:68
SHORT fasp_amg_setup_rs(AMG_data *mgl, AMG_param *param)
Setup phase of Ruge and Stuben's classic AMG.
Definition: PreAMGSetupRS.c:52
SHORT fasp_amg_setup_sa(AMG_data *mgl, AMG_param *param)
Set up phase of smoothed aggregation AMG.
Definition: PreAMGSetupSA.c:63
SHORT fasp_amg_setup_ua(AMG_data *mgl, AMG_param *param)
Set up phase of unsmoothed aggregation AMG.
Definition: PreAMGSetupUA.c:55
void fasp_precond_amg(REAL *r, REAL *z, void *data)
AMG preconditioner.
Definition: PreCSR.c:416
void fasp_precond_ilu(REAL *r, REAL *z, void *data)
ILU preconditioner.
Definition: PreCSR.c:198
void fasp_precond_swz(REAL *r, REAL *z, void *data)
get z from r by Schwarz
Definition: PreCSR.c:371
void fasp_precond_namli(REAL *r, REAL *z, void *data)
Nonlinear AMLI AMG preconditioner.
Definition: PreCSR.c:515
void fasp_precond_amg_nk(REAL *r, REAL *z, void *data)
AMG with extra near kernel solve as preconditioner.
Definition: PreCSR.c:548
void fasp_precond_diag(REAL *r, REAL *z, void *data)
Diagonal preconditioner z=inv(D)*r.
Definition: PreCSR.c:172
void fasp_precond_famg(REAL *r, REAL *z, void *data)
Full AMG preconditioner.
Definition: PreCSR.c:449
void fasp_precond_amli(REAL *r, REAL *z, void *data)
AMLI AMG preconditioner.
Definition: PreCSR.c:482
void fasp_swz_data_free(SWZ_data *swzdata)
Free SWZ_data data memeory space.
Definition: PreDataInit.c:341
AMG_data * fasp_amg_data_create(SHORT max_levels)
Create and initialize AMG_data for classical and SA AMG.
Definition: PreDataInit.c:64
void fasp_amg_data_free(AMG_data *mgl, AMG_param *param)
Free AMG_data data memeory space.
Definition: PreDataInit.c:101
void fasp_ilu_data_free(ILU_data *iludata)
Create ILU_data sturcture.
Definition: PreDataInit.c:300
INT fasp_solver_dcsr_itsolver(dCSRmat *A, dvector *b, dvector *x, precond *pc, ITS_param *itparam)
Solve Ax=b by preconditioned Krylov methods for CSR matrices.
Definition: SolCSR.c:56
INT fasp_solver_dcsr_krylov_amg_nk(dCSRmat *A, dvector *b, dvector *x, ITS_param *itparam, AMG_param *amgparam, dCSRmat *A_nk, dCSRmat *P_nk, dCSRmat *R_nk)
Solve Ax=b by AMG preconditioned Krylov methods with an extra near kernel solve.
Definition: SolCSR.c:759
INT fasp_solver_dcsr_krylov_ilu_M(dCSRmat *A, dvector *b, dvector *x, ITS_param *itparam, ILU_param *iluparam, dCSRmat *M)
Solve Ax=b by ILUs preconditioned Krylov methods: ILU of M as preconditioner.
Definition: SolCSR.c:676
INT fasp_solver_dcsr_krylov_ilu(dCSRmat *A, dvector *b, dvector *x, ITS_param *itparam, ILU_param *iluparam)
Solve Ax=b by ILUs preconditioned Krylov methods.
Definition: SolCSR.c:593
INT fasp_solver_dcsr_krylov_diag(dCSRmat *A, dvector *b, dvector *x, ITS_param *itparam)
Solve Ax=b by diagonal preconditioned Krylov methods.
Definition: SolCSR.c:343
INT fasp_solver_dcsr_krylov_s(dCSRmat *A, dvector *b, dvector *x, ITS_param *itparam)
Solve Ax=b by standard Krylov methods with safe-net for CSR matrices.
Definition: SolCSR.c:294
INT fasp_solver_dcsr_krylov(dCSRmat *A, dvector *b, dvector *x, ITS_param *itparam)
Solve Ax=b by standard Krylov methods for CSR matrices.
Definition: SolCSR.c:245
INT fasp_solver_dcsr_krylov_swz(dCSRmat *A, dvector *b, dvector *x, ITS_param *itparam, SWZ_param *schparam)
Solve Ax=b by overlapping Schwarz Krylov methods.
Definition: SolCSR.c:405
INT fasp_solver_dcsr_itsolver_s(dCSRmat *A, dvector *b, dvector *x, precond *pc, ITS_param *itparam)
Solve Ax=b by preconditioned Krylov methods with safe-net for CSR matrices.
Definition: SolCSR.c:158
INT fasp_solver_dcsr_krylov_amg(dCSRmat *A, dvector *b, dvector *x, ITS_param *itparam, AMG_param *amgparam)
Solve Ax=b by AMG preconditioned Krylov methods.
Definition: SolCSR.c:483
Main header file for the FASP project.
#define REAL
Definition: fasp.h:67
#define SHORT
FASP integer and floating point numbers.
Definition: fasp.h:63
#define INT
Definition: fasp.h:64
#define AMLI_CYCLE
Definition: fasp_const.h:179
#define NL_AMLI_CYCLE
Definition: fasp_const.h:180
#define SOLVER_GMRES
Definition: fasp_const.h:106
#define FASP_SUCCESS
Definition of return status and error messages.
Definition: fasp_const.h:19
#define PREC_FMG
Definition: fasp_const.h:141
#define ILUtp
Definition: fasp_const.h:150
#define ERROR_SOLVER_TYPE
Definition: fasp_const.h:41
#define SOLVER_VGMRES
Definition: fasp_const.h:107
#define SOLVER_BiCGstab
Definition: fasp_const.h:104
#define ILUt
Definition: fasp_const.h:149
#define SA_AMG
Definition: fasp_const.h:163
#define PRINT_SOME
Definition: fasp_const.h:75
#define SOLVER_MinRes
Definition: fasp_const.h:105
#define SOLVER_GCR
Definition: fasp_const.h:110
#define SOLVER_VFGMRES
Definition: fasp_const.h:108
#define SOLVER_CG
Definition: fasp_const.h:103
#define PRINT_NONE
Print level for all subroutines – not including DEBUG output.
Definition: fasp_const.h:73
#define PRINT_MIN
Definition: fasp_const.h:74
#define SOLVER_GCG
Definition: fasp_const.h:109
#define UA_AMG
Definition: fasp_const.h:164
Data for AMG methods.
Definition: fasp.h:790
dCSRmat A
pointer to the matrix at level level_num
Definition: fasp.h:803
dvector b
pointer to the right-hand side at level level_num
Definition: fasp.h:812
dvector x
pointer to the iterative solution at level level_num
Definition: fasp.h:815
SHORT num_levels
number of levels in use <= max_levels
Definition: fasp.h:798
Parameters for AMG methods.
Definition: fasp.h:447
SHORT AMG_type
type of AMG method
Definition: fasp.h:450
SHORT cycle_type
type of AMG cycle
Definition: fasp.h:468
SHORT max_levels
max number of levels of AMG
Definition: fasp.h:462
Data for ILU setup.
Definition: fasp.h:637
Parameters for ILU.
Definition: fasp.h:396
SHORT ILU_type
ILU type for decomposition.
Definition: fasp.h:402
Parameters for iterative solvers.
Definition: fasp.h:379
SHORT itsolver_type
Definition: fasp.h:382
SHORT print_level
Definition: fasp.h:381
SHORT precond_type
Definition: fasp.h:384
REAL tol
Definition: fasp.h:388
INT restart
Definition: fasp.h:386
SHORT stop_type
Definition: fasp.h:385
INT maxit
Definition: fasp.h:387
Data for Schwarz methods.
Definition: fasp.h:712
dCSRmat A
pointer to the original coefficient matrix
Definition: fasp.h:717
Parameters for Schwarz method.
Definition: fasp.h:422
INT SWZ_mmsize
maximal size of blocks
Definition: fasp.h:434
INT SWZ_maxlvl
maximal level for constructing the blocks
Definition: fasp.h:431
INT SWZ_blksolver
type of Schwarz block solver
Definition: fasp.h:437
SHORT SWZ_type
type for Schwarz method
Definition: fasp.h:428
Sparse matrix of REAL type in CSR format.
Definition: fasp.h:143
INT col
column of matrix A, n
Definition: fasp.h:149
INT row
row number of matrix A, m
Definition: fasp.h:146
INT nnz
number of nonzero entries
Definition: fasp.h:152
Vector with n entries of REAL type.
Definition: fasp.h:346
INT row
number of rows
Definition: fasp.h:349
Data for preconditioners.
Definition: fasp.h:880
AMG_data * mgl_data
AMG preconditioner data.
Definition: fasp.h:940
dCSRmat * A_nk
Matrix data for near kernel.
Definition: fasp.h:951
dCSRmat * R_nk
Restriction for near kernel.
Definition: fasp.h:957
SHORT max_levels
max number of AMG levels
Definition: fasp.h:892
dCSRmat * P_nk
Prolongation for near kernel.
Definition: fasp.h:954
Preconditioner data and action.
Definition: fasp.h:1081
void * data
data for preconditioner, void pointer
Definition: fasp.h:1084
void(* fct)(REAL *, REAL *, void *)
action for preconditioner, void function pointer
Definition: fasp.h:1087