On 8 November 2010 14:38, Joon <groups.and.li...@gmail.com> wrote: > Oh I see. So I guess in invA = solve(Ax, I) and then x = dot(invA, b) case, > there are more places where numerical errors occur, than just x = solve(Ax, > b) case.
That's the heart of the matter, but one can be more specific. You can think of a matrix by how it acts on vectors. Taking the inverse amounts to solving Ax=b for all the standard basis vectors (0,...,0,1,0,...,0); multiplying by the inverse amounts to expressing your vector in terms of these, finding where they go, and adding them together. But it can happen that when you break your vector up like that, the images of the components are large but almost cancel. This sort of near-cancellation amplifies numerical errors tremendously. In comparison, solving directly, if you're using a stable algorithm, is able to avoid ever constructing these nearly-cancelling combinations explicitly. The standard reason for trying to construct an inverse is that you want to solve equations for many vectors with the same matrix. But most solution methods are implemented as a matrix factorization followed by a single cheap operation, so if this is your goal, it's better to simply keep the matrix factorization around. Anne _______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
