Forward and back substitutions

You can use either the slash operator ('/') or olsqr to solve a system of linear equations without using the inverse. The slash operator uses an LU based solver and as you may guess, olsqr uses the QR decomposition. Neither ever computes an inverse.

x = b/A;

or

x = olsqr(b, A);

If you want to solve x*A = b, you should be able to do this:

x = (b' / A')';

Great help, thanks! I knew the slash operator, but not the other one.

I played with arbitrary matrices in Gauss and found the following equivalent solutions.

Assume the system Ax=b where A=PP', i.e., P=chol(A)'. Then x can be calculated in the following ways:

x=inv(A)*b=b/A=(b/P)/P'=utrisol(ltrisol(b,P),P')

Likewise, when we have xA=b, equivalent ways of calculating x are

x=b*inv(A)=(b'/A')'=((b'/P)/P')'=utrisol(ltrisol(b',P),P')'

Note that in the first case b is a column vector, while in the second case it's a row vector.

Now my question is: among the above equivalent ways which is the numerically most stable and less subject to round-off errors?

The solution using the inverse is the least numerically stable. Using the QR decomposition (olsqr) should be the most numerically stable. The LU based solve (using the slash operator) should be in between. Let me see if I can get someone else to comment further.

Thanks for the reply and I'm grateful for further clarifications. I just tried simple and small moment matrices where numerical issues were irrelevant. However, in a filter where covariance matrices are propagated, the above alternatives may lead to different solutions. Of course, smart scaling may help, but it would be nice to know more about the numerical properties of the different solutions.