28#include "fasp_functs.h"
77 const INT MIN_ITER = 0;
84 const REAL cr_max = 0.99;
85 const REAL cr_min = 0.174;
91 REAL r_norm, r_normb, gamma, t;
96 REAL r_norm_old = 0.0;
98 INT restart_max = restart;
101 INT Restart = restart;
102 INT Restart1 = Restart + 1;
103 unsigned LONG worksize = (Restart+4)*(Restart+n)+1-n;
106 REAL *c = NULL, *s = NULL, *rs = NULL;
107 REAL *norms = NULL, *r = NULL, *w = NULL;
109 REAL **p = NULL, **hh = NULL;
112 if ( PrtLvl >
PRINT_NONE ) printf(
"\nCalling VGMRes solver (CSR) ...\n");
115 printf(
"### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
116 printf(
"### DEBUG: maxit = %d, tol = %.4le\n", MaxIt, tol);
123 while ( (work == NULL) && (Restart > 5) ) {
124 Restart = Restart - 5;
125 worksize = (Restart+4)*(Restart+n)+1-n;
127 Restart1 = Restart + 1;
130 if ( work == NULL ) {
131 printf(
"### ERROR: No enough memory! [%s:%d]\n", __FILE__, __LINE__ );
135 if ( PrtLvl >
PRINT_MIN && Restart < restart ) {
136 printf(
"### WARNING: vGMRES restart number set to %d!\n", Restart);
143 r = work; w = r + n; rs = w + n; c = rs + Restart1; s = c + Restart;
145 for ( i = 0; i < Restart1; i++ ) p[i] = s + Restart + i*n;
147 for ( i = 0; i < Restart1; i++ ) hh[i] = p[Restart] + n + i*Restart;
159 relres = r_norm/absres0;
168 relres = r_normb/absres0;
173 relres = absres0/normu;
176 printf(
"### ERROR: Unknown stopping type! [%s]\n", __FUNCTION__);
181 if ( relres < tol || absres0 < 1e-12*tol )
goto FINISHED;
190 while ( iter < MaxIt ) {
192 rs[0] = r_norm_old = r_norm;
201 if ( cr > cr_max || iter == 0 ) {
202 Restart = restart_max;
204 else if ( cr < cr_min ) {
208 if ( Restart - d > restart_min ) {
212 Restart = restart_max;
218 while ( i < Restart && iter < MaxIt ) {
231 for (j = 0; j < i; j ++) {
242 for (j = 1; j < i; ++j) {
244 hh[j-1][i-1] = s[j-1]*hh[j][i-1] + c[j-1]*t;
245 hh[j][i-1] = -s[j-1]*t + c[j-1]*hh[j][i-1];
247 t= hh[i][i-1]*hh[i][i-1];
248 t+= hh[i-1][i-1]*hh[i-1][i-1];
251 if (gamma == 0.0) gamma = epsmac;
252 c[i-1] = hh[i-1][i-1] / gamma;
253 s[i-1] = hh[i][i-1] / gamma;
254 rs[i] = -s[i-1]*rs[i-1];
255 rs[i-1] = c[i-1]*rs[i-1];
256 hh[i-1][i-1] = s[i-1]*hh[i][i-1] + c[i-1]*hh[i-1][i-1];
258 absres = r_norm = fabs(rs[i]);
260 relres = absres/absres0;
262 norms[iter] = relres;
265 fasp_itinfo(PrtLvl, StopType, iter, relres, absres,
266 norms[iter]/norms[iter-1]);
269 if ( relres < tol && iter >= MIN_ITER )
break;
274 rs[i-1] = rs[i-1] / hh[i-1][i-1];
275 for (k = i-2; k >= 0; k --) {
277 for (j = k+1; j < i; j ++) t -= hh[k][j]*rs[j];
280 rs[k] = t / hh[k][k];
298 if ( relres < tol && iter >= MIN_ITER ) {
300 REAL computed_relres = relres;
308 switch ( StopType ) {
311 relres = absres/absres0;
319 relres = absres/absres0;
324 relres = absres/normu;
328 norms[iter] = relres;
330 if ( relres < tol ) {
339 ITS_COMPRES(computed_relres); ITS_REALRES(relres);
345 for ( j = i; j > 0; j-- ) {
346 rs[j-1] = -s[j-1]*rs[j];
347 rs[j] = c[j-1]*rs[j];
362 cr = r_norm / r_norm_old;
367 if ( PrtLvl >
PRINT_NONE ) ITS_FINAL(iter,MaxIt,relres);
378 printf(
"### DEBUG: [--End--] %s ...\n", __FUNCTION__);
420 const SHORT StopType,
424 const INT MIN_ITER = 0;
431 const REAL cr_max = 0.99;
432 const REAL cr_min = 0.174;
438 REAL r_norm, r_normb, gamma, t;
443 REAL r_norm_old = 0.0;
445 INT restart_max = restart;
448 INT Restart = restart;
449 INT Restart1 = Restart + 1;
450 unsigned LONG worksize = (Restart+4)*(Restart+n)+1-n;
453 REAL *c = NULL, *s = NULL, *rs = NULL;
454 REAL *norms = NULL, *r = NULL, *w = NULL;
456 REAL **p = NULL, **hh = NULL;
459 if ( PrtLvl >
PRINT_NONE ) printf(
"\nCalling VGMRes solver (BSR) ...\n");
462 printf(
"### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
463 printf(
"### DEBUG: maxit = %d, tol = %.4le\n", MaxIt, tol);
470 while ( (work == NULL) && (Restart > 5) ) {
471 Restart = Restart - 5;
472 worksize = (Restart+4)*(Restart+n)+1-n;
474 Restart1 = Restart + 1;
477 if ( work == NULL ) {
478 printf(
"### ERROR: No enough memory! [%s:%d]\n", __FILE__, __LINE__ );
482 if ( PrtLvl >
PRINT_MIN && Restart < restart ) {
483 printf(
"### WARNING: vGMRES restart number set to %d!\n", Restart);
490 r = work; w = r + n; rs = w + n; c = rs + Restart1; s = c + Restart;
492 for ( i = 0; i < Restart1; i++ ) p[i] = s + Restart + i*n;
494 for ( i = 0; i < Restart1; i++ ) hh[i] = p[Restart] + n + i*Restart;
506 relres = r_norm/absres0;
515 relres = r_normb/absres0;
520 relres = absres0/normu;
523 printf(
"### ERROR: Unknown stopping type! [%s]\n", __FUNCTION__);
528 if ( relres < tol || absres0 < 1e-12*tol )
goto FINISHED;
537 while ( iter < MaxIt ) {
539 rs[0] = r_norm_old = r_norm;
548 if ( cr > cr_max || iter == 0 ) {
549 Restart = restart_max;
551 else if ( cr < cr_min ) {
555 if ( Restart - d > restart_min ) {
559 Restart = restart_max;
565 while ( i < Restart && iter < MaxIt ) {
578 for (j = 0; j < i; j ++) {
589 for (j = 1; j < i; ++j) {
591 hh[j-1][i-1] = s[j-1]*hh[j][i-1] + c[j-1]*t;
592 hh[j][i-1] = -s[j-1]*t + c[j-1]*hh[j][i-1];
594 t= hh[i][i-1]*hh[i][i-1];
595 t+= hh[i-1][i-1]*hh[i-1][i-1];
598 if (gamma == 0.0) gamma = epsmac;
599 c[i-1] = hh[i-1][i-1] / gamma;
600 s[i-1] = hh[i][i-1] / gamma;
601 rs[i] = -s[i-1]*rs[i-1];
602 rs[i-1] = c[i-1]*rs[i-1];
603 hh[i-1][i-1] = s[i-1]*hh[i][i-1] + c[i-1]*hh[i-1][i-1];
605 absres = r_norm = fabs(rs[i]);
607 relres = absres/absres0;
609 norms[iter] = relres;
612 fasp_itinfo(PrtLvl, StopType, iter, relres, absres,
613 norms[iter]/norms[iter-1]);
616 if ( relres < tol && iter >= MIN_ITER )
break;
621 rs[i-1] = rs[i-1] / hh[i-1][i-1];
622 for (k = i-2; k >= 0; k --) {
624 for (j = k+1; j < i; j ++) t -= hh[k][j]*rs[j];
627 rs[k] = t / hh[k][k];
645 if ( relres < tol && iter >= MIN_ITER ) {
647 REAL computed_relres = relres;
655 switch ( StopType ) {
658 relres = absres/absres0;
666 relres = absres/absres0;
671 relres = absres/normu;
675 norms[iter] = relres;
677 if ( relres < tol ) {
686 ITS_COMPRES(computed_relres); ITS_REALRES(relres);
692 for ( j = i; j > 0; j-- ) {
693 rs[j-1] = -s[j-1]*rs[j];
694 rs[j] = c[j-1]*rs[j];
709 cr = r_norm / r_norm_old;
714 if ( PrtLvl >
PRINT_NONE ) ITS_FINAL(iter,MaxIt,relres);
725 printf(
"### DEBUG: [--End--] %s ...\n", __FUNCTION__);
764 const SHORT StopType,
768 const INT MIN_ITER = 0;
775 const REAL cr_max = 0.99;
776 const REAL cr_min = 0.174;
782 REAL r_norm, r_normb, gamma, t;
787 REAL r_norm_old = 0.0;
789 INT restart_max = restart;
792 INT Restart = restart;
793 INT Restart1 = Restart + 1;
794 unsigned LONG worksize = (Restart+4)*(Restart+n)+1-n;
797 REAL *c = NULL, *s = NULL, *rs = NULL;
798 REAL *norms = NULL, *r = NULL, *w = NULL;
800 REAL **p = NULL, **hh = NULL;
803 if ( PrtLvl >
PRINT_NONE ) printf(
"\nCalling VGMRes solver (BLC) ...\n");
806 printf(
"### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
807 printf(
"### DEBUG: maxit = %d, tol = %.4le\n", MaxIt, tol);
814 while ( (work == NULL) && (Restart > 5) ) {
815 Restart = Restart - 5;
816 worksize = (Restart+4)*(Restart+n)+1-n;
818 Restart1 = Restart + 1;
821 if ( work == NULL ) {
822 printf(
"### ERROR: No enough memory! [%s:%d]\n", __FILE__, __LINE__ );
826 if ( PrtLvl >
PRINT_MIN && Restart < restart ) {
827 printf(
"### WARNING: vGMRES restart number set to %d!\n", Restart);
834 r = work; w = r + n; rs = w + n; c = rs + Restart1; s = c + Restart;
836 for ( i = 0; i < Restart1; i++ ) p[i] = s + Restart + i*n;
838 for ( i = 0; i < Restart1; i++ ) hh[i] = p[Restart] + n + i*Restart;
850 relres = r_norm/absres0;
859 relres = r_normb/absres0;
864 relres = absres0/normu;
867 printf(
"### ERROR: Unknown stopping type! [%s]\n", __FUNCTION__);
872 if ( relres < tol || absres0 < 1e-12*tol )
goto FINISHED;
881 while ( iter < MaxIt ) {
883 rs[0] = r_norm_old = r_norm;
892 if ( cr > cr_max || iter == 0 ) {
893 Restart = restart_max;
895 else if ( cr < cr_min ) {
899 if ( Restart - d > restart_min ) {
903 Restart = restart_max;
909 while ( i < Restart && iter < MaxIt ) {
922 for (j = 0; j < i; j ++) {
933 for (j = 1; j < i; ++j) {
935 hh[j-1][i-1] = s[j-1]*hh[j][i-1] + c[j-1]*t;
936 hh[j][i-1] = -s[j-1]*t + c[j-1]*hh[j][i-1];
938 t= hh[i][i-1]*hh[i][i-1];
939 t+= hh[i-1][i-1]*hh[i-1][i-1];
942 if (gamma == 0.0) gamma = epsmac;
943 c[i-1] = hh[i-1][i-1] / gamma;
944 s[i-1] = hh[i][i-1] / gamma;
945 rs[i] = -s[i-1]*rs[i-1];
946 rs[i-1] = c[i-1]*rs[i-1];
947 hh[i-1][i-1] = s[i-1]*hh[i][i-1] + c[i-1]*hh[i-1][i-1];
949 absres = r_norm = fabs(rs[i]);
951 relres = absres/absres0;
953 norms[iter] = relres;
956 fasp_itinfo(PrtLvl, StopType, iter, relres, absres,
957 norms[iter]/norms[iter-1]);
960 if ( relres < tol && iter >= MIN_ITER )
break;
965 rs[i-1] = rs[i-1] / hh[i-1][i-1];
966 for (k = i-2; k >= 0; k --) {
968 for (j = k+1; j < i; j ++) t -= hh[k][j]*rs[j];
971 rs[k] = t / hh[k][k];
989 if ( relres < tol && iter >= MIN_ITER ) {
991 REAL computed_relres = relres;
999 switch ( StopType ) {
1002 relres = absres/absres0;
1010 relres = absres/absres0;
1015 relres = absres/normu;
1019 norms[iter] = relres;
1021 if ( relres < tol ) {
1030 ITS_COMPRES(computed_relres); ITS_REALRES(relres);
1036 for ( j = i; j > 0; j-- ) {
1037 rs[j-1] = -s[j-1]*rs[j];
1038 rs[j] = c[j-1]*rs[j];
1053 cr = r_norm / r_norm_old;
1058 if ( PrtLvl >
PRINT_NONE ) ITS_FINAL(iter,MaxIt,relres);
1069 printf(
"### DEBUG: [--End--] %s ...\n", __FUNCTION__);
1110 const SHORT restart,
1111 const SHORT StopType,
1115 const INT MIN_ITER = 0;
1122 const REAL cr_max = 0.99;
1123 const REAL cr_min = 0.174;
1129 REAL r_norm, r_normb, gamma, t;
1134 REAL r_norm_old = 0.0;
1136 INT restart_max = restart;
1137 INT restart_min = 3;
1139 INT Restart = restart;
1140 INT Restart1 = Restart + 1;
1141 unsigned LONG worksize = (Restart+4)*(Restart+n)+1-n;
1144 REAL *c = NULL, *s = NULL, *rs = NULL;
1145 REAL *norms = NULL, *r = NULL, *w = NULL;
1147 REAL **p = NULL, **hh = NULL;
1150 if ( PrtLvl >
PRINT_NONE ) printf(
"\nCalling VGMRes solver (STR) ...\n");
1153 printf(
"### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
1154 printf(
"### DEBUG: maxit = %d, tol = %.4le\n", MaxIt, tol);
1161 while ( (work == NULL) && (Restart > 5) ) {
1162 Restart = Restart - 5;
1163 worksize = (Restart+4)*(Restart+n)+1-n;
1165 Restart1 = Restart + 1;
1168 if ( work == NULL ) {
1169 printf(
"### ERROR: No enough memory! [%s:%d]\n", __FILE__, __LINE__ );
1173 if ( PrtLvl >
PRINT_MIN && Restart < restart ) {
1174 printf(
"### WARNING: vGMRES restart number set to %d!\n", Restart);
1181 r = work; w = r + n; rs = w + n; c = rs + Restart1; s = c + Restart;
1183 for ( i = 0; i < Restart1; i++ ) p[i] = s + Restart + i*n;
1185 for ( i = 0; i < Restart1; i++ ) hh[i] = p[Restart] + n + i*Restart;
1197 relres = r_norm/absres0;
1206 relres = r_normb/absres0;
1211 relres = absres0/normu;
1214 printf(
"### ERROR: Unknown stopping type! [%s]\n", __FUNCTION__);
1219 if ( relres < tol || absres0 < 1e-12*tol )
goto FINISHED;
1228 while ( iter < MaxIt ) {
1230 rs[0] = r_norm_old = r_norm;
1239 if ( cr > cr_max || iter == 0 ) {
1240 Restart = restart_max;
1242 else if ( cr < cr_min ) {
1246 if ( Restart - d > restart_min ) {
1250 Restart = restart_max;
1256 while ( i < Restart && iter < MaxIt ) {
1269 for (j = 0; j < i; j ++) {
1280 for (j = 1; j < i; ++j) {
1282 hh[j-1][i-1] = s[j-1]*hh[j][i-1] + c[j-1]*t;
1283 hh[j][i-1] = -s[j-1]*t + c[j-1]*hh[j][i-1];
1285 t= hh[i][i-1]*hh[i][i-1];
1286 t+= hh[i-1][i-1]*hh[i-1][i-1];
1289 if (gamma == 0.0) gamma = epsmac;
1290 c[i-1] = hh[i-1][i-1] / gamma;
1291 s[i-1] = hh[i][i-1] / gamma;
1292 rs[i] = -s[i-1]*rs[i-1];
1293 rs[i-1] = c[i-1]*rs[i-1];
1294 hh[i-1][i-1] = s[i-1]*hh[i][i-1] + c[i-1]*hh[i-1][i-1];
1296 absres = r_norm = fabs(rs[i]);
1298 relres = absres/absres0;
1300 norms[iter] = relres;
1303 fasp_itinfo(PrtLvl, StopType, iter, relres, absres,
1304 norms[iter]/norms[iter-1]);
1307 if ( relres < tol && iter >= MIN_ITER )
break;
1312 rs[i-1] = rs[i-1] / hh[i-1][i-1];
1313 for (k = i-2; k >= 0; k --) {
1315 for (j = k+1; j < i; j ++) t -= hh[k][j]*rs[j];
1318 rs[k] = t / hh[k][k];
1336 if ( relres < tol && iter >= MIN_ITER ) {
1338 REAL computed_relres = relres;
1346 switch ( StopType ) {
1349 relres = absres/absres0;
1357 relres = absres/absres0;
1362 relres = absres/normu;
1366 norms[iter] = relres;
1368 if ( relres < tol ) {
1377 ITS_COMPRES(computed_relres); ITS_REALRES(relres);
1383 for ( j = i; j > 0; j-- ) {
1384 rs[j-1] = -s[j-1]*rs[j];
1385 rs[j] = c[j-1]*rs[j];
1400 cr = r_norm / r_norm_old;
1405 if ( PrtLvl >
PRINT_NONE ) ITS_FINAL(iter,MaxIt,relres);
1416 printf(
"### DEBUG: [--End--] %s ...\n", __FUNCTION__);
1458 const SHORT StopType,
1462 const INT min_iter = 0;
1468 const REAL cr_max = 0.99;
1469 const REAL cr_min = 0.174;
1476 REAL r_norm, b_norm, den_norm;
1477 REAL epsilon, gamma, t;
1479 REAL *c = NULL, *s = NULL, *rs = NULL;
1480 REAL *norms = NULL, *r = NULL, *w = NULL;
1481 REAL **p = NULL, **hh = NULL;
1485 REAL r_norm_old = 0.0;
1487 INT restart_max = restart;
1488 INT restart_min = 3;
1490 INT Restart = restart;
1491 INT Restart1 = Restart + 1;
1492 unsigned LONG worksize = (restart+4)*(restart+n)+1-n;
1495 if ( PrtLvl >
PRINT_NONE ) printf(
"\nCalling VGMRes solver (MatFree) ...\n");
1498 printf(
"### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
1499 printf(
"### DEBUG: maxit = %d, tol = %.4le\n", MaxIt, tol);
1506 while ( (work == NULL) && (Restart > 5) ) {
1507 Restart = Restart - 5;
1508 worksize = (Restart+4)*(Restart+n)+1-n;
1510 Restart1 = Restart + 1;
1513 if ( work == NULL ) {
1514 printf(
"### ERROR: No enough memory! [%s:%d]\n", __FILE__, __LINE__ );
1518 if ( PrtLvl >
PRINT_MIN && Restart < restart ) {
1519 printf(
"### WARNING: vGMRES restart number set to %d!\n", Restart);
1526 r = work; w = r + n; rs = w + n; c = rs + Restart1; s = c + Restart;
1527 for (i = 0; i < Restart1; i ++) p[i] = s + Restart + i*n;
1528 for (i = 0; i < Restart1; i ++) hh[i] = p[Restart] + n + i*Restart;
1539 ITS_PUTNORM(
"right-hand side", b_norm);
1540 ITS_PUTNORM(
"residual", r_norm);
1543 if (b_norm > 0.0) den_norm = b_norm;
1544 else den_norm = r_norm;
1546 epsilon = tol*den_norm;
1549 while (iter < MaxIt) {
1551 r_norm_old = r_norm;
1552 if (r_norm == 0.0) {
1564 if (cr > cr_max || iter == 0) {
1565 Restart = restart_max;
1567 else if (cr < cr_min) {
1571 if (Restart - d > restart_min) {
1575 Restart = restart_max;
1579 if (r_norm <= epsilon && iter >= min_iter) {
1584 if (r_norm <= epsilon) {
1599 while (i < Restart && iter < MaxIt) {
1612 for (j = 0; j < i; j ++) {
1624 for (j = 1; j < i; ++j) {
1626 hh[j-1][i-1] = s[j-1]*hh[j][i-1] + c[j-1]*t;
1627 hh[j][i-1] = -s[j-1]*t + c[j-1]*hh[j][i-1];
1629 t= hh[i][i-1]*hh[i][i-1];
1630 t+= hh[i-1][i-1]*hh[i-1][i-1];
1632 if (gamma == 0.0) gamma = epsmac;
1633 c[i-1] = hh[i-1][i-1] / gamma;
1634 s[i-1] = hh[i][i-1] / gamma;
1635 rs[i] = -s[i-1]*rs[i-1];
1636 rs[i-1] = c[i-1]*rs[i-1];
1637 hh[i-1][i-1] = s[i-1]*hh[i][i-1] + c[i-1]*hh[i-1][i-1];
1638 r_norm = fabs(rs[i]);
1640 norms[iter] = r_norm;
1643 fasp_itinfo(PrtLvl,StopType,iter,norms[iter]/b_norm,
1644 norms[iter],norms[iter]/norms[iter-1]);
1647 fasp_itinfo(PrtLvl,StopType,iter,norms[iter],norms[iter],
1648 norms[iter]/norms[iter-1]);
1652 if (r_norm <= epsilon && iter >= min_iter)
break;
1658 rs[i-1] = rs[i-1] / hh[i-1][i-1];
1659 for (k = i-2; k >= 0; k --) {
1661 for (j = k+1; j < i; j ++) t -= hh[k][j]*rs[j];
1664 rs[k] = t / hh[k][k];
1679 if (r_norm <= epsilon && iter >= min_iter) {
1684 if (r_norm <= epsilon) {
1694 for (j = i; j > 0; j--) {
1695 rs[j-1] = -s[j-1]*rs[j];
1696 rs[j] = c[j-1]*rs[j];
1711 cr = r_norm / r_norm_old;
1715 if (PrtLvl >
PRINT_NONE) ITS_FINAL(iter,MaxIt,r_norm);
1726 printf(
"### DEBUG: [--End--] %s ...\n", __FUNCTION__);
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.
void fasp_chkerr(const SHORT status, const char *fctname)
Check error status and print out error messages before quit.
REAL fasp_blas_darray_dotprod(const INT n, const REAL *x, const REAL *y)
Inner product of two arraies x and y.
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_ax(const INT n, const REAL a, REAL *x)
x = a*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_dbsr_aAxpy(const REAL alpha, const dBSRmat *A, const REAL *x, REAL *y)
Compute y := alpha*A*x + y.
void fasp_blas_dbsr_mxv(const dBSRmat *A, const REAL *x, REAL *y)
Compute y := A*x.
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_pvgmres(dSTRmat *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...
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...
INT fasp_solver_dblc_pvgmres(dBLCmat *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...
INT fasp_solver_pvgmres(mxv_matfree *mf, 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...
INT fasp_solver_dbsr_pvgmres(dBSRmat *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...
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_MAXIT
#define STOP_REL_RES
Definition of iterative solver stopping criteria types.
#define BIGREAL
Some global constants.
#define PRINT_NONE
Print level for all subroutines – not including DEBUG output.
Block REAL CSR matrix format.
Block sparse row storage matrix of REAL type.
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