Solve the linear system of equations A * x = b
by means of the Bi-Conjugate Gradient iterative method.
The input arguments are:
Afun
such that Afun (x, "notransp") = A * x
and Afun (x, "transp") = A' * x
. Additional parameters to Afun
may be passed after x0. b - A * x
. The iteration stops if norm (b - A * x)
≤ tol * norm (b)
. If tol is omitted or empty, then a tolerance of 1e-6 is used. M = M1 * M2
. Both M1 and M2 can be passed as a matrix or as a function handle or inline function g
such that g (x, "notransp") = M1 \ x
or g (x, "notransp") = M2 \ x
and g (x, "transp") = M1' \ x
or g (x, "transp") = M2' \ x
. If M1 is omitted or empty, then preconditioning is not applied. The preconditioned system is theoretically equivalent to applying the bicg
method to the linear system inv (M1) * A * inv (M2) * y = inv
(M1) * b
and inv (M2') * A' * inv (M1') * z =
inv (M2') * b
and then setting x = inv (M2) * y
. Any arguments which follow x0 are treated as parameters, and passed in an appropriate manner to any of the functions (Afun or Mfun) or that have been given to bicg
.
The output parameters are:
A * x = b
. If the algorithm did not converge, then x is the iteration which has the minimum residual. eps * norm (x,2)
. length (resvec) - 1
. Consider a trivial problem with a tridiagonal matrix
n = 20; A = toeplitz (sparse ([1, 1], [1, 2], [2, 1] * n ^ 2, 1, n)) + ... toeplitz (sparse (1, 2, -1, 1, n) * n / 2, ... sparse (1, 2, 1, 1, n) * n / 2); b = A * ones (n, 1); restart = 5; [M1, M2] = ilu (A); # in this tridiag case, it corresponds to lu (A) M = M1 * M2; Afun = @(x, string) strcmp (string, "notransp") * (A * x) + ... strcmp (string, "transp") * (A' * x); Mfun = @(x, string) strcmp (string, "notransp") * (M \ x) + ... strcmp (string, "transp") * (M' \ x); M1fun = @(x, string) strcmp (string, "notransp") * (M1 \ x) + ... strcmp (string, "transp") * (M1' \ x); M2fun = @(x, string) strcmp (string, "notransp") * (M2 \ x) + ... strcmp (string, "transp") * (M2' \ x);
EXAMPLE 1: simplest usage of bicg
x = bicg (A, b)
EXAMPLE 2: bicg
with a function that computes A*x
and A'*x
x = bicg (Afun, b, [], n)
EXAMPLE 3: bicg
with a preconditioner matrix M
x = bicg (A, b, 1e-6, n, M)
EXAMPLE 4: bicg
with a function as preconditioner
x = bicg (Afun, b, 1e-6, n, Mfun)
EXAMPLE 5: bicg
with preconditioner matrices M1 and M2
x = bicg (A, b, 1e-6, n, M1, M2)
EXAMPLE 6: bicg
with functions as preconditioners
x = bicg (Afun, b, 1e-6, n, M1fun, M2fun)
EXAMPLE 7: bicg
with as input a function requiring an argument
function y = Ap (A, x, string, z) ## compute A^z * x or (A^z)' * x y = x; if (strcmp (string, "notransp")) for i = 1:z y = A * y; endfor elseif (strcmp (string, "transp")) for i = 1:z y = A' * y; endfor endif endfunction Apfun = @(x, string, p) Ap (A, x, string, p); x = bicg (Apfun, b, [], [], [], [], [], 2);
References:
Solve A x = b
using the stabilizied Bi-conjugate gradient iterative method.
The input parameters are:
Afun
such that Afun(x) = A * x
. Additional parameters to Afun
are passed after x0. b - A * x
. The iteration stops if norm (b - A * x)
≤ tol * norm (b)
. If tol is omitted or empty, then a tolerance of 1e-6 is used. min (20, numel (b))
is used. M = M1 * M2
. Both M1 and M2 can be passed as a matrix or as a function handle or inline function g
such that g(x) = M1 \ x
or g(x) = M2 \ x
. The technique used is the right preconditioning, i.e., it is solved A * inv (M) * y = b
and then x = inv (M) * y
. zeros (size (b))
is used. The arguments which follow x0 are treated as parameters, and passed in a proper way to any of the functions (A or M) which are passed to bicstab
.
The output parameters are:
(A*x-b) / norm(b)
. iter + 0.5
(length(resvec) - 1) / 2
is possible to see the total number of (total) iterations performed. Let us consider a trivial problem with a tridiagonal matrix
n = 20; A = toeplitz (sparse ([1, 1], [1, 2], [2, 1] * n ^ 2, 1, n)) + ... toeplitz (sparse (1, 2, -1, 1, n) * n / 2, ... sparse (1, 2, 1, 1, n) * n / 2); b = A * ones (n, 1); restart = 5; [M1, M2] = ilu (A); # in this tridiag case, it corresponds to lu (A) M = M1 * M2; Afun = @(x) A * x; Mfun = @(x) M \ x; M1fun = @(x) M1 \ x; M2fun = @(x) M2 \ x;
EXAMPLE 1: simplest usage of bicgstab
x = bicgstab (A, b, [], n)
EXAMPLE 2: bicgstab
with a function which computes A * x
x = bicgstab (Afun, b, [], n)
EXAMPLE 3: bicgstab
with a preconditioner matrix M
x = bicgstab (A, b, [], 1e-06, n, M)
EXAMPLE 4: bicgstab
with a function as preconditioner
x = bicgstab (Afun, b, 1e-6, n, Mfun)
EXAMPLE 5: bicgstab
with preconditioner matrices M1 and M2
x = bicgstab (A, b, [], 1e-6, n, M1, M2)
EXAMPLE 6: bicgstab
with functions as preconditioners
x = bicgstab (Afun, b, 1e-6, n, M1fun, M2fun)
EXAMPLE 7: bicgstab
with as input a function requiring an argument
function y = Ap (A, x, z) # compute A^z * x y = x; for i = 1:z y = A * y; endfor endfunction Apfun = @(x, string, p) Ap (A, x, string, p); x = bicgstab (Apfun, b, [], [], [], [], [], 2);
EXAMPLE 8: explicit example to show that bicgstab
uses a right preconditioner
[M1, M2] = ilu (A + 0.1 * eye (n)); # factorization of A perturbed M = M1 * M2; ## reference solution computed by bicgstab after one iteration [x_ref, fl] = bicgstab (A, b, [], 1, M) ## right preconditioning [y, fl] = bicgstab (A / M, b, [], 1) x = M \ y # compare x and x_ref
References:
Solve A x = b
, where A is a square matrix, using the Conjugate Gradients Squared method.
The input arguments are:
Afun
such that Afun(x) = A * x
. Additional parameters to Afun
are passed after x0. min (20, numel (b))
is used. M = M1 * M2
. Both M1 and M2 can be passed as a matrix or as a function handle or inline function g
such that g(x) = M1 \ x
or g(x) = M2 \ x
. If M1 is empty or not passed then no preconditioners are applied. The technique used is the right preconditioning, i.e., it is solved A*inv(M)*y = b
and then x = inv(M)*y
. zeros (size (b))
is used. The arguments which follow x0 are treated as parameters, and passed in a proper way to any of the functions (A or P) which are passed to cgs
.
The output parameters are:
(A*x-b) / norm(b)
. length(resvec) - 1
is possible to see the total number of iterations performed. Let us consider a trivial problem with a tridiagonal matrix
n = 20; A = toeplitz (sparse ([1, 1], [1, 2], [2, 1] * n ^ 2, 1, n)) + ... toeplitz (sparse (1, 2, -1, 1, n) * n / 2, ... sparse (1, 2, 1, 1, n) * n / 2); b = A * ones (n, 1); restart = 5; [M1, M2] = ilu (A); # in this tridiag case it corresponds to chol (A)' M = M1 * M2; Afun = @(x) A * x; Mfun = @(x) M \ x; M1fun = @(x) M1 \ x; M2fun = @(x) M2 \ x;
EXAMPLE 1: simplest usage of cgs
x = cgs (A, b, [], n)
EXAMPLE 2: cgs
with a function which computes A * x
x = cgs (Afun, b, [], n)
EXAMPLE 3: cgs
with a preconditioner matrix M
x = cgs (A, b, [], 1e-06, n, M)
EXAMPLE 4: cgs
with a function as preconditioner
x = cgs (Afun, b, 1e-6, n, Mfun)
EXAMPLE 5: cgs
with preconditioner matrices M1 and M2
x = cgs (A, b, [], 1e-6, n, M1, M2)
EXAMPLE 6: cgs
with functions as preconditioners
x = cgs (Afun, b, 1e-6, n, M1fun, M2fun)
EXAMPLE 7: cgs
with as input a function requiring an argument
function y = Ap (A, x, z) # compute A^z * x y = x; for i = 1:z y = A * y; endfor endfunction Apfun = @(x, string, p) Ap (A, x, string, p); x = cgs (Apfun, b, [], [], [], [], [], 2);
EXAMPLE 8: explicit example to show that cgs
uses a right preconditioner
[M1, M2] = ilu (A + 0.3 * eye (n)); # factorization of A perturbed M = M1 * M2; ## reference solution computed by cgs after one iteration [x_ref, fl] = cgs (A, b, [], 1, M) ## right preconditioning [y, fl] = cgs (A / M, b, [], 1) x = M \ y # compare x and x_ref
References:
Solve A x = b
using the Preconditioned GMRES iterative method with restart, a.k.a. PGMRES(restart).
The input arguments are:
Afun
such that Afun(x) = A * x
. Additional parameters to Afun
are passed after x0. inv (M) * (b - a * x)
. The iteration stops if norm (inv (M) * (b - a * x))
≤ tol * norm (inv (M) * B)
. If tol is omitted or empty, then a tolerance of 1e-6 is used. min (10, N / restart)
is used. Note that, if restart is empty, then maxit is the maximum number of iterations. If restart and maxit are not empty, then the maximum number of iterations is restart * maxit
. If both restart and maxit are empty, then the maximum number of iterations is set to min (10, N)
. M = M1 * M2
. Both M1 and M2 can be passed as a matrix, function handle, or inline function g
such that g(x) = M1 \ x
or g(x) = M2 \ x
. If M1 is [] or not given, then the preconditioner is not applied. The technique used is the left-preconditioning, i.e., it is solved inv(M) * A * x = inv(M) * b
instead of A * x = b
. zeros (size (b))
is used. The arguments which follow x0 are treated as parameters, and passed in a proper way to any of the functions (A or M or M1 or M2) which are passed to gmres
.
The outputs are:
consecutive iterations is less than eps)
iter(2) = restart
. If restart is empty or N, then 1 ≤ iter(2) ≤ maxit. To be more clear, the approximation x is computed at the iteration (iter(1) - 1) * restart + iter(2)
. Since the output x corresponds to the minimal preconditioned residual solution, the total number of iterations that the method performed is given by length (resvec) - 1
.
norm (A * x0 - b)
. Let us consider a trivial problem with a tridiagonal matrix
n = 20; A = toeplitz (sparse ([1, 1], [1, 2], [2, 1] * n ^ 2, 1, n)) + ... toeplitz (sparse (1, 2, -1, 1, n) * n / 2, ... sparse (1, 2, 1, 1, n) * n / 2); b = A * ones (n, 1); restart = 5; [M1, M2] = ilu (A); # in this tridiag case, it corresponds to lu (A) M = M1 * M2; Afun = @(x) A * x; Mfun = @(x) M \ x; M1fun = @(x) M1 \ x; M2fun = @(x) M2 \ x;
EXAMPLE 1: simplest usage of gmres
x = gmres (A, b, [], [], n)
EXAMPLE 2: gmres
with a function which computes A * x
x = gmres (Afun, b, [], [], n)
EXAMPLE 3: usage of gmres
with the restart
x = gmres (A, b, restart);
EXAMPLE 4: gmres
with a preconditioner matrix M with and without restart
x = gmres (A, b, [], 1e-06, n, M) x = gmres (A, b, restart, 1e-06, n, M)
EXAMPLE 5: gmres
with a function as preconditioner
x = gmres (Afun, b, [], 1e-6, n, Mfun)
EXAMPLE 6: gmres
with preconditioner matrices M1 and M2
x = gmres (A, b, [], 1e-6, n, M1, M2)
EXAMPLE 7: gmres
with functions as preconditioners
x = gmres (Afun, b, 1e-6, n, M1fun, M2fun)
EXAMPLE 8: gmres
with as input a function requiring an argument
function y = Ap (A, x, p) # compute A^p * x y = x; for i = 1:p y = A * y; endfor endfunction Apfun = @(x, p) Ap (A, x, p); x = gmres (Apfun, b, [], [], [], [], [], [], 2);
EXAMPLE 9: explicit example to show that gmres
uses a left preconditioner
[M1, M2] = ilu (A + 0.1 * eye (n)); # factorization of A perturbed M = M1 * M2; ## reference solution computed by gmres after two iterations [x_ref, fl] = gmres (A, b, [], [], 1, M) ## left preconditioning [x, fl] = gmres (M \ A, M \ b, [], [], 1) x # compare x and x_ref
References:
Solve A x = b
using the Quasi-Minimal Residual iterative method (without look-ahead).
min (20, numel (b))
is used. zeros (size (b))
is used. A can be passed as a matrix or as a function handle or inline function f
such that f(x, "notransp") = A*x
and f(x, "transp") = A'*x
.
The preconditioner P is given as P = M1 * M2
. Both M1 and M2 can be passed as a matrix or as a function handle or inline function g
such that g(x, "notransp") = M1 \ x
or g(x, "notransp") = M2 \ x
and g(x, "transp") = M1' \ x
or g(x, "transp") = M2' \ x
.
If called with more than one output parameter
(the value 2 is unused but skipped for compatibility).
References:
Solve A x = b
using the Transpose-Tree qmr method, based on the cgs.
The input parameters are:
Afun
such that Afun(x) = A * x
. Additional parameters to Afun
are passed after x0. min (20, numel (b))
is used. To be compatible, since the method as different behaviors in the iteration number is odd or even, is considered as iteration in tfqmr
the entire odd-even cycle. That is, to make an entire iteration, the algorithm performs two sub-iterations: the odd one and the even one. M = M1 * M2
. Both M1 and M2 can be passed as a matrix or as a function handle or inline function g
such that g(x) = M1 \ x
or g(x) = M2 \ x
. The technique used is the right-preconditioning, i.e., it is solved A*inv(M)*y = b
and then x = inv(M)*y
, instead of A x = b
. zeros (size (b))
is used. The arguments which follow x0 are treated as parameters, and passed in a proper way to any of the functions (A or M) which are passed to tfqmr
.
The output parameters are:
(A*x-b) / norm (b)
. norm (b - A x0)
). Doing length (resvec) - 1
is possible to see the total number of iterations performed. Let us consider a trivial problem with a tridiagonal matrix
n = 20; A = toeplitz (sparse ([1, 1], [1, 2], [2, 1] * n ^ 2, 1, n)) + ... toeplitz (sparse (1, 2, -1, 1, n) * n / 2, ... sparse (1, 2, 1, 1, n) * n / 2); b = A * ones (n, 1); restart = 5; [M1, M2] = ilu (A); # in this tridiag case it corresponds to chol (A)' M = M1 * M2; Afun = @(x) A * x; Mfun = @(x) M \ x; M1fun = @(x) M1 \ x; M2fun = @(x) M2 \ x;
EXAMPLE 1: simplest usage of tfqmr
x = tfqmr (A, b, [], n)
EXAMPLE 2: tfqmr
with a function which computes A * x
x = tfqmr (Afun, b, [], n)
EXAMPLE 3: tfqmr
with a preconditioner matrix M
x = tfqmr (A, b, [], 1e-06, n, M)
EXAMPLE 4: tfqmr
with a function as preconditioner
x = tfqmr (Afun, b, 1e-6, n, Mfun)
EXAMPLE 5: tfqmr
with preconditioner matrices M1 and M2
x = tfqmr (A, b, [], 1e-6, n, M1, M2)
EXAMPLE 6: tfmqr
with functions as preconditioners
x = tfqmr (Afun, b, 1e-6, n, M1fun, M2fun)
EXAMPLE 7: tfqmr
with as input a function requiring an argument
function y = Ap (A, x, z) # compute A^z * x y = x; for i = 1:z y = A * y; endfor endfunction Apfun = @(x, string, p) Ap (A, x, string, p); x = tfqmr (Apfun, b, [], [], [], [], [], 2);
EXAMPLE 8: explicit example to show that tfqmr
uses a right preconditioner
[M1, M2] = ilu (A + 0.3 * eye (n)); # factorization of A perturbed M = M1 * M2; ## reference solution computed by tfqmr after one iteration [x_ref, fl] = tfqmr (A, b, [], 1, M) ## right preconditioning [y, fl] = tfqmr (A / M, b, [], 1) x = M \ y # compare x and x_ref
References:
© 1996–2018 John W. Eaton
Permission is granted to make and distribute verbatim copies of this manual provided the copyright notice and this permission notice are preserved on all copies.
Permission is granted to copy and distribute modified versions of this manual under the conditions for verbatim copying, provided that the entire resulting derived work is distributed under the terms of a permission notice identical to this one.Permission is granted to copy and distribute translations of this manual into another language, under the above conditions for modified versions.
https://octave.org/doc/interpreter/Specialized-Solvers.html