On a related note, the following is in the July 16 (Vol 15 #28) NA-Digest:
From: Mike Sussman sussmanm@math.pitt.edu
Date: July 09, 2015
Subject: MATLAB inv() much less stable in new releases than old
Does anyone know why the folks at the Mathworks have chosen to change the matrix algorithm for the inv function so that it is much less stable than before? Between MATLAB releases R2012a and R2014a the following sequence of instructions produces drastically different
results:
A=hilb(10);
Ainv=inv(A);
x=ones(10,1);
b=Ax;
norm(x-Ainvb)/norm(x)
Using release R2012a, this produces 2.9743e-4, somewhat better than would be expected from the condition number of 1.6025e+13. If, instead, Ainv is computed as Ainv=A\eye(10), the result is 0.0051, larger than with inv(A) but still consistent with the condition number. In contrast, using release R2014a (and more recent versions) one gets 1.4290e+5, substantially larger than one would expect from the condition number.
If I recall correctly, this change occurred in or before R2013a, but I can no longer verify that. In any case I have been unable to find any mention of a change in the release notes or on the internet.
Can anyone help me understand why such a change was made? Is the new algorithm better in some way for some class of problem?
Of course, we all know it’s naughty to use inv, but sometimes one gives into temptation, or it really, truly is necessary to explicitly compute the inverse.
By the way, congratulations to Stephen R. Becker, Emmanuel J. Candès, Michael C. Grant (mcg) for winning the 2015 Beale-Orchard-Hays Prize for TFOCS. This makes two in a row for mcg, who, along with Stephen Boyd, won the Beale-Orchard-Hays Prize in 2012 (when it was msot recently previously awarded) for CVX.