Consider this (MATLAB R2011a):
a = 1e10;
>> b = inv(a)*inv(a)
b =
1.0000e-020
>> c = inv(a*a)
c =
1.0000e-020
>> b==c
ans =
0
>> format hex
>> b
b =
3bc79ca10c924224
>> c
c =
3bc79ca10c924223
When MATLAB calculates the intermediate quantities inv(a)
, or a*a
(whether a
is a scalar or a matrix), it by default stores them as the closest double precision floating point number - which is not exact. So when these slightly inaccurate intermediate results are used in subsequent calculations, there will be round off error.
Instead of comparing floating point numbers for direct equality, such as inv(A*B*C) == inv(C)*inv(B)*inv(A)
, it's often better to compare the absolute difference to a threshold, such as abs(inv(A*B*C) - inv(C)*inv(B)*inv(A)) < thresh
. Here thresh
can be an arbitrary small number, or some expression involving eps
, which gives you the smallest difference between two numbers at the precision at which you're working.
The format
command only controls the display of results at the command line, not the way in which results are internally stored. In particular, format rat
does not make MATLAB do calculations symbolically. For this, you might take a look at the Symbolic Math Toolbox. format hex
is often even more useful than format long
for diagnosing floating point precision issues such as the one you've come across.