Mercurial > octave-nkf
annotate doc/interpreter/sparse.txi @ 16816:12005245b645
doc: Periodic grammarcheck of documentation.
* doc/interpreter/basics.txi, doc/interpreter/container.txi,
doc/interpreter/contrib.txi, doc/interpreter/diagperm.txi,
doc/interpreter/errors.txi, doc/interpreter/install.txi,
doc/interpreter/sparse.txi, libinterp/corefcn/ellipj.cc,
libinterp/corefcn/mappers.cc, libinterp/corefcn/regexp.cc,
libinterp/corefcn/tril.cc, libinterp/dldfcn/__init_fltk__.cc,
libinterp/dldfcn/fftw.cc, libinterp/interpfcn/dirfns.cc,
libinterp/interpfcn/input.cc, libinterp/octave-value/ov-fcn-inline.cc,
libinterp/octave-value/ov-struct.cc, scripts/@ftp/cd.m,
scripts/general/interp1.m, scripts/general/num2str.m, scripts/image/ind2rgb.m,
scripts/image/rgb2ind.m, scripts/io/importdata.m, scripts/io/textread.m,
scripts/java/javamem.m, scripts/linear-algebra/condest.m,
scripts/linear-algebra/onenormest.m, scripts/miscellaneous/error_ids.m,
scripts/miscellaneous/getfield.m, scripts/miscellaneous/setfield.m,
scripts/plot/area.m, scripts/plot/pcolor.m, scripts/plot/stairs.m,
scripts/set/powerset.m, scripts/sparse/bicg.m, scripts/sparse/bicgstab.m,
scripts/sparse/cgs.m, scripts/specfun/ellipke.m,
scripts/special-matrix/gallery.m, scripts/strings/strjoin.m,
scripts/strings/strsplit.m, scripts/testfun/__have_feature__.m,
scripts/testfun/__printf_assert__.m, scripts/testfun/__prog_output_assert__.m,
scripts/testfun/__run_test_suite__.m: grammarcheck documentation.
author | Rik <rik@octave.org> |
---|---|
date | Sat, 22 Jun 2013 19:47:32 -0700 |
parents | b157ba28f123 |
children | a4969508008e |
rev | line source |
---|---|
14138
72c96de7a403
maint: update copyright notices for 2012
John W. Eaton <jwe@octave.org>
parents:
14116
diff
changeset
|
1 @c Copyright (C) 2004-2012 David Bateman |
7018 | 2 @c |
3 @c This file is part of Octave. | |
4 @c | |
5 @c Octave is free software; you can redistribute it and/or modify it | |
6 @c under the terms of the GNU General Public License as published by the | |
7 @c Free Software Foundation; either version 3 of the License, or (at | |
8 @c your option) any later version. | |
9 @c | |
10 @c Octave is distributed in the hope that it will be useful, but WITHOUT | |
11 @c ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or | |
12 @c FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License | |
13 @c for more details. | |
14 @c | |
15 @c You should have received a copy of the GNU General Public License | |
16 @c along with Octave; see the file COPYING. If not, see | |
17 @c <http://www.gnu.org/licenses/>. | |
5164 | 18 |
5648 | 19 @ifhtml |
20 @set htmltex | |
21 @end ifhtml | |
22 @iftex | |
23 @set htmltex | |
24 @end iftex | |
25 | |
5164 | 26 @node Sparse Matrices |
27 @chapter Sparse Matrices | |
28 | |
29 @menu | |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
30 * Basics:: Creation and Manipulation of Sparse Matrices |
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
31 * Sparse Linear Algebra:: Linear Algebra on Sparse Matrices |
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
32 * Iterative Techniques:: Iterative Techniques |
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
33 * Real Life Example:: Using Sparse Matrices |
5164 | 34 @end menu |
35 | |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
36 @node Basics |
5164 | 37 @section The Creation and Manipulation of Sparse Matrices |
38 | |
39 The size of mathematical problems that can be treated at any particular | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
40 time is generally limited by the available computing resources. Both, |
5164 | 41 the speed of the computer and its available memory place limitation on |
42 the problem size. | |
43 | |
5506 | 44 There are many classes of mathematical problems which give rise to |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
45 matrices, where a large number of the elements are zero. In this case |
5164 | 46 it makes sense to have a special matrix type to handle this class of |
47 problems where only the non-zero elements of the matrix are | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
48 stored. Not only does this reduce the amount of memory to store the |
5164 | 49 matrix, but it also means that operations on this type of matrix can |
50 take advantage of the a-priori knowledge of the positions of the | |
51 non-zero elements to accelerate their calculations. | |
52 | |
53 A matrix type that stores only the non-zero elements is generally called | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
54 sparse. It is the purpose of this document to discuss the basics of the |
5164 | 55 storage and creation of sparse matrices and the fundamental operations |
56 on them. | |
57 | |
58 @menu | |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
59 * Storage of Sparse Matrices:: |
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
60 * Creating Sparse Matrices:: |
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
61 * Information:: |
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
62 * Operators and Functions:: |
5164 | 63 @end menu |
64 | |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
65 @node Storage of Sparse Matrices |
5164 | 66 @subsection Storage of Sparse Matrices |
67 | |
68 It is not strictly speaking necessary for the user to understand how | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
69 sparse matrices are stored. However, such an understanding will help |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
70 to get an understanding of the size of sparse matrices. Understanding |
5164 | 71 the storage technique is also necessary for those users wishing to |
72 create their own oct-files. | |
73 | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
74 There are many different means of storing sparse matrix data. What all |
5648 | 75 of the methods have in common is that they attempt to reduce the complexity |
5164 | 76 and storage given a-priori knowledge of the particular class of problems |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
77 that will be solved. A good summary of the available techniques for storing |
10791
3140cb7a05a1
Add spellchecker scripts for Octave and run spellcheck of documentation
Rik <octave@nomad.inbox5.com>
parents:
10711
diff
changeset
|
78 sparse matrix is given by Saad @footnote{Y. Saad "SPARSKIT: A basic toolkit |
5164 | 79 for sparse matrix computation", 1994, |
6620 | 80 @url{http://www-users.cs.umn.edu/~saad/software/SPARSKIT/paper.ps}}. |
5164 | 81 With full matrices, knowledge of the point of an element of the matrix |
82 within the matrix is implied by its position in the computers memory. | |
83 However, this is not the case for sparse matrices, and so the positions | |
84 of the non-zero elements of the matrix must equally be stored. | |
85 | |
86 An obvious way to do this is by storing the elements of the matrix as | |
5506 | 87 triplets, with two elements being their position in the array |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
88 (rows and column) and the third being the data itself. This is conceptually |
5164 | 89 easy to grasp, but requires more storage than is strictly needed. |
90 | |
5648 | 91 The storage technique used within Octave is the compressed column |
16816
12005245b645
doc: Periodic grammarcheck of documentation.
Rik <rik@octave.org>
parents:
16792
diff
changeset
|
92 format. It is similar to the Yale format. |
16600
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
93 @footnote{@url{http://en.wikipedia.org/wiki/Sparse_matrix#Yale_format}} |
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
94 In this format the position of each element in a row and the data are |
16816
12005245b645
doc: Periodic grammarcheck of documentation.
Rik <rik@octave.org>
parents:
16792
diff
changeset
|
95 stored as previously. However, if we assume that all elements in the |
16600
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
96 same column are stored adjacent in the computers memory, then we only |
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
97 need to store information on the number of non-zero elements in each |
16816
12005245b645
doc: Periodic grammarcheck of documentation.
Rik <rik@octave.org>
parents:
16792
diff
changeset
|
98 column, rather than their positions. Thus assuming that the matrix has |
16600
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
99 more non-zero elements than there are columns in the matrix, we win in |
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
100 terms of the amount of memory used. |
5164 | 101 |
102 In fact, the column index contains one more element than the number of | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
103 columns, with the first element always being zero. The advantage of |
7001 | 104 this is a simplification in the code, in that there is no special case |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
105 for the first or last columns. A short example, demonstrating this in |
5164 | 106 C is. |
107 | |
108 @example | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
109 @group |
5164 | 110 for (j = 0; j < nc; j++) |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
111 for (i = cidx(j); i < cidx(j+1); i++) |
5648 | 112 printf ("non-zero element (%i,%i) is %d\n", |
10599
d0e0bb2ebebb
Remove tabs in .txi files causing problems with pdf formatting.
Rik <octave@nomad.inbox5.com>
parents:
9906
diff
changeset
|
113 ridx(i), j, data(i)); |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
114 @end group |
5164 | 115 @end example |
116 | |
117 A clear understanding might be had by considering an example of how the | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
118 above applies to an example matrix. Consider the matrix |
5164 | 119 |
120 @example | |
121 @group | |
122 1 2 0 0 | |
123 0 0 0 3 | |
124 0 0 0 4 | |
125 @end group | |
126 @end example | |
127 | |
128 The non-zero elements of this matrix are | |
129 | |
130 @example | |
131 @group | |
132 (1, 1) @result{} 1 | |
133 (1, 2) @result{} 2 | |
134 (2, 4) @result{} 3 | |
135 (3, 4) @result{} 4 | |
136 @end group | |
137 @end example | |
138 | |
139 This will be stored as three vectors @var{cidx}, @var{ridx} and @var{data}, | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
140 representing the column indexing, row indexing and data respectively. The |
5164 | 141 contents of these three vectors for the above matrix will be |
142 | |
143 @example | |
144 @group | |
5506 | 145 @var{cidx} = [0, 1, 2, 2, 4] |
5164 | 146 @var{ridx} = [0, 0, 1, 2] |
147 @var{data} = [1, 2, 3, 4] | |
148 @end group | |
149 @end example | |
150 | |
151 Note that this is the representation of these elements with the first row | |
5648 | 152 and column assumed to start at zero, while in Octave itself the row and |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
153 column indexing starts at one. Thus the number of elements in the |
5164 | 154 @var{i}-th column is given by @code{@var{cidx} (@var{i} + 1) - |
155 @var{cidx} (@var{i})}. | |
156 | |
5648 | 157 Although Octave uses a compressed column format, it should be noted |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
158 that compressed row formats are equally possible. However, in the |
5648 | 159 context of mixed operations between mixed sparse and dense matrices, |
160 it makes sense that the elements of the sparse matrices are in the | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
161 same order as the dense matrices. Octave stores dense matrices in |
5648 | 162 column major ordering, and so sparse matrices are equally stored in |
163 this manner. | |
5164 | 164 |
5324 | 165 A further constraint on the sparse matrix storage used by Octave is that |
5164 | 166 all elements in the rows are stored in increasing order of their row |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
167 index, which makes certain operations faster. However, it imposes |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
168 the need to sort the elements on the creation of sparse matrices. Having |
7001 | 169 disordered elements is potentially an advantage in that it makes operations |
5164 | 170 such as concatenating two sparse matrices together easier and faster, however |
171 it adds complexity and speed problems elsewhere. | |
172 | |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
173 @node Creating Sparse Matrices |
5164 | 174 @subsection Creating Sparse Matrices |
175 | |
176 There are several means to create sparse matrix. | |
177 | |
178 @table @asis | |
179 @item Returned from a function | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
180 There are many functions that directly return sparse matrices. These include |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7505
diff
changeset
|
181 @dfn{speye}, @dfn{sprand}, @dfn{diag}, etc. |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
182 |
5164 | 183 @item Constructed from matrices or vectors |
184 The function @dfn{sparse} allows a sparse matrix to be constructed from | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
185 three vectors representing the row, column and data. Alternatively, the |
5164 | 186 function @dfn{spconvert} uses a three column matrix format to allow easy |
187 importation of data from elsewhere. | |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
188 |
5164 | 189 @item Created and then filled |
190 The function @dfn{sparse} or @dfn{spalloc} can be used to create an empty | |
191 matrix that is then filled by the user | |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
192 |
5164 | 193 @item From a user binary program |
194 The user can directly create the sparse matrix within an oct-file. | |
195 @end table | |
196 | |
197 There are several basic functions to return specific sparse | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
198 matrices. For example the sparse identity matrix, is a matrix that is |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
199 often needed. It therefore has its own function to create it as |
5164 | 200 @code{speye (@var{n})} or @code{speye (@var{r}, @var{c})}, which |
201 creates an @var{n}-by-@var{n} or @var{r}-by-@var{c} sparse identity | |
202 matrix. | |
203 | |
204 Another typical sparse matrix that is often needed is a random distribution | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
205 of random elements. The functions @dfn{sprand} and @dfn{sprandn} perform |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
206 this for uniform and normal random distributions of elements. They have exactly |
5506 | 207 the same calling convention, where @code{sprand (@var{r}, @var{c}, @var{d})}, |
208 creates an @var{r}-by-@var{c} sparse matrix with a density of filled | |
5164 | 209 elements of @var{d}. |
210 | |
7001 | 211 Other functions of interest that directly create sparse matrices, are |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7505
diff
changeset
|
212 @dfn{diag} or its generalization @dfn{spdiags}, that can take the |
5164 | 213 definition of the diagonals of the matrix and create the sparse matrix |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
214 that corresponds to this. For example, |
5164 | 215 |
216 @example | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
217 s = diag (sparse (randn (1,n)), -1); |
5164 | 218 @end example |
219 | |
10846
a4f482e66b65
Grammarcheck more of the documentation.
Rik <octave@nomad.inbox5.com>
parents:
10828
diff
changeset
|
220 @noindent |
5164 | 221 creates a sparse (@var{n}+1)-by-(@var{n}+1) sparse matrix with a single |
222 diagonal defined. | |
223 | |
6620 | 224 @DOCSTRING(spdiags) |
225 | |
226 @DOCSTRING(speye) | |
227 | |
228 @DOCSTRING(spones) | |
229 | |
230 @DOCSTRING(sprand) | |
231 | |
232 @DOCSTRING(sprandn) | |
233 | |
234 @DOCSTRING(sprandsym) | |
235 | |
5164 | 236 The recommended way for the user to create a sparse matrix, is to create |
5648 | 237 two vectors containing the row and column index of the data and a third |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
238 vector of the same size containing the data to be stored. For example, |
5164 | 239 |
240 @example | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
241 @group |
6421 | 242 ri = ci = d = []; |
243 for j = 1:c | |
13198
0a1774f1a70d
Update example in sparse.txi to use new calling form of randperm
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
11593
diff
changeset
|
244 ri = [ri; randperm(r,n)']; |
6421 | 245 ci = [ci; j*ones(n,1)]; |
246 d = [d; rand(n,1)]; | |
247 endfor | |
248 s = sparse (ri, ci, d, r, c); | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
249 @end group |
5164 | 250 @end example |
251 | |
10846
a4f482e66b65
Grammarcheck more of the documentation.
Rik <octave@nomad.inbox5.com>
parents:
10828
diff
changeset
|
252 @noindent |
6421 | 253 creates an @var{r}-by-@var{c} sparse matrix with a random distribution |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
254 of @var{n} (<@var{r}) elements per column. The elements of the vectors |
6421 | 255 do not need to be sorted in any particular order as Octave will sort |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
256 them prior to storing the data. However, pre-sorting the data will |
6421 | 257 make the creation of the sparse matrix faster. |
5164 | 258 |
259 The function @dfn{spconvert} takes a three or four column real matrix. | |
260 The first two columns represent the row and column index respectively and | |
261 the third and four columns, the real and imaginary parts of the sparse | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
262 matrix. The matrix can contain zero elements and the elements can be |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
263 sorted in any order. Adding zero elements is a convenient way to define |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
264 the size of the sparse matrix. For example: |
5164 | 265 |
266 @example | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
267 @group |
5164 | 268 s = spconvert ([1 2 3 4; 1 3 4 4; 1 2 3 0]') |
269 @result{} Compressed Column Sparse (rows=4, cols=4, nnz=3) | |
270 (1 , 1) -> 1 | |
271 (2 , 3) -> 2 | |
272 (3 , 4) -> 3 | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
273 @end group |
5164 | 274 @end example |
275 | |
276 An example of creating and filling a matrix might be | |
277 | |
278 @example | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
279 @group |
5164 | 280 k = 5; |
281 nz = r * k; | |
282 s = spalloc (r, c, nz) | |
283 for j = 1:c | |
284 idx = randperm (r); | |
5648 | 285 s (:, j) = [zeros(r - k, 1); ... |
286 rand(k, 1)] (idx); | |
5164 | 287 endfor |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
288 @end group |
5164 | 289 @end example |
290 | |
5324 | 291 It should be noted, that due to the way that the Octave |
5164 | 292 assignment functions are written that the assignment will reallocate |
5506 | 293 the memory used by the sparse matrix at each iteration of the above loop. |
294 Therefore the @dfn{spalloc} function ignores the @var{nz} argument and | |
10791
3140cb7a05a1
Add spellchecker scripts for Octave and run spellcheck of documentation
Rik <octave@nomad.inbox5.com>
parents:
10711
diff
changeset
|
295 does not pre-assign the memory for the matrix. Therefore, it is vitally |
5648 | 296 important that code using to above structure should be vectorized |
297 as much as possible to minimize the number of assignments and reduce the | |
5164 | 298 number of memory allocations. |
299 | |
6620 | 300 @DOCSTRING(full) |
301 | |
302 @DOCSTRING(spalloc) | |
303 | |
304 @DOCSTRING(sparse) | |
305 | |
306 @DOCSTRING(spconvert) | |
307 | |
8106
8a42498edb30
Clarify doc for sparse function
David Bateman <dbateman@free.fr>
parents:
7681
diff
changeset
|
308 The above problem of memory reallocation can be avoided in |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
309 oct-files. However, the construction of a sparse matrix from an oct-file |
8828 | 310 is more complex than can be discussed here, and |
8106
8a42498edb30
Clarify doc for sparse function
David Bateman <dbateman@free.fr>
parents:
7681
diff
changeset
|
311 you are referred to chapter @ref{Dynamically Linked Functions}, to have |
8a42498edb30
Clarify doc for sparse function
David Bateman <dbateman@free.fr>
parents:
7681
diff
changeset
|
312 a full description of the techniques involved. |
5164 | 313 |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
314 @node Information |
5648 | 315 @subsection Finding out Information about Sparse Matrices |
316 | |
317 There are a number of functions that allow information concerning | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
318 sparse matrices to be obtained. The most basic of these is |
5648 | 319 @dfn{issparse} that identifies whether a particular Octave object is |
320 in fact a sparse matrix. | |
321 | |
322 Another very basic function is @dfn{nnz} that returns the number of | |
323 non-zero entries there are in a sparse matrix, while the function | |
324 @dfn{nzmax} returns the amount of storage allocated to the sparse | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
325 matrix. Note that Octave tends to crop unused memory at the first |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
326 opportunity for sparse objects. There are some cases of user created |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
327 sparse objects where the value returned by @dfn{nzmax} will not be |
5648 | 328 the same as @dfn{nnz}, but in general they will give the same |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
329 result. The function @dfn{spstats} returns some basic statistics on |
5648 | 330 the columns of a sparse matrix including the number of elements, the |
331 mean and the variance of each column. | |
332 | |
6620 | 333 @DOCSTRING(issparse) |
334 | |
335 @DOCSTRING(nnz) | |
336 | |
337 @DOCSTRING(nonzeros) | |
338 | |
339 @DOCSTRING(nzmax) | |
340 | |
341 @DOCSTRING(spstats) | |
342 | |
5648 | 343 When solving linear equations involving sparse matrices Octave |
344 determines the means to solve the equation based on the type of the | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
345 matrix as discussed in @ref{Sparse Linear Algebra}. Octave probes the |
5648 | 346 matrix type when the div (/) or ldiv (\) operator is first used with |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
347 the matrix and then caches the type. However the @dfn{matrix_type} |
5648 | 348 function can be used to determine the type of the sparse matrix prior |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
349 to use of the div or ldiv operators. For example, |
5648 | 350 |
351 @example | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
352 @group |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
353 a = tril (sprandn (1024, 1024, 0.02), -1) ... |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
354 + speye (1024); |
5648 | 355 matrix_type (a); |
356 ans = Lower | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
357 @end group |
5648 | 358 @end example |
359 | |
10846
a4f482e66b65
Grammarcheck more of the documentation.
Rik <octave@nomad.inbox5.com>
parents:
10828
diff
changeset
|
360 @noindent |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
361 shows that Octave correctly determines the matrix type for lower |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
362 triangular matrices. @dfn{matrix_type} can also be used to force |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
363 the type of a matrix to be a particular type. For example: |
5648 | 364 |
365 @example | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
366 @group |
5648 | 367 a = matrix_type (tril (sprandn (1024, ... |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
368 1024, 0.02), -1) + speye (1024), "Lower"); |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
369 @end group |
5648 | 370 @end example |
371 | |
372 This allows the cost of determining the matrix type to be | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
373 avoided. However, incorrectly defining the matrix type will result in |
5648 | 374 incorrect results from solutions of linear equations, and so it is |
375 entirely the responsibility of the user to correctly identify the | |
376 matrix type | |
377 | |
378 There are several graphical means of finding out information about | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
379 sparse matrices. The first is the @dfn{spy} command, which displays |
5648 | 380 the structure of the non-zero elements of the |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
381 matrix. @xref{fig:spmatrix}, for an example of the use of |
5704 | 382 @dfn{spy}. More advanced graphical information can be obtained with the |
5648 | 383 @dfn{treeplot}, @dfn{etreeplot} and @dfn{gplot} commands. |
384 | |
385 @float Figure,fig:spmatrix | |
9088
77e71f3da3d6
Fix documentation image printing under new development code
Rik <rdrider0-list@yahoo.com>
parents:
9066
diff
changeset
|
386 @center @image{spmatrix,4in} |
5648 | 387 @caption{Structure of simple sparse matrix.} |
388 @end float | |
389 | |
390 One use of sparse matrices is in graph theory, where the | |
7001 | 391 interconnections between nodes are represented as an adjacency |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
392 matrix. That is, if the i-th node in a graph is connected to the j-th |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
393 node. Then the ij-th node (and in the case of undirected graphs the |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
394 ji-th node) of the sparse adjacency matrix is non-zero. If each node |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
395 is then associated with a set of coordinates, then the @dfn{gplot} |
5648 | 396 command can be used to graphically display the interconnections |
397 between nodes. | |
398 | |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
399 As a trivial example of the use of @dfn{gplot} consider the example, |
5648 | 400 |
401 @example | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
402 @group |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
403 A = sparse ([2,6,1,3,2,4,3,5,4,6,1,5], |
5648 | 404 [1,1,2,2,3,3,4,4,5,5,6,6],1,6,6); |
405 xy = [0,4,8,6,4,2;5,0,5,7,5,7]'; | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
406 gplot (A,xy) |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
407 @end group |
5648 | 408 @end example |
409 | |
10846
a4f482e66b65
Grammarcheck more of the documentation.
Rik <octave@nomad.inbox5.com>
parents:
10828
diff
changeset
|
410 @noindent |
5648 | 411 which creates an adjacency matrix @code{A} where node 1 is connected |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
412 to nodes 2 and 6, node 2 with nodes 1 and 3, etc. The coordinates of |
5648 | 413 the nodes are given in the n-by-2 matrix @code{xy}. |
414 @ifset htmltex | |
415 @xref{fig:gplot}. | |
416 | |
417 @float Figure,fig:gplot | |
9088
77e71f3da3d6
Fix documentation image printing under new development code
Rik <rdrider0-list@yahoo.com>
parents:
9066
diff
changeset
|
418 @center @image{gplot,4in} |
5648 | 419 @caption{Simple use of the @dfn{gplot} command.} |
420 @end float | |
421 @end ifset | |
422 | |
11593
1577c6f80926
Use non-breaking spaces between certain adjectives and their nouns in docstrings.
Rik <octave@nomad.inbox5.com>
parents:
11573
diff
changeset
|
423 The dependencies between the nodes of a Cholesky@tie{}factorization can be |
5648 | 424 calculated in linear time without explicitly needing to calculate the |
11593
1577c6f80926
Use non-breaking spaces between certain adjectives and their nouns in docstrings.
Rik <octave@nomad.inbox5.com>
parents:
11573
diff
changeset
|
425 Cholesky@tie{}factorization by the @code{etree} command. This command |
5648 | 426 returns the elimination tree of the matrix and can be displayed |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
427 graphically by the command @code{treeplot (etree (A))} if @code{A} is |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
428 symmetric or @code{treeplot (etree (A+A'))} otherwise. |
5648 | 429 |
6620 | 430 @DOCSTRING(spy) |
431 | |
432 @DOCSTRING(etree) | |
433 | |
434 @DOCSTRING(etreeplot) | |
435 | |
436 @DOCSTRING(gplot) | |
437 | |
438 @DOCSTRING(treeplot) | |
439 | |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
440 @DOCSTRING(treelayout) |
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
441 |
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
442 @node Operators and Functions |
5164 | 443 @subsection Basic Operators and Functions on Sparse Matrices |
444 | |
445 @menu | |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
446 * Sparse Functions:: |
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
447 * Return Types of Operators and Functions:: |
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
448 * Mathematical Considerations:: |
5164 | 449 @end menu |
450 | |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
451 @node Sparse Functions |
5648 | 452 @subsubsection Sparse Functions |
453 | |
11396
7b563cf94d8d
Remove documentation on deprecated function dispatch
Rik <octave@nomad.inbox5.com>
parents:
11392
diff
changeset
|
454 Many Octave functions have been overloaded to work with either sparse or full |
7b563cf94d8d
Remove documentation on deprecated function dispatch
Rik <octave@nomad.inbox5.com>
parents:
11392
diff
changeset
|
455 matrices. There is no difference in calling convention when using an |
7b563cf94d8d
Remove documentation on deprecated function dispatch
Rik <octave@nomad.inbox5.com>
parents:
11392
diff
changeset
|
456 overloaded function with a sparse matrix, however, there is also no access to |
7b563cf94d8d
Remove documentation on deprecated function dispatch
Rik <octave@nomad.inbox5.com>
parents:
11392
diff
changeset
|
457 potentially sparse-specific features. At any time the sparse matrix specific |
7b563cf94d8d
Remove documentation on deprecated function dispatch
Rik <octave@nomad.inbox5.com>
parents:
11392
diff
changeset
|
458 version of a function can be used by explicitly calling its function name. |
5648 | 459 |
11396
7b563cf94d8d
Remove documentation on deprecated function dispatch
Rik <octave@nomad.inbox5.com>
parents:
11392
diff
changeset
|
460 The table below lists all of the sparse functions of Octave. Note that the |
7b563cf94d8d
Remove documentation on deprecated function dispatch
Rik <octave@nomad.inbox5.com>
parents:
11392
diff
changeset
|
461 names of the specific sparse forms of the functions are typically the same as |
7b563cf94d8d
Remove documentation on deprecated function dispatch
Rik <octave@nomad.inbox5.com>
parents:
11392
diff
changeset
|
462 the general versions with a @dfn{sp} prefix. In the table below, and in the |
7b563cf94d8d
Remove documentation on deprecated function dispatch
Rik <octave@nomad.inbox5.com>
parents:
11392
diff
changeset
|
463 rest of this article, the specific sparse versions of functions are used. |
6620 | 464 |
465 @c Table includes in comments the missing sparse functions | |
5648 | 466 |
467 @table @asis | |
468 @item Generate sparse matrices: | |
469 @dfn{spalloc}, @dfn{spdiags}, @dfn{speye}, @dfn{sprand}, | |
470 @dfn{sprandn}, @dfn{sprandsym} | |
471 | |
472 @item Sparse matrix conversion: | |
7505
f5005d9510f4
Remove dispatched sparse functions and treat in the generic versions of the functions
David Bateman <dbateman@free.fr>
parents:
7455
diff
changeset
|
473 @dfn{full}, @dfn{sparse}, @dfn{spconvert} |
5648 | 474 |
475 @item Manipulate sparse matrices | |
476 @dfn{issparse}, @dfn{nnz}, @dfn{nonzeros}, @dfn{nzmax}, | |
6620 | 477 @dfn{spfun}, @dfn{spones}, @dfn{spy} |
5164 | 478 |
5648 | 479 @item Graph Theory: |
480 @dfn{etree}, @dfn{etreeplot}, @dfn{gplot}, | |
6620 | 481 @dfn{treeplot} |
482 @c @dfn{treelayout} | |
5648 | 483 |
484 @item Sparse matrix reordering: | |
7619 | 485 @dfn{amd}, @dfn{ccolamd}, @dfn{colamd}, @dfn{colperm}, @dfn{csymamd}, |
6620 | 486 @dfn{dmperm}, @dfn{symamd}, @dfn{randperm}, @dfn{symrcm} |
5648 | 487 |
488 @item Linear algebra: | |
8417
654bcfb937bf
Add the eigs and svds functions
David Bateman <dbateman@free.fr>
parents:
8347
diff
changeset
|
489 @dfn{condest}, @dfn{eigs}, @dfn{matrix_type}, @dfn{normest}, @dfn{sprank}, |
654bcfb937bf
Add the eigs and svds functions
David Bateman <dbateman@free.fr>
parents:
8347
diff
changeset
|
490 @dfn{spaugment}, @dfn{svds} |
5648 | 491 |
492 @item Iterative techniques: | |
6620 | 493 @dfn{luinc}, @dfn{pcg}, @dfn{pcr} |
494 @c @dfn{bicg}, @dfn{bicgstab}, @dfn{cholinc}, @dfn{cgs}, @dfn{gmres}, | |
495 @c @dfn{lsqr}, @dfn{minres}, @dfn{qmr}, @dfn{symmlq} | |
5648 | 496 |
497 @item Miscellaneous: | |
7515
f3c00dc0912b
Eliminate the rest of the dispatched sparse functions
David Bateman <dbateman@free.fr>
parents:
7505
diff
changeset
|
498 @dfn{spparms}, @dfn{symbfact}, @dfn{spstats} |
5648 | 499 @end table |
500 | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
501 In addition all of the standard Octave mapper functions (i.e., basic |
9088
77e71f3da3d6
Fix documentation image printing under new development code
Rik <rdrider0-list@yahoo.com>
parents:
9066
diff
changeset
|
502 math functions that take a single argument) such as @dfn{abs}, etc. |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
503 can accept sparse matrices. The reader is referred to the documentation |
5648 | 504 supplied with these functions within Octave itself for further |
505 details. | |
5164 | 506 |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
507 @node Return Types of Operators and Functions |
5164 | 508 @subsubsection The Return Types of Operators and Functions |
509 | |
5506 | 510 The two basic reasons to use sparse matrices are to reduce the memory |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
511 usage and to not have to do calculations on zero elements. The two are |
5164 | 512 closely related in that the computation time on a sparse matrix operator |
5506 | 513 or function is roughly linear with the number of non-zero elements. |
5164 | 514 |
515 Therefore, there is a certain density of non-zero elements of a matrix | |
516 where it no longer makes sense to store it as a sparse matrix, but rather | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
517 as a full matrix. For this reason operators and functions that have a |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
518 high probability of returning a full matrix will always return one. For |
5164 | 519 example adding a scalar constant to a sparse matrix will almost always |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
520 make it a full matrix, and so the example, |
5164 | 521 |
522 @example | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
523 @group |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
524 speye (3) + 0 |
5164 | 525 @result{} 1 0 0 |
526 0 1 0 | |
527 0 0 1 | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
528 @end group |
5164 | 529 @end example |
530 | |
10846
a4f482e66b65
Grammarcheck more of the documentation.
Rik <octave@nomad.inbox5.com>
parents:
10828
diff
changeset
|
531 @noindent |
7330 | 532 returns a full matrix as can be seen. |
533 | |
534 | |
535 Additionally, if @code{sparse_auto_mutate} is true, all sparse functions | |
536 test the amount of memory occupied by the sparse matrix to see if the | |
537 amount of storage used is larger than the amount used by the full | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
538 equivalent. Therefore @code{speye (2) * 1} will return a full matrix as |
5164 | 539 the memory used is smaller for the full version than the sparse version. |
540 | |
541 As all of the mixed operators and functions between full and sparse | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
542 matrices exist, in general this does not cause any problems. However, |
5164 | 543 one area where it does cause a problem is where a sparse matrix is |
544 promoted to a full matrix, where subsequent operations would resparsify | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
545 the matrix. Such cases are rare, but can be artificially created, for |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
546 example @code{(fliplr (speye (3)) + speye (3)) - speye (3)} gives a full |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
547 matrix when it should give a sparse one. In general, where such cases |
5164 | 548 occur, they impose only a small memory penalty. |
549 | |
5648 | 550 There is however one known case where this behavior of Octave's |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
551 sparse matrices will cause a problem. That is in the handling of the |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
552 @dfn{diag} function. Whether @dfn{diag} returns a sparse or full matrix |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
553 depending on the type of its input arguments. So |
5164 | 554 |
555 @example | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
556 a = diag (sparse ([1,2,3]), -1); |
5164 | 557 @end example |
558 | |
10846
a4f482e66b65
Grammarcheck more of the documentation.
Rik <octave@nomad.inbox5.com>
parents:
10828
diff
changeset
|
559 @noindent |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
560 should return a sparse matrix. To ensure this actually happens, the |
5164 | 561 @dfn{sparse} function, and other functions based on it like @dfn{speye}, |
562 always returns a sparse matrix, even if the memory used will be larger | |
563 than its full representation. | |
564 | |
7330 | 565 @DOCSTRING(sparse_auto_mutate) |
566 | |
567 Note that the @code{sparse_auto_mutate} option is incompatible with | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
568 @sc{matlab}, and so it is off by default. |
7330 | 569 |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
570 @node Mathematical Considerations |
5164 | 571 @subsubsection Mathematical Considerations |
572 | |
573 The attempt has been made to make sparse matrices behave in exactly the | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
574 same manner as there full counterparts. However, there are certain differences |
5164 | 575 and especially differences with other products sparse implementations. |
576 | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
577 Firstly, the "./" and ".^" operators must be used with care. Consider what |
5164 | 578 the examples |
579 | |
580 @example | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
581 @group |
5164 | 582 s = speye (4); |
583 a1 = s .^ 2; | |
584 a2 = s .^ s; | |
585 a3 = s .^ -2; | |
586 a4 = s ./ 2; | |
587 a5 = 2 ./ s; | |
588 a6 = s ./ s; | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
589 @end group |
5164 | 590 @end example |
591 | |
10846
a4f482e66b65
Grammarcheck more of the documentation.
Rik <octave@nomad.inbox5.com>
parents:
10828
diff
changeset
|
592 @noindent |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
593 will give. The first example of @var{s} raised to the power of 2 causes |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
594 no problems. However @var{s} raised element-wise to itself involves a |
6431 | 595 large number of terms @code{0 .^ 0} which is 1. There @code{@var{s} .^ |
5164 | 596 @var{s}} is a full matrix. |
597 | |
8347
fa78cb8d8a5c
corrections for typos
Brian Gough<bjg@network-theory.co.uk>
parents:
8286
diff
changeset
|
598 Likewise @code{@var{s} .^ -2} involves terms like @code{0 .^ -2} which |
5164 | 599 is infinity, and so @code{@var{s} .^ -2} is equally a full matrix. |
600 | |
601 For the "./" operator @code{@var{s} ./ 2} has no problems, but | |
602 @code{2 ./ @var{s}} involves a large number of infinity terms as well | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
603 and is equally a full matrix. The case of @code{@var{s} ./ @var{s}} |
5164 | 604 involves terms like @code{0 ./ 0} which is a @code{NaN} and so this |
605 is equally a full matrix with the zero elements of @var{s} filled with | |
606 @code{NaN} values. | |
607 | |
5648 | 608 The above behavior is consistent with full matrices, but is not |
5164 | 609 consistent with sparse implementations in other products. |
610 | |
611 A particular problem of sparse matrices comes about due to the fact that | |
612 as the zeros are not stored, the sign-bit of these zeros is equally not | |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
613 stored. In certain cases the sign-bit of zero is important. For example: |
5164 | 614 |
615 @example | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
616 @group |
5164 | 617 a = 0 ./ [-1, 1; 1, -1]; |
618 b = 1 ./ a | |
619 @result{} -Inf Inf | |
620 Inf -Inf | |
621 c = 1 ./ sparse (a) | |
622 @result{} Inf Inf | |
623 Inf Inf | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
624 @end group |
5164 | 625 @end example |
626 | |
5648 | 627 To correct this behavior would mean that zero elements with a negative |
5164 | 628 sign-bit would need to be stored in the matrix to ensure that their |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
629 sign-bit was respected. This is not done at this time, for reasons of |
6750 | 630 efficiency, and so the user is warned that calculations where the sign-bit |
5164 | 631 of zero is important must not be done using sparse matrices. |
632 | |
5648 | 633 In general any function or operator used on a sparse matrix will |
634 result in a sparse matrix with the same or a larger number of non-zero | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
635 elements than the original matrix. This is particularly true for the |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
636 important case of sparse matrix factorizations. The usual way to |
5648 | 637 address this is to reorder the matrix, such that its factorization is |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
638 sparser than the factorization of the original matrix. That is the |
5648 | 639 factorization of @code{L * U = P * S * Q} has sparser terms @code{L} |
640 and @code{U} than the equivalent factorization @code{L * U = S}. | |
641 | |
642 Several functions are available to reorder depending on the type of the | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
643 matrix to be factorized. If the matrix is symmetric positive-definite, |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
644 then @dfn{symamd} or @dfn{csymamd} should be used. Otherwise |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
645 @dfn{amd}, @dfn{colamd} or @dfn{ccolamd} should be used. For completeness |
5648 | 646 the reordering functions @dfn{colperm} and @dfn{randperm} are |
647 also available. | |
648 | |
8829 | 649 @xref{fig:simplematrix}, for an example of the structure of a simple |
5648 | 650 positive definite matrix. |
5506 | 651 |
5648 | 652 @float Figure,fig:simplematrix |
9088
77e71f3da3d6
Fix documentation image printing under new development code
Rik <rdrider0-list@yahoo.com>
parents:
9066
diff
changeset
|
653 @center @image{spmatrix,4in} |
5648 | 654 @caption{Structure of simple sparse matrix.} |
655 @end float | |
5506 | 656 |
11593
1577c6f80926
Use non-breaking spaces between certain adjectives and their nouns in docstrings.
Rik <octave@nomad.inbox5.com>
parents:
11573
diff
changeset
|
657 The standard Cholesky@tie{}factorization of this matrix can be |
5648 | 658 obtained by the same command that would be used for a full |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
659 matrix. This can be visualized with the command |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
660 @code{r = chol (A); spy (r);}. |
5652 | 661 @xref{fig:simplechol}. |
662 The original matrix had | |
5648 | 663 @ifinfo |
664 @ifnothtml | |
665 43 | |
666 @end ifnothtml | |
667 @end ifinfo | |
668 @ifset htmltex | |
669 598 | |
670 @end ifset | |
11593
1577c6f80926
Use non-breaking spaces between certain adjectives and their nouns in docstrings.
Rik <octave@nomad.inbox5.com>
parents:
11573
diff
changeset
|
671 non-zero terms, while this Cholesky@tie{}factorization has |
5648 | 672 @ifinfo |
673 @ifnothtml | |
674 71, | |
675 @end ifnothtml | |
676 @end ifinfo | |
677 @ifset htmltex | |
678 10200, | |
679 @end ifset | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
680 with only half of the symmetric matrix being stored. This |
5648 | 681 is a significant level of fill in, and although not an issue |
682 for such a small test case, can represents a large overhead | |
683 in working with other sparse matrices. | |
5164 | 684 |
5648 | 685 The appropriate sparsity preserving permutation of the original |
686 matrix is given by @dfn{symamd} and the factorization using this | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
687 reordering can be visualized using the command @code{q = symamd (A); |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
688 r = chol (A(q,q)); spy (r)}. This gives |
5648 | 689 @ifinfo |
690 @ifnothtml | |
691 29 | |
692 @end ifnothtml | |
693 @end ifinfo | |
694 @ifset htmltex | |
695 399 | |
696 @end ifset | |
697 non-zero terms which is a significant improvement. | |
5164 | 698 |
11593
1577c6f80926
Use non-breaking spaces between certain adjectives and their nouns in docstrings.
Rik <octave@nomad.inbox5.com>
parents:
11573
diff
changeset
|
699 The Cholesky@tie{}factorization itself can be used to determine the |
5648 | 700 appropriate sparsity preserving reordering of the matrix during the |
701 factorization, In that case this might be obtained with three return | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
702 arguments as @code{[r, p, q] = chol (A); spy (r)}. |
5164 | 703 |
5648 | 704 @float Figure,fig:simplechol |
9088
77e71f3da3d6
Fix documentation image printing under new development code
Rik <rdrider0-list@yahoo.com>
parents:
9066
diff
changeset
|
705 @center @image{spchol,4in} |
11593
1577c6f80926
Use non-breaking spaces between certain adjectives and their nouns in docstrings.
Rik <octave@nomad.inbox5.com>
parents:
11573
diff
changeset
|
706 @caption{Structure of the un-permuted Cholesky@tie{}factorization of the above matrix.} |
5648 | 707 @end float |
5164 | 708 |
5648 | 709 @float Figure,fig:simplecholperm |
9088
77e71f3da3d6
Fix documentation image printing under new development code
Rik <rdrider0-list@yahoo.com>
parents:
9066
diff
changeset
|
710 @center @image{spcholperm,4in} |
11593
1577c6f80926
Use non-breaking spaces between certain adjectives and their nouns in docstrings.
Rik <octave@nomad.inbox5.com>
parents:
11573
diff
changeset
|
711 @caption{Structure of the permuted Cholesky@tie{}factorization of the above matrix.} |
5648 | 712 @end float |
5164 | 713 |
5648 | 714 In the case of an asymmetric matrix, the appropriate sparsity |
715 preserving permutation is @dfn{colamd} and the factorization using | |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
716 this reordering can be visualized using the command |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
717 @code{q = colamd (A); [l, u, p] = lu (A(:,q)); spy (l+u)}. |
5164 | 718 |
5648 | 719 Finally, Octave implicitly reorders the matrix when using the div (/) |
720 and ldiv (\) operators, and so no the user does not need to explicitly | |
721 reorder the matrix to maximize performance. | |
722 | |
7619 | 723 @DOCSTRING(amd) |
724 | |
6620 | 725 @DOCSTRING(ccolamd) |
726 | |
727 @DOCSTRING(colamd) | |
728 | |
729 @DOCSTRING(colperm) | |
730 | |
731 @DOCSTRING(csymamd) | |
732 | |
733 @DOCSTRING(dmperm) | |
734 | |
735 @DOCSTRING(symamd) | |
736 | |
737 @DOCSTRING(symrcm) | |
738 | |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
739 @node Sparse Linear Algebra |
5164 | 740 @section Linear Algebra on Sparse Matrices |
741 | |
8488
cdb4788879b3
[docs] poly-morphic => polymorphic
Brian Gough <bjg@gnu.org>
parents:
8417
diff
changeset
|
742 Octave includes a polymorphic solver for sparse matrices, where |
5164 | 743 the exact solver used to factorize the matrix, depends on the properties |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
744 of the sparse matrix itself. Generally, the cost of determining the matrix type |
5322 | 745 is small relative to the cost of factorizing the matrix itself, but in any |
746 case the matrix type is cached once it is calculated, so that it is not | |
747 re-determined each time it is used in a linear equation. | |
5164 | 748 |
749 The selection tree for how the linear equation is solve is | |
750 | |
751 @enumerate 1 | |
5648 | 752 @item If the matrix is diagonal, solve directly and goto 8 |
5164 | 753 |
754 @item If the matrix is a permuted diagonal, solve directly taking into | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
755 account the permutations. Goto 8 |
5164 | 756 |
5648 | 757 @item If the matrix is square, banded and if the band density is less |
758 than that given by @code{spparms ("bandden")} continue, else goto 4. | |
5164 | 759 |
760 @enumerate a | |
761 @item If the matrix is tridiagonal and the right-hand side is not sparse | |
5648 | 762 continue, else goto 3b. |
5164 | 763 |
764 @enumerate | |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
765 @item If the matrix is Hermitian, with a positive real diagonal, attempt |
11593
1577c6f80926
Use non-breaking spaces between certain adjectives and their nouns in docstrings.
Rik <octave@nomad.inbox5.com>
parents:
11573
diff
changeset
|
766 Cholesky@tie{}factorization using @sc{lapack} xPTSV. |
5164 | 767 |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
768 @item If the above failed or the matrix is not Hermitian with a positive |
5164 | 769 real diagonal use Gaussian elimination with pivoting using |
9088
77e71f3da3d6
Fix documentation image printing under new development code
Rik <rdrider0-list@yahoo.com>
parents:
9066
diff
changeset
|
770 @sc{lapack} xGTSV, and goto 8. |
5164 | 771 @end enumerate |
772 | |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
773 @item If the matrix is Hermitian with a positive real diagonal, attempt |
11593
1577c6f80926
Use non-breaking spaces between certain adjectives and their nouns in docstrings.
Rik <octave@nomad.inbox5.com>
parents:
11573
diff
changeset
|
774 Cholesky@tie{}factorization using @sc{lapack} xPBTRF. |
5164 | 775 |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
776 @item if the above failed or the matrix is not Hermitian with a positive |
5164 | 777 real diagonal use Gaussian elimination with pivoting using |
9088
77e71f3da3d6
Fix documentation image printing under new development code
Rik <rdrider0-list@yahoo.com>
parents:
9066
diff
changeset
|
778 @sc{lapack} xGBTRF, and goto 8. |
5164 | 779 @end enumerate |
780 | |
781 @item If the matrix is upper or lower triangular perform a sparse forward | |
5648 | 782 or backward substitution, and goto 8 |
5164 | 783 |
11573
6f8ffe2c6f76
Grammarcheck txi files for 3.4 release.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
784 @item If the matrix is an upper triangular matrix with column permutations |
5322 | 785 or lower triangular matrix with row permutations, perform a sparse forward |
5648 | 786 or backward substitution, and goto 8 |
5164 | 787 |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
788 @item If the matrix is square, Hermitian with a real positive diagonal, attempt |
11593
1577c6f80926
Use non-breaking spaces between certain adjectives and their nouns in docstrings.
Rik <octave@nomad.inbox5.com>
parents:
11573
diff
changeset
|
789 sparse Cholesky@tie{}factorization using @sc{cholmod}. |
5164 | 790 |
11593
1577c6f80926
Use non-breaking spaces between certain adjectives and their nouns in docstrings.
Rik <octave@nomad.inbox5.com>
parents:
11573
diff
changeset
|
791 @item If the sparse Cholesky@tie{}factorization failed or the matrix is not |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
792 Hermitian with a real positive diagonal, and the matrix is square, factorize |
10711
fbd7843974fa
Periodic grammar check of documentation files to ensure common format.
Rik <octave@nomad.inbox5.com>
parents:
10668
diff
changeset
|
793 using @sc{umfpack}. |
5164 | 794 |
795 @item If the matrix is not square, or any of the previous solvers flags | |
5648 | 796 a singular or near singular matrix, find a minimum norm solution using |
10828
322f43e0e170
Grammarcheck .txi documentation files.
Rik <octave@nomad.inbox5.com>
parents:
10791
diff
changeset
|
797 @sc{cxsparse}@footnote{The @sc{cholmod}, @sc{umfpack} and @sc{cxsparse} packages were |
7096 | 798 written by Tim Davis and are available at |
16792
b157ba28f123
doc: Use Texinfo @url command consistently throughout documentation.
Rik <rik@octave.org>
parents:
16601
diff
changeset
|
799 @url{http://www.cise.ufl.edu/research/sparse/}}. |
5164 | 800 @end enumerate |
801 | |
16600
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
802 The band density is defined as the number of non-zero values in the band |
16816
12005245b645
doc: Periodic grammarcheck of documentation.
Rik <rik@octave.org>
parents:
16792
diff
changeset
|
803 divided by the total number of values in the full band. The banded |
16600
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
804 matrix solvers can be entirely disabled by using @dfn{spparms} to set |
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
805 @code{bandden} to 1 (i.e., @code{spparms ("bandden", 1)}). |
5164 | 806 |
16600
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
807 The QR@tie{}solver factorizes the problem with a Dulmage-Mendelsohn |
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
808 decomposition, to separate the problem into blocks that can be treated |
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
809 as over-determined, multiple well determined blocks, and a final |
16816
12005245b645
doc: Periodic grammarcheck of documentation.
Rik <rik@octave.org>
parents:
16792
diff
changeset
|
810 over-determined block. For matrices with blocks of strongly connected |
16600
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
811 nodes this is a big win as LU@tie{}decomposition can be used for many |
16816
12005245b645
doc: Periodic grammarcheck of documentation.
Rik <rik@octave.org>
parents:
16792
diff
changeset
|
812 blocks. It also significantly improves the chance of finding a solution |
16600
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
813 to over-determined problems rather than just returning a vector of |
f2f5dd09e97d
doc: fix some minor sparse documentation oversights
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
14138
diff
changeset
|
814 @dfn{NaN}'s. |
5681 | 815 |
816 All of the solvers above, can calculate an estimate of the condition | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
817 number. This can be used to detect numerical stability problems in the |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
818 solution and force a minimum norm solution to be used. However, for |
5681 | 819 narrow banded, triangular or diagonal matrices, the cost of |
820 calculating the condition number is significant, and can in fact | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
821 exceed the cost of factoring the matrix. Therefore the condition |
6939 | 822 number is not calculated in these cases, and Octave relies on simpler |
9088
77e71f3da3d6
Fix documentation image printing under new development code
Rik <rdrider0-list@yahoo.com>
parents:
9066
diff
changeset
|
823 techniques to detect singular matrices or the underlying @sc{lapack} code in |
5681 | 824 the case of banded matrices. |
5164 | 825 |
5322 | 826 The user can force the type of the matrix with the @code{matrix_type} |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
827 function. This overcomes the cost of discovering the type of the matrix. |
8347
fa78cb8d8a5c
corrections for typos
Brian Gough<bjg@network-theory.co.uk>
parents:
8286
diff
changeset
|
828 However, it should be noted that identifying the type of the matrix incorrectly |
5322 | 829 will lead to unpredictable results, and so @code{matrix_type} should be |
5506 | 830 used with care. |
5322 | 831 |
6620 | 832 @DOCSTRING(normest) |
833 | |
8286
6f2d95255911
fix @seealso references to point to existing anchors
Thorsten Meyer <thorsten.meyier@gmx.de>
parents:
8106
diff
changeset
|
834 @DOCSTRING(onenormest) |
6f2d95255911
fix @seealso references to point to existing anchors
Thorsten Meyer <thorsten.meyier@gmx.de>
parents:
8106
diff
changeset
|
835 |
7189 | 836 @DOCSTRING(condest) |
837 | |
6620 | 838 @DOCSTRING(spparms) |
839 | |
840 @DOCSTRING(sprank) | |
841 | |
842 @DOCSTRING(symbfact) | |
843 | |
7681
b1c1133641ee
Add the spaugment function
David Bateman <dbateman@free.fr>
parents:
7619
diff
changeset
|
844 For non square matrices, the user can also utilize the @code{spaugment} |
b1c1133641ee
Add the spaugment function
David Bateman <dbateman@free.fr>
parents:
7619
diff
changeset
|
845 function to find a least squares solution to a linear equation. |
b1c1133641ee
Add the spaugment function
David Bateman <dbateman@free.fr>
parents:
7619
diff
changeset
|
846 |
b1c1133641ee
Add the spaugment function
David Bateman <dbateman@free.fr>
parents:
7619
diff
changeset
|
847 @DOCSTRING(spaugment) |
b1c1133641ee
Add the spaugment function
David Bateman <dbateman@free.fr>
parents:
7619
diff
changeset
|
848 |
8417
654bcfb937bf
Add the eigs and svds functions
David Bateman <dbateman@free.fr>
parents:
8347
diff
changeset
|
849 Finally, the function @code{eigs} can be used to calculate a limited |
654bcfb937bf
Add the eigs and svds functions
David Bateman <dbateman@free.fr>
parents:
8347
diff
changeset
|
850 number of eigenvalues and eigenvectors based on a selection criteria |
654bcfb937bf
Add the eigs and svds functions
David Bateman <dbateman@free.fr>
parents:
8347
diff
changeset
|
851 and likewise for @code{svds} which calculates a limited number of |
654bcfb937bf
Add the eigs and svds functions
David Bateman <dbateman@free.fr>
parents:
8347
diff
changeset
|
852 singular values and vectors. |
654bcfb937bf
Add the eigs and svds functions
David Bateman <dbateman@free.fr>
parents:
8347
diff
changeset
|
853 |
654bcfb937bf
Add the eigs and svds functions
David Bateman <dbateman@free.fr>
parents:
8347
diff
changeset
|
854 @DOCSTRING(eigs) |
654bcfb937bf
Add the eigs and svds functions
David Bateman <dbateman@free.fr>
parents:
8347
diff
changeset
|
855 |
654bcfb937bf
Add the eigs and svds functions
David Bateman <dbateman@free.fr>
parents:
8347
diff
changeset
|
856 @DOCSTRING(svds) |
654bcfb937bf
Add the eigs and svds functions
David Bateman <dbateman@free.fr>
parents:
8347
diff
changeset
|
857 |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
858 @node Iterative Techniques |
5164 | 859 @section Iterative Techniques applied to sparse matrices |
860 | |
6620 | 861 The left division @code{\} and right division @code{/} operators, |
862 discussed in the previous section, use direct solvers to resolve a | |
863 linear equation of the form @code{@var{x} = @var{A} \ @var{b}} or | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
864 @code{@var{x} = @var{b} / @var{A}}. Octave equally includes a number of |
6620 | 865 functions to solve sparse linear equations using iterative techniques. |
866 | |
867 @DOCSTRING(pcg) | |
868 | |
869 @DOCSTRING(pcr) | |
5837 | 870 |
6620 | 871 The speed with which an iterative solver converges to a solution can be |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
872 accelerated with the use of a pre-conditioning matrix @var{M}. In this |
6620 | 873 case the linear equation @code{@var{M}^-1 * @var{x} = @var{M}^-1 * |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
874 @var{A} \ @var{b}} is solved instead. Typical pre-conditioning matrices |
6620 | 875 are partial factorizations of the original matrix. |
5648 | 876 |
6620 | 877 @DOCSTRING(luinc) |
878 | |
8817
03b7f618ab3d
include docstrings for new functions in the manual
John W. Eaton <jwe@octave.org>
parents:
8488
diff
changeset
|
879 @node Real Life Example |
5648 | 880 @section Real Life Example of the use of Sparse Matrices |
881 | |
882 A common application for sparse matrices is in the solution of Finite | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
883 Element Models. Finite element models allow numerical solution of |
5648 | 884 partial differential equations that do not have closed form solutions, |
885 typically because of the complex shape of the domain. | |
886 | |
887 In order to motivate this application, we consider the boundary value | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
888 Laplace equation. This system can model scalar potential fields, such |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
889 as heat or electrical potential. Given a medium |
5648 | 890 @tex |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
891 $\Omega$ with boundary $\partial\Omega$. At all points on the $\partial\Omega$ |
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
892 the boundary conditions are known, and we wish to calculate the potential in |
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
893 $\Omega$. |
5648 | 894 @end tex |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
895 @ifnottex |
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
896 Omega with boundary dOmega. At all points on the dOmega |
5648 | 897 the boundary conditions are known, and we wish to calculate the potential in |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
898 Omega. |
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
899 @end ifnottex |
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
900 Boundary conditions may specify the potential (Dirichlet |
5648 | 901 boundary condition), its normal derivative across the boundary |
902 (Neumann boundary condition), or a weighted sum of the potential and | |
903 its derivative (Cauchy boundary condition). | |
904 | |
905 In a thermal model, we want to calculate the temperature in | |
906 @tex | |
907 $\Omega$ | |
908 @end tex | |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
909 @ifnottex |
5648 | 910 Omega |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
911 @end ifnottex |
5648 | 912 and know the boundary temperature (Dirichlet condition) |
913 or heat flux (from which we can calculate the Neumann condition | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
914 by dividing by the thermal conductivity at the boundary). Similarly, |
5648 | 915 in an electrical model, we want to calculate the voltage in |
916 @tex | |
917 $\Omega$ | |
918 @end tex | |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
919 @ifnottex |
5648 | 920 Omega |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
921 @end ifnottex |
5648 | 922 and know the boundary voltage (Dirichlet) or current |
923 (Neumann condition after diving by the electrical conductivity). | |
924 In an electrical model, it is common for much of the boundary | |
925 to be electrically isolated; this is a Neumann boundary condition | |
926 with the current equal to zero. | |
927 | |
928 The simplest finite element models will divide | |
929 @tex | |
930 $\Omega$ | |
931 @end tex | |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
932 @ifnottex |
5648 | 933 Omega |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
934 @end ifnottex |
5648 | 935 into simplexes (triangles in 2D, pyramids in 3D). |
936 @ifset htmltex | |
11573
6f8ffe2c6f76
Grammarcheck txi files for 3.4 release.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
937 We take as a 3-D example a cylindrical liquid filled tank with a small |
5648 | 938 non-conductive ball from the EIDORS project@footnote{EIDORS - Electrical |
939 Impedance Tomography and Diffuse optical Tomography Reconstruction Software | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
940 @url{http://eidors3d.sourceforge.net}}. This is model is designed to reflect |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
941 an application of electrical impedance tomography, where current patterns |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
942 are applied to such a tank in order to image the internal conductivity |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
943 distribution. In order to describe the FEM geometry, we have a matrix of |
5648 | 944 vertices @code{nodes} and simplices @code{elems}. |
945 @end ifset | |
946 | |
11573
6f8ffe2c6f76
Grammarcheck txi files for 3.4 release.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
947 The following example creates a simple rectangular 2-D electrically |
5648 | 948 conductive medium with 10 V and 20 V imposed on opposite sides |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
949 (Dirichlet boundary conditions). All other edges are electrically |
5648 | 950 isolated. |
951 | |
952 @example | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
953 @group |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
954 node_y = [1;1.2;1.5;1.8;2]*ones(1,11); |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
955 node_x = ones(5,1)*[1,1.05,1.1,1.2, ... |
5648 | 956 1.3,1.5,1.7,1.8,1.9,1.95,2]; |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
957 nodes = [node_x(:), node_y(:)]; |
5648 | 958 |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
959 [h,w] = size (node_x); |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
960 elems = []; |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
961 for idx = 1:w-1 |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
962 widx = (idx-1)*h; |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
963 elems = [elems; ... |
5648 | 964 widx+[(1:h-1);(2:h);h+(1:h-1)]'; ... |
965 widx+[(2:h);h+(2:h);h+(1:h-1)]' ]; | |
966 endfor | |
967 | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
968 E = size (elems,1); # No. of simplices |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
969 N = size (nodes,1); # No. of vertices |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
970 D = size (elems,2); # dimensions+1 |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
971 @end group |
5648 | 972 @end example |
973 | |
974 This creates a N-by-2 matrix @code{nodes} and a E-by-3 matrix | |
975 @code{elems} with values, which define finite element triangles: | |
5164 | 976 |
5648 | 977 @example |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
978 @group |
5648 | 979 nodes(1:7,:)' |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
980 1.00 1.00 1.00 1.00 1.00 1.05 1.05 @dots{} |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
981 1.00 1.20 1.50 1.80 2.00 1.00 1.20 @dots{} |
5648 | 982 |
983 elems(1:7,:)' | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
984 1 2 3 4 2 3 4 @dots{} |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
985 2 3 4 5 7 8 9 @dots{} |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
986 6 7 8 9 6 7 8 @dots{} |
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
987 @end group |
5648 | 988 @end example |
989 | |
990 Using a first order FEM, we approximate the electrical conductivity | |
991 distribution in | |
992 @tex | |
993 $\Omega$ | |
994 @end tex | |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
995 @ifnottex |
5648 | 996 Omega |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
997 @end ifnottex |
5648 | 998 as constant on each simplex (represented by the vector @code{conductivity}). |
999 Based on the finite element geometry, we first calculate a system (or | |
1000 stiffness) matrix for each simplex (represented as 3-by-3 elements on the | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
1001 diagonal of the element-wise system matrix @code{SE}. Based on @code{SE} |
5648 | 1002 and a N-by-DE connectivity matrix @code{C}, representing the connections |
7001 | 1003 between simplices and vertices, the global connectivity matrix @code{S} is |
5648 | 1004 calculated. |
1005 | |
1006 @example | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1007 ## Element conductivity |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1008 conductivity = [1*ones(1,16), ... |
5648 | 1009 2*ones(1,48), 1*ones(1,16)]; |
1010 | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1011 ## Connectivity matrix |
5648 | 1012 C = sparse ((1:D*E), reshape (elems', ... |
1013 D*E, 1), 1, D*E, N); | |
1014 | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1015 ## Calculate system matrix |
5648 | 1016 Siidx = floor ([0:D*E-1]'/D) * D * ... |
1017 ones(1,D) + ones(D*E,1)*(1:D) ; | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1018 Sjidx = [1:D*E]'*ones (1,D); |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1019 Sdata = zeros (D*E,D); |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1020 dfact = factorial (D-1); |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1021 for j = 1:E |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1022 a = inv ([ones(D,1), ... |
5648 | 1023 nodes(elems(j,:), :)]); |
1024 const = conductivity(j) * 2 / ... | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1025 dfact / abs (det (a)); |
5648 | 1026 Sdata(D*(j-1)+(1:D),:) = const * ... |
1027 a(2:D,:)' * a(2:D,:); | |
1028 endfor | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1029 ## Element-wise system matrix |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1030 SE = sparse(Siidx,Sjidx,Sdata); |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1031 ## Global system matrix |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1032 S = C'* SE *C; |
5648 | 1033 @end example |
1034 | |
1035 The system matrix acts like the conductivity | |
1036 @tex | |
1037 $S$ | |
1038 @end tex | |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
1039 @ifnottex |
5648 | 1040 @code{S} |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
1041 @end ifnottex |
5648 | 1042 in Ohm's law |
1043 @tex | |
1044 $SV = I$. | |
1045 @end tex | |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
1046 @ifnottex |
5648 | 1047 @code{S * V = I}. |
10668
72585f1ca7a2
Replace @ifinfo with @ifnottex.
Rik <octave@nomad.inbox5.com>
parents:
10599
diff
changeset
|
1048 @end ifnottex |
5648 | 1049 Based on the Dirichlet and Neumann boundary conditions, we are able to |
1050 solve for the voltages at each vertex @code{V}. | |
1051 | |
1052 @example | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1053 ## Dirichlet boundary conditions |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1054 D_nodes = [1:5, 51:55]; |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1055 D_value = [10*ones(1,5), 20*ones(1,5)]; |
5648 | 1056 |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1057 V = zeros (N,1); |
5648 | 1058 V(D_nodes) = D_value; |
1059 idx = 1:N; # vertices without Dirichlet | |
1060 # boundary condns | |
1061 idx(D_nodes) = []; | |
1062 | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1063 ## Neumann boundary conditions. Note that |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1064 ## N_value must be normalized by the |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1065 ## boundary length and element conductivity |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1066 N_nodes = []; |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1067 N_value = []; |
5648 | 1068 |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1069 Q = zeros (N,1); |
5648 | 1070 Q(N_nodes) = N_value; |
1071 | |
1072 V(idx) = S(idx,idx) \ ( Q(idx) - ... | |
1073 S(idx,D_nodes) * V(D_nodes)); | |
1074 @end example | |
1075 | |
1076 Finally, in order to display the solution, we show each solved voltage | |
1077 value in the z-axis for each simplex vertex. | |
1078 @ifset htmltex | |
1079 @xref{fig:femmodel}. | |
1080 @end ifset | |
1081 | |
1082 @example | |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
1083 @group |
5648 | 1084 elemx = elems(:,[1,2,3,1])'; |
1085 xelems = reshape (nodes(elemx, 1), 4, E); | |
1086 yelems = reshape (nodes(elemx, 2), 4, E); | |
1087 velems = reshape (V(elemx), 4, E); | |
14856
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1088 plot3 (xelems,yelems,velems,"k"); |
c3fd61c59e9c
maint: Use Octave coding conventions for cuddling parentheses in doc directory
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
1089 print "grid.eps"; |
9066
be150a172010
Cleanup documentation for diagperm.texi, sparse.texi
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
1090 @end group |
5648 | 1091 @end example |
1092 | |
1093 | |
1094 @ifset htmltex | |
1095 @float Figure,fig:femmodel | |
9088
77e71f3da3d6
Fix documentation image printing under new development code
Rik <rdrider0-list@yahoo.com>
parents:
9066
diff
changeset
|
1096 @center @image{grid,4in} |
5648 | 1097 @caption{Example finite element model the showing triangular elements. |
1098 The height of each vertex corresponds to the solution value.} | |
1099 @end float | |
1100 @end ifset |