26#include "fasp_functs.h"
77 INT iter = 0, stag = 1, more_step = 1, restart_step = 1;
80 REAL normu2, normuu, normp, infnormu, factor;
81 REAL alpha, alpha0, alpha1, temp2;
85 REAL *p0=work, *p1=work+m, *p2=p1+m, *z0=p2+m, *z1=z0+m;
86 REAL *t0=z1+m, *t1=t0+m, *t=t1+m, *tp=t+m, *tz=tp+m, *r=tz+m;
89 if ( PrtLvl >
PRINT_NONE ) printf(
"\nCalling MinRes solver (CSR) ...\n");
92 printf(
"### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
93 printf(
"### DEBUG: maxit = %d, tol = %.4le\n", MaxIt, tol);
110 switch ( StopType ) {
114 relres = absres0/normr0;
119 relres = absres0/normr0;
124 relres = absres0/normu2;
127 printf(
"### ERROR: Unknown stopping type! [%s]\n", __FUNCTION__);
132 if ( relres < tol || absres0 < 1e-12*tol )
goto FINISHED;
135 fasp_itinfo(PrtLvl,StopType,iter,relres,absres0,0.0);
164 while ( iter++ < MaxIt ) {
219 switch ( StopType ) {
222 absres = sqrt(temp2);
223 relres = absres/normr0;
231 absres = sqrt(temp2);
232 relres = absres/normr0;
236 absres = sqrt(temp2);
237 relres = absres/normu2;
242 factor = absres/absres0;
245 fasp_itinfo(PrtLvl,StopType,iter,relres,absres,factor);
247 if ( factor > 0.9 ) {
251 if (infnormu <= sol_inf_tol) {
259 normuu =
ABS(alpha)*(normuu/normu2);
261 if ( normuu < maxdiff ) {
263 if ( stag < MaxStag ) {
265 ITS_DIFFRES(normuu,relres);
277 absres = sqrt(temp2);
278 relres = absres/normr0;
286 absres = sqrt(temp2);
287 relres = absres/normr0;
291 absres = sqrt(temp2);
292 relres = absres/normu2;
296 if ( PrtLvl >=
PRINT_MORE ) ITS_REALRES(relres);
301 if ( stag >= MaxStag ) {
348 if ( relres < tol ) {
350 if ( PrtLvl >=
PRINT_MORE ) ITS_COMPRES(relres);
360 absres = sqrt(temp2);
361 relres = absres/normr0;
369 absres = sqrt(temp2);
370 relres = absres/normr0;
374 absres = sqrt(temp2);
375 relres = absres/normu2;
379 if ( PrtLvl >=
PRINT_MORE ) ITS_REALRES(relres);
382 if ( relres < tol )
break;
384 if ( more_step >= MaxRestartStep ) {
436 if ( PrtLvl >
PRINT_NONE ) ITS_FINAL(iter,MaxIt,relres);
442 printf(
"### DEBUG: [--End--] %s ...\n", __FUNCTION__);
481 const SHORT StopType,
490 INT iter = 0, stag = 1, more_step = 1, restart_step = 1;
493 REAL normu2, normuu, normp, infnormu, factor;
494 REAL alpha, alpha0, alpha1, temp2;
498 REAL *p0=work, *p1=work+m, *p2=p1+m, *z0=p2+m, *z1=z0+m;
499 REAL *t0=z1+m, *t1=t0+m, *t=t1+m, *tp=t+m, *tz=tp+m, *r=tz+m;
502 if ( PrtLvl >
PRINT_NONE ) printf(
"\nCalling MinRes solver (BLC) ...\n");
505 printf(
"### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
506 printf(
"### DEBUG: maxit = %d, tol = %.4le\n", MaxIt, tol);
523 switch ( StopType ) {
527 relres = absres0/normr0;
532 relres = absres0/normr0;
537 relres = absres0/normu2;
540 printf(
"### ERROR: Unknown stopping type! [%s]\n", __FUNCTION__);
545 if ( relres < tol || absres0 < 1e-12*tol )
goto FINISHED;
548 fasp_itinfo(PrtLvl,StopType,iter,relres,absres0,0.0);
577 while ( iter++ < MaxIt ) {
632 switch ( StopType ) {
635 absres = sqrt(temp2);
636 relres = absres/normr0;
644 absres = sqrt(temp2);
645 relres = absres/normr0;
649 absres = sqrt(temp2);
650 relres = absres/normu2;
655 factor = absres/absres0;
658 fasp_itinfo(PrtLvl,StopType,iter,relres,absres,factor);
660 if ( factor > 0.9 ) {
664 if (infnormu <= sol_inf_tol) {
672 normuu =
ABS(alpha)*(normuu/normu2);
674 if ( normuu < maxdiff ) {
676 if ( stag < MaxStag ) {
678 ITS_DIFFRES(normuu,relres);
690 absres = sqrt(temp2);
691 relres = absres/normr0;
699 absres = sqrt(temp2);
700 relres = absres/normr0;
704 absres = sqrt(temp2);
705 relres = absres/normu2;
709 if ( PrtLvl >=
PRINT_MORE ) ITS_REALRES(relres);
714 if ( stag >= MaxStag ) {
761 if ( relres < tol ) {
763 if ( PrtLvl >=
PRINT_MORE ) ITS_COMPRES(relres);
773 absres = sqrt(temp2);
774 relres = absres/normr0;
782 absres = sqrt(temp2);
783 relres = absres/normr0;
787 absres = sqrt(temp2);
788 relres = absres/normu2;
792 if ( PrtLvl >=
PRINT_MORE ) ITS_REALRES(relres);
795 if ( relres < tol )
break;
797 if ( more_step >= MaxRestartStep ) {
849 if ( PrtLvl >
PRINT_NONE ) ITS_FINAL(iter,MaxIt,relres);
855 printf(
"### DEBUG: [--End--] %s ...\n", __FUNCTION__);
891 const SHORT StopType,
900 INT iter = 0, stag = 1, more_step = 1, restart_step = 1;
903 REAL normu2, normuu, normp, infnormu, factor;
904 REAL alpha, alpha0, alpha1, temp2;
908 REAL *p0=work, *p1=work+m, *p2=p1+m, *z0=p2+m, *z1=z0+m;
909 REAL *t0=z1+m, *t1=t0+m, *t=t1+m, *tp=t+m, *tz=tp+m, *r=tz+m;
912 if ( PrtLvl >
PRINT_NONE ) printf(
"\nCalling MinRes solver (STR) ...\n");
915 printf(
"### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
916 printf(
"### DEBUG: maxit = %d, tol = %.4le\n", MaxIt, tol);
933 switch ( StopType ) {
937 relres = absres0/normr0;
942 relres = absres0/normr0;
947 relres = absres0/normu2;
950 printf(
"### ERROR: Unknown stopping type! [%s]\n", __FUNCTION__);
955 if ( relres < tol || absres0 < 1e-12*tol )
goto FINISHED;
958 fasp_itinfo(PrtLvl,StopType,iter,relres,absres0,0.0);
987 while ( iter++ < MaxIt ) {
1022 normp = sqrt(normp);
1042 switch ( StopType ) {
1045 absres = sqrt(temp2);
1046 relres = absres/normr0;
1054 absres = sqrt(temp2);
1055 relres = absres/normr0;
1059 absres = sqrt(temp2);
1060 relres = absres/normu2;
1065 factor = absres/absres0;
1068 fasp_itinfo(PrtLvl,StopType,iter,relres,absres,factor);
1070 if ( factor > 0.9 ) {
1074 if (infnormu <= sol_inf_tol) {
1082 normuu =
ABS(alpha)*(normuu/normu2);
1084 if ( normuu < maxdiff ) {
1086 if ( stag < MaxStag ) {
1088 ITS_DIFFRES(normuu,relres);
1100 absres = sqrt(temp2);
1101 relres = absres/normr0;
1109 absres = sqrt(temp2);
1110 relres = absres/normr0;
1114 absres = sqrt(temp2);
1115 relres = absres/normu2;
1119 if ( PrtLvl >=
PRINT_MORE ) ITS_REALRES(relres);
1124 if ( stag >= MaxStag ) {
1150 normp = sqrt(normp);
1170 if ( relres < tol ) {
1172 if ( PrtLvl >=
PRINT_MORE ) ITS_COMPRES(relres);
1182 absres = sqrt(temp2);
1183 relres = absres/normr0;
1191 absres = sqrt(temp2);
1192 relres = absres/normr0;
1196 absres = sqrt(temp2);
1197 relres = absres/normu2;
1201 if ( PrtLvl >=
PRINT_MORE ) ITS_REALRES(relres);
1204 if ( relres < tol )
break;
1206 if ( more_step >= MaxRestartStep ) {
1234 normp = sqrt(normp);
1258 if ( PrtLvl >
PRINT_NONE ) ITS_FINAL(iter,MaxIt,relres);
1264 printf(
"### DEBUG: [--End--] %s ...\n", __FUNCTION__);
1302 const SHORT StopType,
1311 INT iter = 0, stag = 1, more_step = 1, restart_step = 1;
1314 REAL normu2, normuu, normp, infnormu, factor;
1315 REAL alpha, alpha0, alpha1, temp2;
1319 REAL *p0=work, *p1=work+m, *p2=p1+m, *z0=p2+m, *z1=z0+m;
1320 REAL *t0=z1+m, *t1=t0+m, *t=t1+m, *tp=t+m, *tz=tp+m, *r=tz+m;
1323 if ( PrtLvl >
PRINT_NONE ) printf(
"\nCalling MinRes solver (MatFree) ...\n");
1326 printf(
"### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
1327 printf(
"### DEBUG: maxit = %d, tol = %.4le\n", MaxIt, tol);
1331 stag=1; more_step=1; restart_step=1;
1351 relres=absres0/normr0;
1356 relres=absres0/normu2;
1361 relres=absres0/normr0;
1366 if ( relres < tol || absres0 < 1e-12*tol )
goto FINISHED;
1394 while( iter++ < MaxIt) {
1459 relres=sqrt(temp2)/normr0;
1462 relres=sqrt(temp2)/normu2;
1465 relres=sqrt(temp2)/normr0;
1470 factor=absres/absres0;
1473 fasp_itinfo(PrtLvl,StopType,iter,relres,absres,factor);
1477 if ( infnormu <= sol_inf_tol ) {
1484 normuu=
ABS(alpha)*(normuu/normu2);
1487 if (normuu<maxdiff) {
1488 if ( stag < MaxStag ) {
1490 ITS_DIFFRES(normuu,relres);
1502 relres=sqrt(temp2)/normr0;
1510 relres=sqrt(temp2)/normr0;
1513 relres=sqrt(temp2)/normu2;
1517 if ( PrtLvl >=
PRINT_MORE ) ITS_REALRES(relres);
1522 if ( stag >= MaxStag ) {
1568 if ( relres < tol ) {
1569 if ( PrtLvl >=
PRINT_MORE ) ITS_COMPRES(relres);
1578 relres=sqrt(temp2)/normr0;
1586 relres=sqrt(temp2)/normr0;
1589 relres=sqrt(temp2)/normu2;
1593 if ( PrtLvl >=
PRINT_MORE ) ITS_REALRES(relres);
1596 if ( relres < tol )
break;
1598 if ( more_step >= MaxRestartStep ) {
1604 if ( more_step < MaxRestartStep ) {
1653 if ( PrtLvl >
PRINT_NONE ) ITS_FINAL(iter,MaxIt,relres);
1659 printf(
"### DEBUG: [--End--] %s ...\n", __FUNCTION__);
void fasp_darray_set(const INT n, REAL *x, const REAL val)
Set initial value for an array to be x=val.
void fasp_darray_cp(const INT n, const REAL *x, REAL *y)
Copy an array to the other y=x.
void fasp_mem_free(void *mem)
Free up previous allocated memory body and set pointer to NULL.
void * fasp_mem_calloc(const unsigned int size, const unsigned int type)
Allocate, initiate, and check memory.
void fasp_itinfo(const INT ptrlvl, const INT stop_type, const INT iter, const REAL relres, const REAL absres, const REAL factor)
Print out iteration information for iterative solvers.
REAL fasp_blas_darray_dotprod(const INT n, const REAL *x, const REAL *y)
Inner product of two arraies x and y.
REAL fasp_blas_darray_norminf(const INT n, const REAL *x)
Linf norm of array x.
void fasp_blas_darray_axpby(const INT n, const REAL a, const REAL *x, const REAL b, REAL *y)
y = a*x + b*y
REAL fasp_blas_darray_norm2(const INT n, const REAL *x)
L2 norm of array x.
void fasp_blas_darray_axpy(const INT n, const REAL a, const REAL *x, REAL *y)
y = a*x + y
void fasp_blas_dblc_mxv(const dBLCmat *A, const REAL *x, REAL *y)
Matrix-vector multiplication y = A*x.
void fasp_blas_dblc_aAxpy(const REAL alpha, const dBLCmat *A, const REAL *x, REAL *y)
Matrix-vector multiplication y = alpha*A*x + y.
void fasp_blas_dcsr_mxv(const dCSRmat *A, const REAL *x, REAL *y)
Matrix-vector multiplication y = A*x.
void fasp_blas_dcsr_aAxpy(const REAL alpha, const dCSRmat *A, const REAL *x, REAL *y)
Matrix-vector multiplication y = alpha*A*x + y.
void fasp_blas_dstr_aAxpy(const REAL alpha, const dSTRmat *A, const REAL *x, REAL *y)
Matrix-vector multiplication y = alpha*A*x + y.
void fasp_blas_dstr_mxv(const dSTRmat *A, const REAL *x, REAL *y)
Matrix-vector multiplication y = A*x.
INT fasp_solver_dstr_pminres(dSTRmat *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.
INT fasp_solver_dblc_pminres(dBLCmat *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.
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.
INT fasp_solver_pminres(mxv_matfree *mf, 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.
Main header file for the FASP project.
#define SHORT
FASP integer and floating point numbers.
#define MAX(a, b)
Definition of max, min, abs.
#define ERROR_SOLVER_STAG
#define ERROR_SOLVER_MAXIT
#define STOP_REL_RES
Definition of iterative solver stopping criteria types.
#define ERROR_SOLVER_SOLSTAG
#define BIGREAL
Some global constants.
#define PRINT_NONE
Print level for all subroutines – not including DEBUG output.
#define ERROR_SOLVER_TOLSMALL
Block REAL CSR matrix format.
Sparse matrix of REAL type in CSR format.
Structure matrix of REAL type.
Vector with n entries of REAL type.
REAL * val
actual vector entries
Matrix-vector multiplication, replace the actual matrix.
void * data
data for MxV, can be a Matrix or something else
void(* fct)(const void *, const REAL *, REAL *)
action for MxV, void function pointer
Preconditioner data and action.
void * data
data for preconditioner, void pointer
void(* fct)(REAL *, REAL *, void *)
action for preconditioner, void function pointer