--- a/doc/interpreter/sparse.txi	Thu Jan 01 18:47:36 2015 -0800
+++ b/doc/interpreter/sparse.txi	Tue Jan 20 08:26:57 2015 -0500
@@ -6,12 +6,12 @@
 @c under the terms of the GNU General Public License as published by the
 @c Free Software Foundation; either version 3 of the License, or (at
 @c your option) any later version.
-@c 
+@c
 @c Octave is distributed in the hope that it will be useful, but WITHOUT
 @c ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
 @c FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
 @c for more details.
-@c 
+@c
 @c You should have received a copy of the GNU General Public License
 @c along with Octave; see the file COPYING.  If not, see
 @c <http://www.gnu.org/licenses/>.
@@ -39,7 +39,7 @@
 The size of mathematical problems that can be treated at any particular
 time is generally limited by the available computing resources.  Both,
 the speed of the computer and its available memory place limitation on
-the problem size. 
+the problem size.
 
 There are many classes of mathematical problems which give rise to
 matrices, where a large number of the elements are zero.  In this case
@@ -68,8 +68,8 @@
 It is not strictly speaking necessary for the user to understand how
 sparse matrices are stored.  However, such an understanding will help
 to get an understanding of the size of sparse matrices.  Understanding
-the storage technique is also necessary for those users wishing to 
-create their own oct-files. 
+the storage technique is also necessary for those users wishing to
+create their own oct-files.
 
 There are many different means of storing sparse matrix data.  What all
 of the methods have in common is that they attempt to reduce the complexity
@@ -79,12 +79,12 @@
 for sparse matrix computation", 1994,
 @url{http://www-users.cs.umn.edu/~saad/software/SPARSKIT/paper.ps}}.
 With full matrices, knowledge of the point of an element of the matrix
-within the matrix is implied by its position in the computers memory. 
+within the matrix is implied by its position in the computers memory.
 However, this is not the case for sparse matrices, and so the positions
-of the non-zero elements of the matrix must equally be stored. 
+of the non-zero elements of the matrix must equally be stored.
 
 An obvious way to do this is by storing the elements of the matrix as
-triplets, with two elements being their position in the array 
+triplets, with two elements being their position in the array
 (rows and column) and the third being the data itself.  This is conceptually
 easy to grasp, but requires more storage than is strictly needed.
 
@@ -109,7 +109,7 @@
 @group
   for (j = 0; j < nc; j++)
     for (i = cidx(j); i < cidx(j+1); i++)
-       printf ("non-zero element (%i,%i) is %d\n", 
+       printf ("non-zero element (%i,%i) is %d\n",
            ridx(i), j, data(i));
 @end group
 @end example
@@ -149,9 +149,9 @@
 @end example
 
 Note that this is the representation of these elements with the first row
-and column assumed to start at zero, while in Octave itself the row and 
-column indexing starts at one.  Thus the number of elements in the 
-@var{i}-th column is given by @code{@var{cidx} (@var{i} + 1) - 
+and column assumed to start at zero, while in Octave itself the row and
+column indexing starts at one.  Thus the number of elements in the
+@var{i}-th column is given by @code{@var{cidx} (@var{i} + 1) -
 @var{cidx} (@var{i})}.
 
 Although Octave uses a compressed column format, it should be noted
@@ -162,7 +162,7 @@
 column major ordering, and so sparse matrices are equally stored in
 this manner.
 
-A further constraint on the sparse matrix storage used by Octave is that 
+A further constraint on the sparse matrix storage used by Octave is that
 all elements in the rows are stored in increasing order of their row
 index, which makes certain operations faster.  However, it imposes
 the need to sort the elements on the creation of sparse matrices.  Having
@@ -181,7 +181,7 @@
 @dfn{speye}, @dfn{sprand}, @dfn{diag}, etc.
 
 @item Constructed from matrices or vectors
-The function @dfn{sparse} allows a sparse matrix to be constructed from 
+The function @dfn{sparse} allows a sparse matrix to be constructed from
 three vectors representing the row, column and data.  Alternatively, the
 function @dfn{spconvert} uses a three column matrix format to allow easy
 importation of data from elsewhere.
@@ -210,7 +210,7 @@
 
 Other functions of interest that directly create sparse matrices, are
 @dfn{diag} or its generalization @dfn{spdiags}, that can take the
-definition of the diagonals of the matrix and create the sparse matrix 
+definition of the diagonals of the matrix and create the sparse matrix
 that corresponds to this.  For example,
 
 @example
@@ -233,7 +233,7 @@
 
 @DOCSTRING(sprandsym)
 
-The recommended way for the user to create a sparse matrix, is to create 
+The recommended way for the user to create a sparse matrix, is to create
 two vectors containing the row and column index of the data and a third
 vector of the same size containing the data to be stored.  For example,
 
@@ -259,7 +259,7 @@
 The function @dfn{spconvert} takes a three or four column real matrix.
 The first two columns represent the row and column index respectively and
 the third and four columns, the real and imaginary parts of the sparse
-matrix.  The matrix can contain zero elements and the elements can be 
+matrix.  The matrix can contain zero elements and the elements can be
 sorted in any order.  Adding zero elements is a convenient way to define
 the size of the sparse matrix.  For example:
 
@@ -290,8 +290,8 @@
 
 It should be noted, that due to the way that the Octave
 assignment functions are written that the assignment will reallocate
-the memory used by the sparse matrix at each iteration of the above loop. 
-Therefore the @dfn{spalloc} function ignores the @var{nz} argument and 
+the memory used by the sparse matrix at each iteration of the above loop.
+Therefore the @dfn{spalloc} function ignores the @var{nz} argument and
 does not pre-assign the memory for the matrix.  Therefore, it is vitally
 important that code using to above structure should be vectorized
 as much as possible to minimize the number of assignments and reduce the
@@ -350,7 +350,7 @@
 @example
 @group
 a = tril (sprandn (1024, 1024, 0.02), -1) ...
-    + speye (1024); 
+    + speye (1024);
 matrix_type (a);
 ans = Lower
 @end group
@@ -410,7 +410,7 @@
 which creates an adjacency matrix @code{A} where node 1 is connected
 to nodes 2 and 6, node 2 with nodes 1 and 3, etc.  The coordinates of
 the nodes are given in the n-by-2 matrix @code{xy}.
-@ifset htmltex 
+@ifset htmltex
 @xref{fig:gplot}.
 
 @float Figure,fig:gplot
@@ -454,7 +454,7 @@
 matrices.  There is no difference in calling convention when using an
 overloaded function with a sparse matrix, however, there is also no access to
 potentially sparse-specific features.  At any time the sparse matrix specific
-version of a function can be used by explicitly calling its function name. 
+version of a function can be used by explicitly calling its function name.
 
 The table below lists all of the sparse functions of Octave.  Note that the
 names of the specific sparse forms of the functions are typically the same as
@@ -465,7 +465,7 @@
 
 @table @asis
 @item Generate sparse matrices:
-  @dfn{spalloc}, @dfn{spdiags}, @dfn{speye}, @dfn{sprand}, 
+  @dfn{spalloc}, @dfn{spdiags}, @dfn{speye}, @dfn{sprand},
   @dfn{sprandn}, @dfn{sprandsym}
 
 @item Sparse matrix conversion:
@@ -476,7 +476,7 @@
   @dfn{spfun}, @dfn{spones}, @dfn{spy}
 
 @item Graph Theory:
-  @dfn{etree}, @dfn{etreeplot}, @dfn{gplot}, 
+  @dfn{etree}, @dfn{etreeplot}, @dfn{gplot},
   @dfn{treeplot}
 @c @dfn{treelayout}
 
@@ -490,7 +490,7 @@
 
 @item Iterative techniques:
   @dfn{luinc}, @dfn{pcg}, @dfn{pcr}
-@c @dfn{bicg}, @dfn{bicgstab}, @dfn{cholinc}, @dfn{cgs}, @dfn{gmres}, 
+@c @dfn{bicg}, @dfn{bicgstab}, @dfn{cholinc}, @dfn{cgs}, @dfn{gmres},
 @c @dfn{lsqr}, @dfn{minres}, @dfn{qmr}, @dfn{symmlq}
 
 @item Miscellaneous:
@@ -506,14 +506,14 @@
 @node Return Types of Operators and Functions
 @subsubsection Return Types of Operators and Functions
 
-The two basic reasons to use sparse matrices are to reduce the memory 
+The two basic reasons to use sparse matrices are to reduce the memory
 usage and to not have to do calculations on zero elements.  The two are
 closely related in that the computation time on a sparse matrix operator
 or function is roughly linear with the number of non-zero elements.
 
-Therefore, there is a certain density of non-zero elements of a matrix 
+Therefore, there is a certain density of non-zero elements of a matrix
 where it no longer makes sense to store it as a sparse matrix, but rather
-as a full matrix.  For this reason operators and functions that have a 
+as a full matrix.  For this reason operators and functions that have a
 high probability of returning a full matrix will always return one.  For
 example adding a scalar constant to a sparse matrix will almost always
 make it a full matrix, and so the example,
@@ -528,7 +528,7 @@
 @end example
 
 @noindent
-returns a full matrix as can be seen. 
+returns a full matrix as can be seen.
 
 
 Additionally, if @code{sparse_auto_mutate} is true, all sparse functions
@@ -537,19 +537,19 @@
 equivalent.  Therefore @code{speye (2) * 1} will return a full matrix as
 the memory used is smaller for the full version than the sparse version.
 
-As all of the mixed operators and functions between full and sparse 
+As all of the mixed operators and functions between full and sparse
 matrices exist, in general this does not cause any problems.  However,
 one area where it does cause a problem is where a sparse matrix is
 promoted to a full matrix, where subsequent operations would resparsify
 the matrix.  Such cases are rare, but can be artificially created, for
 example @code{(fliplr (speye (3)) + speye (3)) - speye (3)} gives a full
-matrix when it should give a sparse one.  In general, where such cases 
+matrix when it should give a sparse one.  In general, where such cases
 occur, they impose only a small memory penalty.
 
 There is however one known case where this behavior of Octave's
 sparse matrices will cause a problem.  That is in the handling of the
 @dfn{diag} function.  Whether @dfn{diag} returns a sparse or full matrix
-depending on the type of its input arguments.  So 
+depending on the type of its input arguments.  So
 
 @example
  a = diag (sparse ([1,2,3]), -1);
@@ -557,8 +557,8 @@
 
 @noindent
 should return a sparse matrix.  To ensure this actually happens, the
-@dfn{sparse} function, and other functions based on it like @dfn{speye}, 
-always returns a sparse matrix, even if the memory used will be larger 
+@dfn{sparse} function, and other functions based on it like @dfn{speye},
+always returns a sparse matrix, even if the memory used will be larger
 than its full representation.
 
 @DOCSTRING(sparse_auto_mutate)
@@ -573,7 +573,7 @@
 same manner as there full counterparts.  However, there are certain differences
 and especially differences with other products sparse implementations.
 
-First, the @qcode{"./"} and @qcode{".^"} operators must be used with care. 
+First, the @qcode{"./"} and @qcode{".^"} operators must be used with care.
 Consider what the examples
 
 @example
@@ -592,19 +592,19 @@
 will give.  The first example of @var{s} raised to the power of 2 causes
 no problems.  However @var{s} raised element-wise to itself involves a
 large number of terms @code{0 .^ 0} which is 1. There @code{@var{s} .^
-@var{s}} is a full matrix. 
+@var{s}} is a full matrix.
 
 Likewise @code{@var{s} .^ -2} involves terms like @code{0 .^ -2} which
 is infinity, and so @code{@var{s} .^ -2} is equally a full matrix.
 
-For the "./" operator @code{@var{s} ./ 2} has no problems, but 
+For the "./" operator @code{@var{s} ./ 2} has no problems, but
 @code{2 ./ @var{s}} involves a large number of infinity terms as well
 and is equally a full matrix.  The case of @code{@var{s} ./ @var{s}}
 involves terms like @code{0 ./ 0} which is a @code{NaN} and so this
 is equally a full matrix with the zero elements of @var{s} filled with
 @code{NaN} values.
 
-The above behavior is consistent with full matrices, but is not 
+The above behavior is consistent with full matrices, but is not
 consistent with sparse implementations in other products.
 
 A particular problem of sparse matrices comes about due to the fact that
@@ -622,9 +622,9 @@
      Inf            Inf
 @end group
 @end example
- 
+
 To correct this behavior would mean that zero elements with a negative
-sign-bit would need to be stored in the matrix to ensure that their 
+sign-bit would need to be stored in the matrix to ensure that their
 sign-bit was respected.  This is not done at this time, for reasons of
 efficiency, and so the user is warned that calculations where the sign-bit
 of zero is important must not be done using sparse matrices.
@@ -645,7 +645,7 @@
 the reordering functions @dfn{colperm} and @dfn{randperm} are
 also available.
 
-@xref{fig:simplematrix}, for an example of the structure of a simple 
+@xref{fig:simplematrix}, for an example of the structure of a simple
 positive definite matrix.
 
 @float Figure,fig:simplematrix
@@ -655,10 +655,10 @@
 
 The standard Cholesky@tie{}factorization of this matrix can be
 obtained by the same command that would be used for a full
-matrix.  This can be visualized with the command 
+matrix.  This can be visualized with the command
 @code{r = chol (A); spy (r);}.
 @xref{fig:simplechol}.
-The original matrix had 
+The original matrix had
 @ifinfo
 @ifnothtml
 43
@@ -678,13 +678,13 @@
 @end ifset
 with only half of the symmetric matrix being stored.  This
 is a significant level of fill in, and although not an issue
-for such a small test case, can represents a large overhead 
+for such a small test case, can represents a large overhead
 in working with other sparse matrices.
 
 The appropriate sparsity preserving permutation of the original
 matrix is given by @dfn{symamd} and the factorization using this
 reordering can be visualized using the command @code{q = symamd (A);
-r = chol (A(q,q)); spy (r)}.  This gives 
+r = chol (A(q,q)); spy (r)}.  This gives
 @ifinfo
 @ifnothtml
 29
@@ -738,7 +738,7 @@
 @node Sparse Linear Algebra
 @section Linear Algebra on Sparse Matrices
 
-Octave includes a polymorphic solver for sparse matrices, where 
+Octave includes a polymorphic solver for sparse matrices, where
 the exact solver used to factorize the matrix, depends on the properties
 of the sparse matrix itself.  Generally, the cost of determining the matrix type
 is small relative to the cost of factorizing the matrix itself, but in any
@@ -757,7 +757,7 @@
 than that given by @code{spparms ("bandden")} continue, else goto 4.
 
 @enumerate a
-@item If the matrix is tridiagonal and the right-hand side is not sparse 
+@item If the matrix is tridiagonal and the right-hand side is not sparse
 continue, else goto 3b.
 
 @enumerate
@@ -765,7 +765,7 @@
       Cholesky@tie{}factorization using @sc{lapack} xPTSV.
 
 @item If the above failed or the matrix is not Hermitian with a positive
-      real diagonal use Gaussian elimination with pivoting using 
+      real diagonal use Gaussian elimination with pivoting using
       @sc{lapack} xGTSV, and goto 8.
 @end enumerate
 
@@ -773,7 +773,7 @@
       Cholesky@tie{}factorization using @sc{lapack} xPBTRF.
 
 @item if the above failed or the matrix is not Hermitian with a positive
-      real diagonal use Gaussian elimination with pivoting using 
+      real diagonal use Gaussian elimination with pivoting using
       @sc{lapack} xGBTRF, and goto 8.
 @end enumerate
 
@@ -781,14 +781,14 @@
 or backward substitution, and goto 8
 
 @item If the matrix is an upper triangular matrix with column permutations
-or lower triangular matrix with row permutations, perform a sparse forward 
+or lower triangular matrix with row permutations, perform a sparse forward
 or backward substitution, and goto 8
 
 @item If the matrix is square, Hermitian with a real positive diagonal, attempt
 sparse Cholesky@tie{}factorization using @sc{cholmod}.
 
 @item If the sparse Cholesky@tie{}factorization failed or the matrix is not
-Hermitian with a real positive diagonal, and the matrix is square, factorize 
+Hermitian with a real positive diagonal, and the matrix is square, factorize
 using @sc{umfpack}.
 
 @item If the matrix is not square, or any of the previous solvers flags
@@ -885,7 +885,7 @@
 
 In order to motivate this application, we consider the boundary value
 Laplace equation.  This system can model scalar potential fields, such
-as heat or electrical potential.  Given a medium 
+as heat or electrical potential.  Given a medium
 @tex
 $\Omega$ with boundary $\partial\Omega$.  At all points on the $\partial\Omega$
 the boundary conditions are known, and we wish to calculate the potential in
@@ -903,17 +903,17 @@
 
 In a thermal model, we want to calculate the temperature in
 @tex
-$\Omega$ 
+$\Omega$
 @end tex
 @ifnottex
 Omega
 @end ifnottex
 and know the boundary temperature (Dirichlet condition)
 or heat flux (from which we can calculate the Neumann condition
-by dividing by the thermal conductivity at the boundary).  Similarly, 
+by dividing by the thermal conductivity at the boundary).  Similarly,
 in an electrical model, we want to calculate the voltage in
 @tex
-$\Omega$ 
+$\Omega$
 @end tex
 @ifnottex
 Omega
@@ -924,27 +924,27 @@
 to be electrically isolated; this is a Neumann boundary condition
 with the current equal to zero.
 
-The simplest finite element models will divide 
+The simplest finite element models will divide
 @tex
-$\Omega$ 
+$\Omega$
 @end tex
 @ifnottex
 Omega
 @end ifnottex
 into simplexes (triangles in 2D, pyramids in 3D).
 @ifset htmltex
-We take as a 3-D example a cylindrical liquid filled tank with a small 
-non-conductive ball from the EIDORS project@footnote{EIDORS - Electrical 
-Impedance Tomography and Diffuse optical Tomography Reconstruction Software 
+We take as a 3-D example a cylindrical liquid filled tank with a small
+non-conductive ball from the EIDORS project@footnote{EIDORS - Electrical
+Impedance Tomography and Diffuse optical Tomography Reconstruction Software
 @url{http://eidors3d.sourceforge.net}}.  This is model is designed to reflect
 an application of electrical impedance tomography, where current patterns
 are applied to such a tank in order to image the internal conductivity
-distribution.  In order to describe the FEM geometry, we have a matrix of 
+distribution.  In order to describe the FEM geometry, we have a matrix of
 vertices @code{nodes} and simplices @code{elems}.
 @end ifset
 
 The following example creates a simple rectangular 2-D electrically
-conductive medium with 10 V and 20 V imposed on opposite sides 
+conductive medium with 10 V and 20 V imposed on opposite sides
 (Dirichlet boundary conditions).  All other edges are electrically
 isolated.
 
@@ -961,7 +961,7 @@
      widx = (idx-1)*h;
      elems = [elems; ...
        widx+[(1:h-1);(2:h);h+(1:h-1)]'; ...
-       widx+[(2:h);h+(2:h);h+(1:h-1)]' ]; 
+       widx+[(2:h);h+(2:h);h+(1:h-1)]' ];
    endfor
 
    E = size (elems,1); # No. of simplices
@@ -986,10 +986,10 @@
 @end group
 @end example
 
-Using a first order FEM, we approximate the electrical conductivity 
-distribution in 
+Using a first order FEM, we approximate the electrical conductivity
+distribution in
 @tex
-$\Omega$ 
+$\Omega$
 @end tex
 @ifnottex
 Omega
@@ -997,8 +997,8 @@
 as constant on each simplex (represented by the vector @code{conductivity}).
 Based on the finite element geometry, we first calculate a system (or
 stiffness) matrix for each simplex (represented as 3-by-3 elements on the
-diagonal of the element-wise system matrix @code{SE}.  Based on @code{SE} 
-and a N-by-DE connectivity matrix @code{C}, representing the connections 
+diagonal of the element-wise system matrix @code{SE}.  Based on @code{SE}
+and a N-by-DE connectivity matrix @code{C}, representing the connections
 between simplices and vertices, the global connectivity matrix @code{S} is
 calculated.
 
@@ -1018,7 +1018,7 @@
   Sdata = zeros (D*E,D);
   dfact = factorial (D-1);
   for j = 1:E
-     a = inv ([ones(D,1), ... 
+     a = inv ([ones(D,1), ...
          nodes(elems(j,:), :)]);
      const = conductivity(j) * 2 / ...
          dfact / abs (det (a));
@@ -1031,31 +1031,31 @@
   S = C'* SE *C;
 @end example
 
-The system matrix acts like the conductivity 
+The system matrix acts like the conductivity
 @tex
-$S$ 
+$S$
 @end tex
 @ifnottex
 @code{S}
 @end ifnottex
-in Ohm's law 
+in Ohm's law
 @tex
-$SV = I$. 
+$SV = I$.
 @end tex
 @ifnottex
 @code{S * V = I}.
 @end ifnottex
-Based on the Dirichlet and Neumann boundary conditions, we are able to 
-solve for the voltages at each vertex @code{V}. 
+Based on the Dirichlet and Neumann boundary conditions, we are able to
+solve for the voltages at each vertex @code{V}.
 
 @example
   ## Dirichlet boundary conditions
-  D_nodes = [1:5, 51:55]; 
-  D_value = [10*ones(1,5), 20*ones(1,5)]; 
+  D_nodes = [1:5, 51:55];
+  D_value = [10*ones(1,5), 20*ones(1,5)];
 
   V = zeros (N,1);
   V(D_nodes) = D_value;
-  idx = 1:N; # vertices without Dirichlet 
+  idx = 1:N; # vertices without Dirichlet
              # boundary condns
   idx(D_nodes) = [];
 
@@ -1072,7 +1072,7 @@
             S(idx,D_nodes) * V(D_nodes));
 @end example
 
-Finally, in order to display the solution, we show each solved voltage 
+Finally, in order to display the solution, we show each solved voltage
 value in the z-axis for each simplex vertex.
 @ifset htmltex
 @xref{fig:femmodel}.
@@ -1084,7 +1084,7 @@
   xelems = reshape (nodes(elemx, 1), 4, E);
   yelems = reshape (nodes(elemx, 2), 4, E);
   velems = reshape (V(elemx), 4, E);
-  plot3 (xelems,yelems,velems,"k"); 
+  plot3 (xelems,yelems,velems,"k");
   print "grid.eps";
 @end group
 @end example
@@ -1093,7 +1093,7 @@
 @ifset htmltex
 @float Figure,fig:femmodel
 @center @image{grid,4in}
-@caption{Example finite element model the showing triangular elements. 
+@caption{Example finite element model the showing triangular elements.
 The height of each vertex corresponds to the solution value.}
 @end float
 @end ifset