DITSOL_PTFQMR (matvec, pcondl, pcondr, mstop, a, ia, x, b, n, ql, iql , qr, iqr, iparam, rparam, iwork, rwork, ierror)
DITSOL_PTFQMR has the standard parameter list for an iterative solver.
The quasi-minimal residual (QMR) method [Freund and Nachtigal 1991] is one of the algorithms proposed as a remedy for the irregular convergence behavior of the bi-conjugate gradient and the conjugate graient squared algorithms. Since these algorithms are not characterized by a minimization property, the residual norm often oscillates wildly. The QMR algorithm, by generating iterates that are defined by a quasi-minimization of the residual norm, results in smooth convergence curves. CXML includes TFQMR, the transpose-free variant of the QMR method, implemented without look-ahead [Freund 1993]. The implementation of the transpose-free quasi-minimal residual method requires the routine MATVEC to provide operations for job= 0. The routines MATVEC, PCONDL (if used), PCONDR (if used) and MSTOP (if used) should be declared external in your calling (sub)program. An upper bound for the two norm of the residual of the system being solved, is obtained during the implementation of the TFQMR method at no extra cost. This is the residual of the system to which the method is applied, which in the left and split preconditioned case is the preconditioned residual, inverse(QL) * r. To obtain the true residual, a non-negligible amount of extra computation would be required. Hence, for this method, only stopping criteria 3 and 4 are allowed. In the unpreconditioned case, the stopping criteria default to 1 and 2, respectively. Thus only istop = 3 and istop = 4 are permitted for both the preconditioned and unpreconditioned case. Additionally, a user-defined MSTOP is allowed, but the vectors r and z, corresponding the the real and preconditioned residuals, respectively, and passed as input parameters to the routine MSTOP, are undefined. CXML provides the following four forms of the method: • Unpreconditioned transpose-free quasi-minimal residual method: This is the transpose-free quasi-minimal residual method applied to A * x = b where A is a general matrix. As no preconditioning is used, both PCONDL and PCONDR are dummy input parameters. For the unpreconditioned transpose-free quasi-minimal residual method, the length of the real work space array, defined by the variable nrwk (IPARAM(4)), should be at least 7*n, where n is the order of the matrix A. The vectors r and z, passed as input arguments to the routine MSTOP, are not defined. • Transpose-free quasi-minimal residual method with left preconditioning: This is the transpose-free quasi-minimal residual method applied to (inverse(QL) * A )* x = (inverse(QL) * b) The routine PCONDL, with job= 0 should evaluate v = inverse(QL) * u The routine PCONDR is not used and is therefore a dummy input parameter. For the transpose-free quasi-minimal residual method, with left preconditioning, the length of the real work space array, defined by the variable nrwk (IPARAM(4)), should be at least 8*n, where n is the order of the matrix A. This does not include the memory requirements of the preconditioner. The vectors r and z, passed as input arguments to the routine MSTOP, are not defined. • Transpose-free quasi-minimal residual method with right preconditioning: This is the transpose-free quasi-minimal residual method applied to ( A * inverse(QR)) * y = b where y = QR * x The routine PCONDR, with job= 0 should evaluate v = inverse(QR) * u The routine PCONDL is not used and is therefore a dummy input parameter. For the transpose-free quasi-minimal residual method, with right preconditioning, the length of the real work space array, defined by the variable nrwk (IPARAM(4)), should be at least 8*n, where n is the order of the matrix A. This does not include the memory requirements of the preconditioner. The vectors r and z, passed as input arguments to the routine MSTOP, are not defined. • Transpose-free quasi-minimal residual method with split preconditioning: This is the transpose-free quasi-minimal residual method applied to (inverse(QL) * A * inverse(QR)) * y = (inverse(QL) * b) where y = QR * x The routine PCONDL, with job= 0 should evaluate v = inverse(QL) * u and the routine PCONDR, with job= 0 should evaluate v = inverse(QR) * u For the transpose-free quasi-minimal residual method, with split preconditioning, the length of the real work space array, defined by the variable nrwk (IPARAM(4)), should be at least 9*n, where n is the order of the matrix A. This does not include the memory requirements of the preconditioner. The vectors r and z, passed as input arguments to the routine MSTOP, are not defined. This routine is available in both serial and parallel versions. The routine names and parameter list are identical for both versions. For information about linking to the serial or to the parallel library, refer to the CXML Reference Manual.