Hi AB,
This is an interesting question. The conventional wisdom is based on the computational complexity, meaning it prioritizes large matrices. For small matrices, things become more complicated (and, as John mentioned, much more painful to measure and understand due to noise).
BLOT: For small matrices, estimating the condition number in mldivide takes a considerable time of the computation. If you're confident in your conditioning, use pagemldivide which doesn't compute the condition number unless a second output is requested, which is faster than inv and mldivide.
We have recently taken an interest in matrices with smaller sizes, as part of the introduction of page* functions. These are functions that take a 3D array representing a stack of matrices, and apply the same operation to each matrix in this stack (each "page" A(:, :, i)); which can be much faster than using a for-loop due to reduction in function call overhead.
Remember that MATLAB's backslash (or mldivide) will warn if the left matrix is ill-conditioned:
magic(10) \ randn(10, 1);
For us to give this warning, we need to estimate the condition number of that matrix, which is done through an iterative algorithm, solving several linear systems. For large matrices, this isn't too expensive compared to the factorization itself ((O(n^2) to solve linear system, O(n^3) for the factorization); but for small matrices, it can become dominant.
In pagemldivide, which is focused on many small matrices, we decided not to warn when one of the matrices is ill-conditioned. Instead, a second output gives the estimated condition number of each matrix. This means we don't need to estimate the condition numbers when the second output isn't requested.
Let's compare inv, mldivide, and pagemldivide timings for your example:
function tim = timeMldivide(n)
A = randn(n);
A = A*A';
s = randn(n,1);
rep = 1e7 / n^2;
tic;
for ii=1:rep
y = inv(A)*s;
end
tInv = toc;
tic;
for ii=1:rep
y = A \ s;
end
tMldivide = toc;
tic;
for ii=1:rep
y = pagemldivide(A, s);
end
tPagemldivide = toc;
if nargout == 0
tInv
tMldivide
tPagemldivide
else
tim = [tInv tMldivide tPagemldivide] ./ rep;
end
end
timeMldivide(20)
We can see that inv is faster than mldivide, and pagemldivide is faster than both. So, if you know your problem to be well-conditioned, pagemldivide is a better choice than inv.
You may wonder: Doesn't inv also warn for ill-conditioned inputs - meaning it also needs to compute the condition number? This is true, but the condition number is defined as cond(A) = norm(A, 1) * norm(inv(A), 1). Since inv(A) is already is computed, it's cheaper to just compute its 1-norm, instead of using an iterative estimator as is done in mldivide.
To finish, let's look at timings for different matrix sizes, to see how these timings behave there:
nlist = [2 3 5 8 12 16 20 32 64 128 256 512];
tim = [];
for n=nlist
tim = [tim; timeMldivide(n)];
end
loglog(nlist, tim);
legend('inv', 'mldivide', 'pagemldivide', Location='northwest')
xlabel('Matrix size')
ylabel('Time (s)')
We can see pagemldivide (which doesn't check condition number) is fastest in all cases, although the difference to mldivide decreases as the matrix size increases. We see inv becomes much slower for large sizes, but is faster than mldivide for sizes 20 and 32 (the specific bounds will depend on the machine, and can change from release to release).
That mldivide is faster for very small matrices I would put down to the special-case treatment in mldivide for very small matrices - we always use LU there, because the cost of structure detection was becoming noticeable for these very small matrices (see the doc for mldivide). In inv, this structure detection is needed, because inv for Hermitian matrices returns an output that is also Hermitian - and this requires a separate algorithm.
