Previous: Functions of a Matrix, Up: Linear Algebra [Contents][Index]
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);
Reference:
Y. Saad, Iterative Methods for Sparse Linear Systems, Second edition, 2003, SIAM.
Solve A x = b using the stabilized 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
Reference:
Y. Saad, Iterative Methods for Sparse Linear Systems, Second edition, 2003, SIAM
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:
Y. Saad, Iterative Methods for Sparse Linear Systems, Second edition, 2003, SIAM
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
Reference:
Y. Saad, Iterative Methods for Sparse Linear Systems, Second edition, 2003, SIAM
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
Reference:
Y. Saad, Iterative Methods for Sparse Linear Systems, Second edition, 2003, SIAM
Previous: Functions of a Matrix, Up: Linear Algebra [Contents][Index]