Mercurial > octave
annotate scripts/sparse/bicg.m @ 27918:b442ec6dda5c
use centralized file for copyright info for individual contributors
* COPYRIGHT.md: New file.
* In most other files, use "Copyright (C) YYYY-YYYY The Octave Project
Developers" instead of tracking individual names in separate source
files. The motivation is to reduce the effort required to update the
notices each year.
Until now, the Octave source files contained copyright notices that
list individual contributors. I adopted these file-scope copyright
notices because that is what everyone was doing 30 years ago in the
days before distributed version control systems. But now, with many
contributors and modern version control systems, having these
file-scope copyright notices causes trouble when we update copyright
years or refactor code.
Over time, the file-scope copyright notices may become outdated as new
contributions are made or code is moved from one file to
another. Sometimes people contribute significant patches but do not
add a line claiming copyright. Other times, people add a copyright
notice for their contribution but then a later refactoring moves part
or all of their contribution to another file and the notice is not
moved with the code. As a practical matter, moving such notices is
difficult -- determining what parts are due to a particular
contributor requires a time-consuming search through the project
history. Even managing the yearly update of copyright years is
problematic. We have some contributors who are no longer
living. Should we update the copyright dates for their contributions
when we release new versions? Probably not, but we do still want to
claim copyright for the project as a whole.
To minimize the difficulty of maintaining the copyright notices, I
would like to change Octave's sources to use what is described here:
https://softwarefreedom.org/resources/2012/ManagingCopyrightInformation.html
in the section "Maintaining centralized copyright notices":
The centralized notice approach consolidates all copyright
notices in a single location, usually a top-level file.
This file should contain all of the copyright notices
provided project contributors, unless the contribution was
clearly insignificant. It may also credit -- without a copyright
notice -- anyone who helped with the project but did not
contribute code or other copyrighted material.
This approach captures less information about contributions
within individual files, recognizing that the DVCS is better
equipped to record those details. As we mentioned before, it
does have one disadvantage as compared to the file-scope
approach: if a single file is separated from the distribution,
the recipient won't see the contributors' copyright notices.
But this can be easily remedied by including a single
copyright notice in each file's header, pointing to the
top-level file:
Copyright YYYY-YYYY The Octave Project Developers
See the COPYRIGHT file at the top-level directory
of this distribution or at https://octave.org/COPYRIGHT.html.
followed by the usual GPL copyright statement.
For more background, see the discussion here:
https://lists.gnu.org/archive/html/octave-maintainers/2020-01/msg00009.html
Most files in the following directories have been skipped intentinally
in this changeset:
doc
libgui/qterminal
liboctave/external
m4
author | John W. Eaton <jwe@octave.org> |
---|---|
date | Mon, 06 Jan 2020 15:38:17 -0500 |
parents | 4e632e942e84 |
children | 1891570abac8 |
rev | line source |
---|---|
27918
b442ec6dda5c
use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents:
27356
diff
changeset
|
1 ## Copyright (C) 2006-2019 The Octave Project Developers |
b442ec6dda5c
use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents:
27356
diff
changeset
|
2 ## |
b442ec6dda5c
use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents:
27356
diff
changeset
|
3 ## See the file COPYRIGHT.md in the top-level directory of this distribution |
b442ec6dda5c
use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents:
27356
diff
changeset
|
4 ## or <https://octave.org/COPYRIGHT.html/>. |
b442ec6dda5c
use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents:
27356
diff
changeset
|
5 ## |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
6 ## |
17795
0a8c35ae5ce1
maint: Fix various problems with GPL copyright statements.
Rik <rik@octave.org>
parents:
17744
diff
changeset
|
7 ## This file is part of Octave. |
0a8c35ae5ce1
maint: Fix various problems with GPL copyright statements.
Rik <rik@octave.org>
parents:
17744
diff
changeset
|
8 ## |
24534
194eb4bd202b
maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents:
23084
diff
changeset
|
9 ## Octave is free software: you can redistribute it and/or modify it |
22755
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22323
diff
changeset
|
10 ## under the terms of the GNU General Public License as published by |
24534
194eb4bd202b
maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents:
23084
diff
changeset
|
11 ## the Free Software Foundation, either version 3 of the License, or |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
12 ## (at your option) any later version. |
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
13 ## |
22755
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22323
diff
changeset
|
14 ## Octave is distributed in the hope that it will be useful, but |
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22323
diff
changeset
|
15 ## WITHOUT ANY WARRANTY; without even the implied warranty of |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
16 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
17 ## GNU General Public License for more details. |
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
18 ## |
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
19 ## You should have received a copy of the GNU General Public License |
22755
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22323
diff
changeset
|
20 ## along with Octave; see the file COPYING. If not, see |
24534
194eb4bd202b
maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents:
23084
diff
changeset
|
21 ## <https://www.gnu.org/licenses/>. |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
22 |
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
23 ## -*- texinfo -*- |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
24 ## @deftypefn {} {@var{x} =} bicg (@var{A}, @var{b}) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
25 ## @deftypefnx {} {@var{x} =} bicg (@var{A}, @var{b}, @var{tol}) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
26 ## @deftypefnx {} {@var{x} =} bicg (@var{A}, @var{b}, @var{tol}, @var{maxit}) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
27 ## @deftypefnx {} {@var{x} =} bicg (@var{A}, @var{b}, @var{tol}, @var{maxit}, @var{M}) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
28 ## @deftypefnx {} {@var{x} =} bicg (@var{A}, @var{b}, @var{tol}, @var{maxit}, @var{M1}, @var{M2}) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
29 ## @deftypefnx {} {@var{x} =} bicg (@var{A}, @var{b}, @var{tol}, @var{maxit}, @var{M}, [], @var{x0}) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
30 ## @deftypefnx {} {@var{x} =} bicg (@var{A}, @var{b}, @var{tol}, @var{maxit}, @var{M1}, @var{M2}, @var{x0}) |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
31 ## @deftypefnx {} {@var{x} =} bicg (@var{A}, @var{b}, @var{tol}, @var{maxit}, @var{M}, [], @var{x0}, @dots{}) |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
32 ## @deftypefnx {} {@var{x} =} bicg (@var{A}, @var{b}, @var{tol}, @var{maxit}, @var{M1}, @var{M2}, @var{x0}, @dots{}) |
20852
516bb87ea72e
2015 Code Sprint: remove class of function from docstring for all m-files.
Rik <rik@octave.org>
parents:
20711
diff
changeset
|
33 ## @deftypefnx {} {[@var{x}, @var{flag}, @var{relres}, @var{iter}, @var{resvec}] =} bicg (@var{A}, @var{b}, @dots{}) |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
34 ## Solve the linear system of equations @w{@code{@var{A} * @var{x} = @var{b}}} |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
35 ## by means of the Bi-Conjugate Gradient iterative method. |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
36 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
37 ## The input arguments are: |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
38 ## |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
39 ## @itemize |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
40 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
41 ## @item @var{A} is the matrix of the linear system and it must be square. |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
42 ## @var{A} can be passed as a matrix, function handle, or inline function |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
43 ## @code{Afun} such that @w{@code{Afun (x, "notransp") = A * x}} and |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
44 ## @w{@code{Afun (x, "transp") = A' * x}}. Additional parameters to |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
45 ## @code{Afun} may be passed after @var{x0}. |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
46 ## |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
47 ## @item @var{b} is the right-hand side vector. It must be a column vector |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
48 ## with the same number of rows as @var{A}. |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
49 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
50 ## @item |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
51 ## @var{tol} is the required relative tolerance for the residual error, |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
52 ## @w{@code{@var{b} - @var{A} * @var{x}}}. The iteration stops if |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
53 ## @w{@code{norm (@var{b} - @var{A} * @var{x})} @leq{} |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
54 ## @w{@code{@var{tol} * norm (@var{b})}}}. |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
55 ## If @var{tol} is omitted or empty, then a tolerance of 1e-6 is used. |
13929
9cae456085c2
Grammarcheck of documentation before 3.6.0 release.
Rik <octave@nomad.inbox5.com>
parents:
13796
diff
changeset
|
56 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
57 ## @item |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
58 ## @var{maxit} is the maximum allowed number of iterations; if @var{maxit} |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
59 ## is omitted or empty then a value of 20 is used. |
13929
9cae456085c2
Grammarcheck of documentation before 3.6.0 release.
Rik <octave@nomad.inbox5.com>
parents:
13796
diff
changeset
|
60 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
61 ## @item |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
62 ## @var{M1}, @var{M2} are the preconditioners. The preconditioner @var{M} is |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
63 ## given as @code{@var{M} = @var{M1} * @var{M2}}. Both @var{M1} and @var{M2} |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
64 ## can be passed as a matrix or as a function handle or inline function |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
65 ## @code{g} such that @w{@code{g (@var{x}, "notransp") = @var{M1} \ @var{x}}} |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
66 ## or @w{@code{g (@var{x}, "notransp") = @var{M2} \ @var{x}}} and |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
67 ## @w{@code{g (@var{x}, "transp") = @var{M1}' \ @var{x}}} or |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
68 ## @w{@code{g (@var{x}, "transp") = @var{M2}' \ @var{x}}}. |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
69 ## If @var{M1} is omitted or empty, then preconditioning is not applied. |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
70 ## The preconditioned system is theoretically equivalent to applying the |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
71 ## @code{bicg} method to the linear system |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
72 ## @code{inv (@var{M1}) * A * inv (@var{M2}) * @var{y} = inv |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
73 ## (@var{M1}) * @var{b}} and |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
74 ## @code{inv (@var{M2'}) * A' * inv (@var{M1'}) * @var{z} = |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
75 ## inv (@var{M2'}) * @var{b}} and then setting |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
76 ## @code{@var{x} = inv (@var{M2}) * @var{y}}. |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
77 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
78 ## @item |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
79 ## @var{x0} is the initial guess. If @var{x0} is omitted or empty then the |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
80 ## function sets @var{x0} to a zero vector by default. |
13091
e5aaba072d2b
maint: style fixes in iterative linear solvers
Carlo de Falco <kingcrimson@tiscali.it>
parents:
13073
diff
changeset
|
81 ## @end itemize |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
82 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
83 ## Any arguments which follow @var{x0} are treated as parameters, and passed in |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
84 ## an appropriate manner to any of the functions (@var{Afun} or @var{Mfun}) or |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
85 ## that have been given to @code{bicg}. |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
86 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
87 ## The output parameters are: |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
88 ## |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
89 ## @itemize |
13091
e5aaba072d2b
maint: style fixes in iterative linear solvers
Carlo de Falco <kingcrimson@tiscali.it>
parents:
13073
diff
changeset
|
90 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
91 ## @item |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
92 ## @var{x} is the computed approximation to the solution of |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
93 ## @w{@code{@var{A} * @var{x} = @var{b}}}. If the algorithm did not converge, |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
94 ## then @var{x} is the iteration which has the minimum residual. |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
95 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
96 ## @item |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
97 ## @var{flag} indicates the exit status: |
14853
72b8b39e12be
doc: Periodic grammarcheck of documentation.
Rik <octave@nomad.inbox5.com>
parents:
14609
diff
changeset
|
98 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
99 ## @itemize |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
100 ## @item 0: The algorithm converged to within the prescribed tolerance. |
13929
9cae456085c2
Grammarcheck of documentation before 3.6.0 release.
Rik <octave@nomad.inbox5.com>
parents:
13796
diff
changeset
|
101 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
102 ## @item 1: The algorithm did not converge and it reached the maximum number of |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
103 ## iterations. |
13929
9cae456085c2
Grammarcheck of documentation before 3.6.0 release.
Rik <octave@nomad.inbox5.com>
parents:
13796
diff
changeset
|
104 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
105 ## @item 2: The preconditioner matrix is singular. |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
106 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
107 ## @item 3: The algorithm stagnated, i.e., the absolute value of the |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
108 ## difference between the current iteration @var{x} and the previous is less |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
109 ## than @code{eps * norm (@var{x},2)}. |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
110 ## |
25200
b3ee0179d7b0
bicg.m: Check for denominators much smaller than numerators bicg.m (bug #53589).
Marco Caliari <marco.caliari@univr.it>
parents:
25166
diff
changeset
|
111 ## @item 4: The algorithm could not continue because intermediate values |
b3ee0179d7b0
bicg.m: Check for denominators much smaller than numerators bicg.m (bug #53589).
Marco Caliari <marco.caliari@univr.it>
parents:
25166
diff
changeset
|
112 ## became too small or too large for reliable computation. |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
113 ## @end itemize |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
114 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
115 ## @item |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
116 ## @var{relres} is the ratio of the final residual to its initial value, |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
117 ## measured in the Euclidean norm. |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
118 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
119 ## @item |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
120 ## @var{iter} is the iteration which @var{x} is computed. |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
121 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
122 ## @item |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
123 ## @var{resvec} is a vector containing the residual at each iteration. |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
124 ## The total number of iterations performed is given by |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
125 ## @code{length (@var{resvec}) - 1}. |
13091
e5aaba072d2b
maint: style fixes in iterative linear solvers
Carlo de Falco <kingcrimson@tiscali.it>
parents:
13073
diff
changeset
|
126 ## @end itemize |
16816
12005245b645
doc: Periodic grammarcheck of documentation.
Rik <rik@octave.org>
parents:
16768
diff
changeset
|
127 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
128 ## Consider a trivial problem with a tridiagonal matrix |
13929
9cae456085c2
Grammarcheck of documentation before 3.6.0 release.
Rik <octave@nomad.inbox5.com>
parents:
13796
diff
changeset
|
129 ## |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
130 ## @example |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
131 ## @group |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
132 ## n = 20; |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
133 ## A = toeplitz (sparse ([1, 1], [1, 2], [2, 1] * n ^ 2, 1, n)) + ... |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
134 ## toeplitz (sparse (1, 2, -1, 1, n) * n / 2, ... |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
135 ## sparse (1, 2, 1, 1, n) * n / 2); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
136 ## b = A * ones (n, 1); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
137 ## restart = 5; |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
138 ## [M1, M2] = ilu (A); # in this tridiag case, it corresponds to lu (A) |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
139 ## M = M1 * M2; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
140 ## Afun = @@(x, string) strcmp (string, "notransp") * (A * x) + ... |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
141 ## strcmp (string, "transp") * (A' * x); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
142 ## Mfun = @@(x, string) strcmp (string, "notransp") * (M \ x) + ... |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
143 ## strcmp (string, "transp") * (M' \ x); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
144 ## M1fun = @@(x, string) strcmp (string, "notransp") * (M1 \ x) + ... |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
145 ## strcmp (string, "transp") * (M1' \ x); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
146 ## M2fun = @@(x, string) strcmp (string, "notransp") * (M2 \ x) + ... |
25579
07c2c42f457e
doc: Miscellaneous documentation fixes all over the manual (bug #54288).
Rik <rik@octave.org>
parents:
25200
diff
changeset
|
147 ## strcmp (string, "transp") * (M2' \ x); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
148 ## @end group |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
149 ## @end example |
13929
9cae456085c2
Grammarcheck of documentation before 3.6.0 release.
Rik <octave@nomad.inbox5.com>
parents:
13796
diff
changeset
|
150 ## |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
151 ## @sc{Example 1:} simplest usage of @code{bicg} |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
152 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
153 ## @example |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
154 ## x = bicg (A, b) |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
155 ## @end example |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
156 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
157 ## @sc{Example 2:} @code{bicg} with a function that computes |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
158 ## @code{@var{A}*@var{x}} and @code{@var{A'}*@var{x}} |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
159 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
160 ## @example |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
161 ## x = bicg (Afun, b, [], n) |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
162 ## @end example |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
163 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
164 ## @sc{Example 3:} @code{bicg} with a preconditioner matrix @var{M} |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
165 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
166 ## @example |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
167 ## x = bicg (A, b, 1e-6, n, M) |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
168 ## @end example |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
169 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
170 ## @sc{Example 4:} @code{bicg} with a function as preconditioner |
13929
9cae456085c2
Grammarcheck of documentation before 3.6.0 release.
Rik <octave@nomad.inbox5.com>
parents:
13796
diff
changeset
|
171 ## |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
172 ## @example |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
173 ## x = bicg (Afun, b, 1e-6, n, Mfun) |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
174 ## @end example |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
175 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
176 ## @sc{Example 5:} @code{bicg} with preconditioner matrices @var{M1} |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
177 ## and @var{M2} |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
178 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
179 ## @example |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
180 ## x = bicg (A, b, 1e-6, n, M1, M2) |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
181 ## @end example |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
182 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
183 ## @sc{Example 6:} @code{bicg} with functions as preconditioners |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
184 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
185 ## @example |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
186 ## x = bicg (Afun, b, 1e-6, n, M1fun, M2fun) |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
187 ## @end example |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
188 ## |
24992
3d5a39077708
avoid some warnings from old versions of makeinfo (bug #53479)
Guillaume Flandin
parents:
24985
diff
changeset
|
189 ## @sc{Example 7:} @code{bicg} with as input a function requiring an argument |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
190 ## |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
191 ## @example |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
192 ## @group |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
193 ## function y = Ap (A, x, string, z) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
194 ## ## compute A^z * x or (A^z)' * x |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
195 ## y = x; |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
196 ## if (strcmp (string, "notransp")) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
197 ## for i = 1:z |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
198 ## y = A * y; |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
199 ## endfor |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
200 ## elseif (strcmp (string, "transp")) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
201 ## for i = 1:z |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
202 ## y = A' * y; |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
203 ## endfor |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
204 ## endif |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
205 ## endfunction |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
206 ## |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
207 ## Apfun = @@(x, string, p) Ap (A, x, string, p); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
208 ## x = bicg (Apfun, b, [], [], [], [], [], 2); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
209 ## @end group |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
210 ## @end example |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
211 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
212 ## References: |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
213 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
214 ## @enumerate |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
215 ## |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
216 ## @item @nospell{Y. Saad}, @cite{Iterative Methods for Sparse Linear Systems}, |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
217 ## Second edition, 2003, SIAM. |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
218 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
219 ## @end enumerate |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
220 ## |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
221 ## @seealso{bicgstab, cgs, gmres, pcg, qmr, tfqmr} |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
222 ## @end deftypefn |
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
223 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
224 function [x_min, flag, relres, iter_min, resvec] = ... |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
225 bicg (A, b, tol = [], maxit = [], M1 = [], M2 = [], x0 = [], varargin) |
13141
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
226 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
227 [Afun, M1fun, M2fun] = __alltohandles__ (A, b, M1, M2, "bicg"); |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
228 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
229 [tol, maxit, x0] = __default__input__ ({1e-06, min(rows(b), 20), ... |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
230 zeros(rows (b),1)}, tol, maxit, x0); |
13141
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
231 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
232 if (columns (b) == 2) |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
233 c = b(:,2); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
234 b = b(:,1); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
235 else |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
236 c = b; |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
237 endif |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
238 norm_b = norm (b, 2); |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
239 |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
240 if (norm_b == 0) # the only (only iff det(A) == 0) solution is x = 0 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
241 if (nargout < 2) |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
242 printf ("The right hand side vector is all zero so bicg\n") |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
243 printf ("returned an all zero solution without iterating.\n") |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
244 endif |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
245 x_min = zeros (numel (b), 1); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
246 flag = 0; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
247 relres = 0; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
248 iter_min = 0; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
249 resvec = 0; |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
250 return; |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
251 endif |
13141
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
252 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
253 x = x_min = x_pr = x0; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
254 iter = iter_min = 0; |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
255 flag = 1; # Default flag is "maximum number of iterations reached" |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
256 resvec = zeros (maxit + 1, 1); |
27356
4e632e942e84
doc: Improve documentation of sparse functions.
Kai T. Ohlhus <k.ohlhus@gmail.com>
parents:
26376
diff
changeset
|
257 r0 = b - Afun (x, "notransp", varargin{:}); # Residual of the system |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
258 s0 = c - Afun (x, "transp", varargin{:}); # Res. of the "dual system" |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
259 resvec(1) = norm (r0, 2); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
260 |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
261 try |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
262 warning ("error", "Octave:singular-matrix", "local") |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
263 prec_r0 = M1fun (r0, "notransp", varargin{:}); # r0 preconditioned |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
264 prec_s0 = s0; |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
265 prec_r0 = M2fun (prec_r0, "notransp", varargin{:}); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
266 prec_s0 = M2fun (prec_s0, "transp", varargin{:}); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
267 prec_s0 = M1fun (prec_s0, "transp", varargin{:}); # s0 preconditioned |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
268 p = prec_r0; # Direction of the system |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
269 q = prec_s0; # Direction of the "dual system" |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
270 catch |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
271 flag = 2; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
272 end_try_catch |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
273 |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
274 while ((flag != 2) && (iter < maxit) && (resvec(iter+1) >= norm_b * tol)) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
275 v = Afun (p, "notransp", varargin{:}); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
276 prod_qv = q' * v; |
25200
b3ee0179d7b0
bicg.m: Check for denominators much smaller than numerators bicg.m (bug #53589).
Marco Caliari <marco.caliari@univr.it>
parents:
25166
diff
changeset
|
277 alpha = (s0' * prec_r0); |
b3ee0179d7b0
bicg.m: Check for denominators much smaller than numerators bicg.m (bug #53589).
Marco Caliari <marco.caliari@univr.it>
parents:
25166
diff
changeset
|
278 if (abs (prod_qv) <= eps * abs (alpha)) |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
279 flag = 4; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
280 break |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
281 endif |
25200
b3ee0179d7b0
bicg.m: Check for denominators much smaller than numerators bicg.m (bug #53589).
Marco Caliari <marco.caliari@univr.it>
parents:
25166
diff
changeset
|
282 alpha ./= prod_qv; |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
283 x += alpha * p; |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
284 prod_rs = (s0' * prec_r0); # Product between r0 and s0 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
285 r0 -= alpha * v; |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
286 s0 -= conj (alpha) * Afun (q, "transp", varargin{:}); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
287 prec_r0 = M1fun (r0, "notransp", varargin{:}); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
288 prec_s0 = s0; |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
289 prec_r0 = M2fun (prec_r0, "notransp", varargin{:}); |
25200
b3ee0179d7b0
bicg.m: Check for denominators much smaller than numerators bicg.m (bug #53589).
Marco Caliari <marco.caliari@univr.it>
parents:
25166
diff
changeset
|
290 beta = s0' * prec_r0; |
b3ee0179d7b0
bicg.m: Check for denominators much smaller than numerators bicg.m (bug #53589).
Marco Caliari <marco.caliari@univr.it>
parents:
25166
diff
changeset
|
291 if (abs (prod_rs) <= abs (beta)) |
b3ee0179d7b0
bicg.m: Check for denominators much smaller than numerators bicg.m (bug #53589).
Marco Caliari <marco.caliari@univr.it>
parents:
25166
diff
changeset
|
292 flag = 4; |
b3ee0179d7b0
bicg.m: Check for denominators much smaller than numerators bicg.m (bug #53589).
Marco Caliari <marco.caliari@univr.it>
parents:
25166
diff
changeset
|
293 break; |
b3ee0179d7b0
bicg.m: Check for denominators much smaller than numerators bicg.m (bug #53589).
Marco Caliari <marco.caliari@univr.it>
parents:
25166
diff
changeset
|
294 endif |
b3ee0179d7b0
bicg.m: Check for denominators much smaller than numerators bicg.m (bug #53589).
Marco Caliari <marco.caliari@univr.it>
parents:
25166
diff
changeset
|
295 beta ./= prod_rs; |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
296 prec_s0 = M2fun (prec_s0, "transp", varargin{:}); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
297 prec_s0 = M1fun (prec_s0, "transp", varargin{:}); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
298 iter += 1; |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
299 resvec(iter+1) = norm (r0); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
300 if (resvec(iter+1) <= resvec(iter_min+1)) |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
301 x_min = x; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
302 iter_min = iter; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
303 endif |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
304 if (norm (x - x_pr) <= norm (x) * eps) |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
305 flag = 3; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
306 break; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
307 endif |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
308 p = prec_r0 + beta*p; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
309 q = prec_s0 + conj (beta) * q; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
310 endwhile |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
311 resvec = resvec(1:iter+1,1); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
312 |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
313 if (flag == 2) |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
314 relres = 1; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
315 else |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
316 relres = resvec(iter_min+1) / norm_b; |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
317 endif |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
318 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
319 if ((flag == 1) && (relres <= tol)) |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
320 flag = 0; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
321 endif |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
322 |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
323 if (nargout < 2) |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
324 switch (flag) |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
325 case 0 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
326 printf ("bicg converged at iteration %i ", iter_min); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
327 printf ("to a solution with relative residual %e\n", relres); |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
328 case 1 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
329 printf ("bicg stopped at iteration %i ", iter); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
330 printf ("without converging to the desired tolerance %e\n", tol); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
331 printf ("because the maximum number of iterations was reached. "); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
332 printf ("The iterate returned (number %i) has ", iter_min); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
333 printf ("relative residual %e\n", relres); |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
334 case 2 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
335 printf ("bicg stopped at iteration %i ", iter); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
336 printf ("without converging to the desired tolerance %e\n", tol); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
337 printf ("because the preconditioner matrix is singular.\n"); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
338 printf ("The iterate returned (number %i) ", iter_min); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
339 printf ("has relative residual %e\n", relres); |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
340 case 3 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
341 printf ("bicg stopped at iteration %i ", iter); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
342 printf ("without converging to the desired tolerance %e\n", tol); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
343 printf ("because the method stagnated.\n"); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
344 printf ("The iterate returned (number %i) ", iter_min); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
345 printf ("has relative residual %e\n", relres); |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
346 case 4 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
347 printf ("bicg stopped at iteration %i ", iter); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
348 printf ("without converging to the desired tolerance %e\n", tol); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
349 printf ("because the method can't continue.\n"); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
350 printf ("The iterate returned (number %i) ", iter_min); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
351 printf ("has relative residual %e\n", relres); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
352 endswitch |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
353 endif |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
354 |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
355 endfunction |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
356 |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
357 |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
358 %!demo |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
359 %! ## simplest use case |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
360 %! n = 20; |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
361 %! A = toeplitz (sparse ([1, 1], [1, 2], [2, 1] * n ^ 2, 1, n)) + ... |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
362 %! toeplitz (sparse (1, 2, -1, 1, n) * n / 2, ... |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
363 %! sparse (1, 2, 1, 1, n) * n / 2); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
364 %! b = A * ones (n, 1); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
365 %! [M1, M2] = ilu (A + 0.1 * eye (n)); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
366 %! M = M1 * M2; |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
367 %! x = bicg (A, b, [], n); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
368 %! function y = Ap (A, x, string, z) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
369 %! ## compute A^z * x or (A^z)' * x |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
370 %! y = x; |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
371 %! if (strcmp (string, "notransp")) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
372 %! for i = 1:z |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
373 %! y = A * y; |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
374 %! endfor |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
375 %! elseif (strcmp (string, "transp")) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
376 %! for i = 1:z |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
377 %! y = A' * y; |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
378 %! endfor |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
379 %! endif |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
380 %! endfunction |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
381 %! |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
382 %! Afun = @(x, string) Ap (A, x, string, 1); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
383 %! x = bicg (Afun, b, [], n); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
384 %! x = bicg (A, b, 1e-6, n, M); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
385 %! x = bicg (A, b, 1e-6, n, M1, M2); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
386 %! function y = Mfun (M, x, string) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
387 %! if (strcmp (string, "notransp")) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
388 %! y = M \ x; |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
389 %! else |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
390 %! y = M' \ x; |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
391 %! endif |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
392 %! endfunction |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
393 %! |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
394 %! M1fun = @(x, string) Mfun (M, x, string); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
395 %! x = bicg (Afun, b, 1e-6, n, M1fun); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
396 %! M1fun = @(x, string) Mfun (M1, x, string); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
397 %! M2fun = @(x, string) Mfun (M2, x, string); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
398 %! x = bicg (Afun, b, 1e-6, n, M1fun, M2fun); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
399 %! Afun = @(x, string, p) Ap (A, x, string, p); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
400 %! ## Solution of A^2 * x = b |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
401 %! x = bicg (Afun, b, [], 2*n, [], [], [], 2); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
402 |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
403 %!test |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
404 %! ## Check that all type of inputs work |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
405 %! A = sparse (toeplitz ([2, 1, 0, 0, 0], [2, -1, 0, 0, 0])); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
406 %! b = A * ones (5, 1); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
407 %! M1 = diag (sqrt (diag (A))); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
408 %! M2 = M1; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
409 %! Afun = @(z, string) strcmp (string, "notransp") * (A * z) + ... |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
410 %! strcmp (string, "transp") * (A' * z); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
411 %! M1_fun = @(z, string) strcmp (string,"notransp") * (M1 \ z) + ... |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
412 %! strcmp (string, "transp") * (M1' \ z); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
413 %! M2_fun = @(z, string) strcmp (string, "notransp") * (M2 \ z) + ... |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
414 %! strcmp (string, "transp") * (M2' \ z); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
415 %! [x, flag] = bicg (A, b); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
416 %! assert (flag, 0); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
417 %! [x, flag] = bicg (A, b, [], [], M1, M2); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
418 %! assert (flag, 0); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
419 %! [x, flag] = bicg (A, b, [], [], M1_fun, M2_fun); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
420 %! assert (flag, 0); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
421 %! [x, flag] = bicg (A, b,[], [], M1_fun, M2); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
422 %! assert (flag, 0); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
423 %! [x, flag] = bicg (A, b,[], [], M1, M2_fun); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
424 %! assert (flag, 0); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
425 %! [x, flag] = bicg (Afun, b); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
426 %! assert (flag, 0); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
427 %! [x, flag] = bicg (Afun, b,[], [], M1, M2); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
428 %! assert (flag, 0); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
429 %! [x, flag] = bicg (Afun, b,[], [], M1_fun, M2); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
430 %! assert (flag, 0); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
431 %! [x, flag] = bicg (Afun, b,[], [], M1, M2_fun); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
432 %! assert (flag, 0); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
433 %! [x, flag] = bicg (Afun, b,[], [], M1_fun, M2_fun); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
434 %! assert (flag, 0); |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
435 |
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
436 %!test |
13141
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
437 %! n = 100; |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
438 %! A = spdiags ([-2*ones(n,1) 4*ones(n,1) -ones(n,1)], -1:1, n, n); |
13141
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
439 %! b = sum (A, 2); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
440 %! tol = 1e-8; |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
441 %! maxit = 15; |
13141
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
442 %! M1 = spdiags ([ones(n,1)/(-2) ones(n,1)],-1:0, n, n); |
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
443 %! M2 = spdiags ([4*ones(n,1) -ones(n,1)], 0:1, n, n); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
444 %! [x, flag, relres, iter, resvec] = bicg (A, b, tol, maxit, M1, M2); |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
445 %! assert (norm (b - A*x) / norm (b), 0, tol); |
13796
53674ceb9133
maint: fix function definition in test
John W. Eaton <jwe@octave.org>
parents:
13141
diff
changeset
|
446 |
53674ceb9133
maint: fix function definition in test
John W. Eaton <jwe@octave.org>
parents:
13141
diff
changeset
|
447 %!function y = afun (x, t, a) |
16933
e39f00a32dc7
maint: Use parentheses around condition for switch(),while(),if() statements.
Rik <rik@octave.org>
parents:
16816
diff
changeset
|
448 %! switch (t) |
17174
c3c1ebfaa7dc
maint: Use common indentation for switch statement.
Rik <rik@octave.org>
parents:
16933
diff
changeset
|
449 %! case "notransp" |
c3c1ebfaa7dc
maint: Use common indentation for switch statement.
Rik <rik@octave.org>
parents:
16933
diff
changeset
|
450 %! y = a * x; |
c3c1ebfaa7dc
maint: Use common indentation for switch statement.
Rik <rik@octave.org>
parents:
16933
diff
changeset
|
451 %! case "transp" |
c3c1ebfaa7dc
maint: Use common indentation for switch statement.
Rik <rik@octave.org>
parents:
16933
diff
changeset
|
452 %! y = a' * x; |
13796
53674ceb9133
maint: fix function definition in test
John W. Eaton <jwe@octave.org>
parents:
13141
diff
changeset
|
453 %! endswitch |
53674ceb9133
maint: fix function definition in test
John W. Eaton <jwe@octave.org>
parents:
13141
diff
changeset
|
454 %!endfunction |
53674ceb9133
maint: fix function definition in test
John W. Eaton <jwe@octave.org>
parents:
13141
diff
changeset
|
455 %! |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
456 %!test |
13141
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
457 %! n = 100; |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
458 %! A = spdiags ([-2*ones(n,1) 4*ones(n,1) -ones(n,1)], -1:1, n, n); |
13141
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
459 %! b = sum (A, 2); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
460 %! tol = 1e-8; |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
461 %! maxit = 15; |
13141
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
462 %! M1 = spdiags ([ones(n,1)/(-2) ones(n,1)],-1:0, n, n); |
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
463 %! M2 = spdiags ([4*ones(n,1) -ones(n,1)], 0:1, n, n); |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
464 %! |
13141
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
465 %! [x, flag, relres, iter, resvec] = bicg (@(x, t) afun (x, t, A), |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
466 %! b, tol, maxit, M1, M2); |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
467 %! assert (x, ones (size (b)), 1e-7); |
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
468 |
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
469 %!test |
13141
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
470 %! n = 100; |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
471 %! tol = 1e-8; |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
472 %! a = sprand (n, n, .1); |
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
473 %! A = a' * a + 100 * eye (n); |
13141
e81ddf9cacd5
maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
13091
diff
changeset
|
474 %! b = sum (A, 2); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
475 %! [x, flag, relres, iter, resvec] = bicg (A, b, tol, [], diag (diag (A))); |
13021
d55d396a9a55
Add an implementation of the biconjugate gradient iterative method
Carlo de Falco <kingcrimson@tiscali.it>
parents:
diff
changeset
|
476 %! assert (x, ones (size (b)), 1e-7); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
477 |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
478 %!test |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
479 %! ## Check that if the preconditioner is singular, the method doesn't work |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
480 %! A = sparse (toeplitz ([2, 1, 0, 0, 0], [2, -1, 0, 0, 0])); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
481 %! b = ones (5,1); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
482 %! M = ones (5); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
483 %! [x, flag] = bicg (A, b, [], [], M); |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
484 %! assert (flag, 2); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
485 |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
486 %!test |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
487 %! ## If A singular, the algorithm doesn't work due to division by zero |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
488 %! A = ones (5); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
489 %! b = [1:5]'; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
490 %! [x, flag] = bicg (A, b); |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
491 %! assert (flag, 4); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
492 |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
493 %!test |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
494 %! ## test for a complex linear system |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
495 %! A = sparse (toeplitz ([2, 1, 0, 0, 0], [2, -1, 0, 0, 0]) + ... |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
496 %! 1i * toeplitz ([2, 1, 0, 0, 0], [2, -1, 0, 0, 0])); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
497 %! b = sum (A, 2); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
498 %! [x, flag] = bicg (A, b); |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
499 %! assert (flag, 0); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
500 |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
501 %!test |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
502 %! A = single (1); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
503 %! b = 1; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
504 %! [x, flag] = bicg (A, b); |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
505 %! assert (class (x), "single"); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
506 |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
507 %!test |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
508 %! A = 1; |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
509 %! b = single (1); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
510 %! [x, flag] = bicg (A, b); |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
511 %! assert (class (x), "single"); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
512 |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
513 %!test |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
514 %! A = single (1); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
515 %! b = single (1); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
516 %! [x, flag] = bicg (A, b); |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
517 %! assert (class (x), "single"); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
518 |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
519 %!test |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
520 %!function y = Afun (x, trans) |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
521 %! A = sparse (toeplitz ([2, 1, 0, 0], [2, -1, 0, 0])); |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
522 %! if (strcmp (trans, "notransp")) |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
523 %! y = A * x; |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
524 %! else |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
525 %! y = A' * x; |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
526 %! endif |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
527 %!endfunction |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
528 %! |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
529 %! [x, flag] = bicg ("Afun", [1; 2; 2; 3]); |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
530 %! assert (x, ones (4, 1), 1e-6); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
531 |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
532 %!test |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
533 %! ## unpreconditioned residual |
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
534 %! A = sparse (toeplitz ([2, 1, 0, 0, 0], [2, -1, 0, 0, 0])); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
535 %! b = sum (A, 2); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
536 %! M = magic (5); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
537 %! [x, flag, relres] = bicg (A, b, [], 2, M); |
26160
b9debf4436aa
bicg.m: Relax BIST test by 1eps to pass on some systems (bug #55132).
Rik <rik@octave.org>
parents:
25579
diff
changeset
|
538 %! assert (norm (b - A * x) / norm (b), 0, relres + eps); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
539 |
24863
009ed89dd0fa
test: make bicg and pcg tests conditional on CHOLMOD and UMFPACK
Mike Miller <mtmiller@octave.org>
parents:
24860
diff
changeset
|
540 ## Preconditioned technique |
009ed89dd0fa
test: make bicg and pcg tests conditional on CHOLMOD and UMFPACK
Mike Miller <mtmiller@octave.org>
parents:
24860
diff
changeset
|
541 %!testif HAVE_UMFPACK |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
542 %! A = sparse (toeplitz ([2, 1, 0, 0, 0], [2, -1, 0, 0, 0])); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
543 %! b = sum (A, 2); |
25159
b8279dd83664
Turn off warning about sparse lu factorization in bicg BIST test (bug #53390).
Marco Caliari <marco.caliari@univr.it>
parents:
25152
diff
changeset
|
544 %! warning ("off", "Octave:lu:sparse_input", "local"); |
24860
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
545 %! [M1, M2] = lu (A + eye (5)); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
546 %! [x, flag] = bicg (A, b, [], 1, M1, M2); |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
547 %! ## b has two columns! |
6266e321ef22
Merge code from SoCiS 2016 "Improve iterative methods" project (http://socis16octave-improveiterativemethods.blogspot.com)
Cristiano Dorigo <cristiano.dorigo@hotmail.it>
parents:
24534
diff
changeset
|
548 %! [y, flag] = bicg (M1 \ A / M2, [M1 \ b, M2' \ b], [], 1); |
25152
ad37f04a2eb7
bicg.m: Overhaul GSOC-improved code to conform to Octave conventions.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
549 %! assert (x, M2 \ y, 8 * eps); |