LAPACK Archives

[Lapack] ZGELS v.s. ZGESV

your question amounts to ask:

what is the difference between solving a square nonsingular system of 
equations Ax=b with an LU with partial pivoting factorization (say LU-GEPP 
for Gaussian Elimination with Partial Pivoting) of the matrix A (ZGESV) or 
a QR factorization using Householder reflections (say QR-HH) of the matrix 

1- The cost of LU-GEPP is 2/3 n^3 while the cost of the QR-HH 
factorization 4/3 n^3 (for square n-by-n matrix). So LU-GEPP (DGESV) is 
expected to be twice as fast as QR-HH (DGELS).

2- QR-HH is unconditionnally backward stable. This means that's for any 
input matrix A, the computed solution _x_ from QR-HH is the exact solution 
of a slightly perturbed problem. In other words, _x_ is such that:

        ( A + Delta_A ) _x_ = ( b + delta_b )
where || Delta_A || <= small_constant * machine_eps * || A ||
and   || delta_b || <= small_constant * machine_eps * || b ||.

Such inequality (backward stability) is the minimum one ask for any 
reasonnable algorithm.
LU-GEPP has such a good behavior for all the matrix that we have 
encountered except a few. So LU-GEPP is not backward stable in theory but 
is in practice. So in practice both methods will provide you with a 
backward stable solution.

3- if the system is close to singular (we say that it is ill-conditionned) 
then you can have a lot of backward stable solutions. Basically the 
relation is

        || _x_ - true_solution || <= condition_number_A * bwd_error

Both LU_GEPP and QR-HH will give you a small bwd_error but if the 
condition number of A is large then the computed solution _x_ can be 
pretty far from the true solution. Moreover there is no reason for the 
computed solution of LU-GEPP to be the same as the computed solution of 
QR-HH and in practice you will observe that they differ by a factor 
machine_precision * condition_number_A.

However you can not tell one solution is better than the other. They are 
both backward stable solutions. So they are both 'acceptable' answer.

4- If the matrix is ill-conditionned, instead of using ZGELS, you can use 
other algorithms. They will filter the small singular values out of your 
initial matrix A, once those are removed your system becomes 
overdetermined and we solve a least squares problem on it. Those routines 
are ZGELSY (QR with column pivoting), ZGELSD (SVD with Divide and 
Conquer), ZGELSS (SVD with QR).

If you want to know more, you might want to look at a dense linear algebra 
book. (For example: Demmel, Trefethen and Bau, or Golub and Van Loan).


On Fri, 24 Aug 2007, Sorkin Dima wrote:

 What is the difference between solving a _square_
system of equations with ZGELS vs. ZGESV.
Can they give different results ?
What happens when matrix is very close to singular ?

Thank you, regards,
Lapack mailing list

<Prev in Thread] Current Thread [Next in Thread>

For additional information you may use the LAPACK/ScaLAPACK Forum.
Or one of the mailing lists, or