Additional FAQ on the compatibility of mldivide for singular, under- and over-determined matrices

[David Bateman]

The attached changeset adds an FAQ based on Jaorslav's e-mail from a few 
months back on the compatibility of mldivide between Octave and Matlab. 
This question has come up a few times and so pointing the users to a 
pre-packaged response seems to me to be a good idea. It also adds a NEWS 
item for Jaroslav's other recent change.

D.

--
David Bateman                                address@hidden
Motorola Labs - Paris                        +33 1 69 35 48 04 (Ph) 
Parc Les Algorithmes, Commune de St Aubin    +33 6 72 01 06 33 (Mob) 
91193 Gif-Sur-Yvette FRANCE                  +33 1 69 35 77 01 (Fax) 

The information contained in this communication has been classified as: 

[x] General Business Information 
[ ] Motorola Internal Use Only 
[ ] Motorola Confidential Proprietary

# HG changeset patch
# User David Bateman <address@hidden>
# Date 1225375730 -3600
# Node ID 6d6d60e970a998dd03d8d4d99f7322e20ad8ed4c
# Parent  cc6c38eb950e6d97d15b68b3c97989603b5ff19c
Minor NEWS/FAQ update

diff --git a/ChangeLog b/ChangeLog
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,7 @@
+2008-10-30  David Bateman  <address@hidden>
+
+       * NEWS: Minor update to document improved indexing code.
+
 2008-10-23  Jason Riedy  <address@hidden>
 
        * configure.in: Take care to use FT2_CFLAGS when checking
diff --git a/NEWS b/NEWS
--- a/NEWS
+++ b/NEWS
@@ -105,4 +105,8 @@
     so users that make use of data generated by Octave with R or
     visa-versa are warned that compatibility might not be assured.
 
+ ** Improved array indexing
+    The underlying code used for indexing of arrays has been completely
+    rewritten and so the indexing of arrays is now significantly faster.
+
 See NEWS.3 for old news.
diff --git a/doc/ChangeLog b/doc/ChangeLog
--- a/doc/ChangeLog
+++ b/doc/ChangeLog
@@ -1,5 +1,9 @@
 2008-10-30  David Bateman  <address@hidden>
 
+       * faq/Octave-FAQ.texi: Document improved indexing and add an faq for
+       the compatibility of mldivide/mrdivide for singualr, under- and 
+       over-determined matrices.
+       
        * interpreter/plot.txi: Document contour groups.
 
 2008-10-29  Thorsten Meyer  <address@hidden>
diff --git a/doc/faq/Octave-FAQ.texi b/doc/faq/Octave-FAQ.texi
--- a/doc/faq/Octave-FAQ.texi
+++ b/doc/faq/Octave-FAQ.texi
@@ -251,6 +251,10 @@
 underlying LAPACK code.
 
 @item Single precision type
+
address@hidden Improved array indexing
+The underlying code used for indexing of arrays has been completely
+rewritten and so the indexing of arrays is now significantly faster.
 @end itemize
 
 
@@ -861,7 +865,7 @@
 @item Block comments
 Block comments denoted by "address@hidden" and "address@hidden" markers are 
supported by
 Octave with some limitations. The major limitation is that block
-comments are not supported with [] or @address@hidden
+comments are not supported within [] or @address@hidden
 
 @item Mat-File format
 There are some differences in the mat v5 file format accepted by
@@ -981,6 +985,80 @@
 any elements are nonzero. Octave however duplicates this behavior for if
 statements containing empty matrices.
 
address@hidden Solvers for singular, under- and over-determined matrices
+
+Matlab's solvers as used by the operators mldivide (\) and mrdivide (/),
+use a different approach than Octave's in the case of singular, under-, 
+or over-determined matrices. In the case of a singular matrix, Matlab
+returns the result given by the LU decomposition, even though the underlying
+solver has flagged the result as erroneous. Octave has made the choice
+of falling back to a minimum norm solution of matrices that have been
+flagged as singular which arguably is a better result for these cases.
+
+In the case of under- or over-determined matrices, Octave continues to
+use a minimum norm solution, whereas Matlab uses an approach that is
+equivalent to
+
address@hidden
address@hidden
+function x = mldivide (A, b)
+  [Q, R, E] = qr(A);
+  x = [A \ b, E(:, 1:m) * (R(:, 1:m) \ (Q' * b))]
+end
address@hidden group
address@hidden example
+
address@hidden
+While this approach is certainly faster and uses less memory than
+Octave's minimum norm approach, this approach seems to be inferior in
+other ways.
+
+A numerical question arises: how big can the null space component become,
+relative to the minimum-norm solution? Can it be nicely bounded, or can it be
+arbitrarily big? Consider this example:
+
address@hidden
address@hidden
+m = 10; 
+n = 10000; 
+A = ones(m, n) + 1e-6 * randn(m,n); 
+b = ones(m, 1) + 1e-6 * randn(m,1); 
+norm(A \ b)
address@hidden group
address@hidden example
+
address@hidden
+while Octave's minimum-norm values are around 3e-2, Matlab's results
+are 50-times larger. For another issue, try this code:
+
address@hidden
address@hidden
+m = 5; 
+n = 100; 
+j = floor(m * rand(1, n)) + 1; 
+b = ones(m, 1);
+A = zeros(m, n);
+A(sub2ind(size(A),j,1:n)) = 1;
+x = A \ b; 
+[dummy,p] = sort(rand(1,n)); 
+y = A(:,p)\b; 
+norm(x(p)-y)
address@hidden group
address@hidden example
+
address@hidden
+It shows that unlike in Octave, mldivide in Matlab is not invariant
+with respect to column permutations. If there are multiple columns of
+the same norm, permuting columns of the matrix gets you different
+result than permuting the solution vector. This will surprise many
+users.
+
+Since the mldivide (\) and mrdivide (/) operators are often part of a more 
+complex expression, where there is no room to react to warnings or flags, it 
+should prefer intelligence (robustness) to speed, and so the Octave developers
+are firmly of the opinion that Octave's approach for singular, under- and
+over-determined matrices is a better choice that Matlab's
+
 @item Octave extensions
 The extensions in Octave over @sc{Matlab} syntax are
 very useful, but might cause issues when sharing with @sc{Matlab} users.