6
$\begingroup$

I want to have an estimation, that my solution has an error, let's say less than 1e-8.

Usually, I stop the Gauss-Seidel algorithm, when the residual is "small enough" and this is already the problem. How do I know when the residual is small enough, because even when the residual is small, the solution may still have too much error. So this is no good method.

What do you use as stopping criterion?

On another website (math-linux.com) I found a stopping criterion: $$ \|r\|/\|b\|\leq \epsilon $$ But again, what theory is behind that?

This, by the way, is the code I used in my last project, just for information, how I did it:

void relax(double epsilon, vector<double> &x, SparseMatrix &A, const vector<double> &f) {
    int maxIter = 100;
    int iter = 0;
    double residual = 1.0;
    double minResidual 0.000001; //I also tried 1e-14;

    while (iter < maxIter && residual >= minResidual) {
        for (int i = 0; i < A.dim; ++i) {
            double ls = A.lineScalar(i, x);
            x[i] = (1.0/A(i,i)) * (f[i] - ls);
        }
        vector<double> temp = A.multiply(x);
        residual = L2(temp - f);
    }
}
Improve this question
$\endgroup$
5
  • $\begingroup$ Seems to me pretty inefficient to test the residual every iteration this way. $\endgroup$
    – Vladimir F Героям слава
    Jul 25, 2012 at 8:29
  • $\begingroup$ Yes, that's true. $\endgroup$
    – vanCompute
    Jul 25, 2012 at 8:35
  • 3
    $\begingroup$ scicomp.stackexchange.com/q/582/119 $\endgroup$
    – Jed Brown
    Jul 25, 2012 at 12:01
  • $\begingroup$ I would like to pose another question for this topic: Can we estimate the lowest singular value of $A$ from below in a reasonable way provided the matrix has additional structure, like strict diagonal dominance? As far as I know, Kaasschieter has described such a procedure for the CGM. $\endgroup$
    – shuhalo
    Jul 25, 2012 at 13:13
  • 2
    $\begingroup$ @Martin: You should pose this as a separate question rather than a comment. $\endgroup$
    – Arnold Neumaier
    Jul 25, 2012 at 13:52

2 Answers 2

Reset to default
7
$\begingroup$

One cannot conclude from a residual how accurate the solution is. Between the best and the worst case in norm, there is a factor of exactly the condition number. More precisely, if the residual norm is r and the error norm is e then $\|A\|^{-1}\le e/r \le \|A^{-1}\|$, and both bounds are attainable. Taking the quotient of the bounds proves the claim.

The order of magnitude of the norm of the inverse can often be estimated from a few steps of the Lanczos/Arnoldi method. Thus one can apply the bounds.

However, generally, model equations are not really exact anyway, and getting the residual smaller than the accuracy with which uncertain data determine the residual is considered the scientifically corrrect procedure. Depending on the data, this may be far earlier that stagnation.

The theory behind this is called backward error analysis. For example, if the only source of uncertainty are measurement errors in the right hand side of 1e-4 in some norm, and your residuum norm is smaller than 1e-4, it cannot be distinguished from the residuum the exact solution of your solved problem would have when used as approximation in a problem with a different rhs of the same accuracy.

Improve this answer
$\endgroup$
9
  • $\begingroup$ As one needs $\|A^{-1}\|$ this approach is not feasible, as one has to know the solution for calculation of the error. It may be better to look, wether the residuum stagnates? $\endgroup$
    – vanCompute
    Jul 25, 2012 at 12:45
  • $\begingroup$ The order of magnitude of the inverse norm can often be estimated from a few steps of the Lanczos/Arnoldi method. - The fact that one cannot conclude from a residual how accurate the solution is is independent of the knowledge of the value of the inverse norm. - Generally, model equations are not really exact anyway, and getting the residual smaller than the accuracy with which uncertain data determine the residual is considered the scientifically corrrect procedure. Depending on the data, this may be far earlier that stagnation. $\endgroup$
    – Arnold Neumaier
    Jul 25, 2012 at 13:03
  • $\begingroup$ Ok, if my matrix A comes from a PDE and the right hand side is a measurement with accuracy 1e-4 (for example), then I would stop ich the residuum is smaller then 1e-4? But what theory is behind that? $\endgroup$
    – vanCompute
    Jul 25, 2012 at 13:26
  • 1
    $\begingroup$ The theory behind it is called backward error analysis. If your residuum is smaller than 1e-4 it cannot be distinguished from the residuum the the exact solution of your particular problem would have when used as approximation in a problem with a different rhs of the same accuracy. $\endgroup$
    – Arnold Neumaier
    Jul 25, 2012 at 13:51
  • 1
    $\begingroup$ "order of magnitude of the inverse norm can often be estimated from a few steps" -- A "few steps" is very unreliable, it's easy to be orders of magnitude off after 100 iterations. Arnoldi is worse than Lanczos, but in either case, you really have to be converged to some modest tolerance to have a useful estimate. In contrast, a few iterations (less than 10) almost always gives you an accurate estimate of the largest singular value of $A$. $\endgroup$
    – Jed Brown
    Jul 25, 2012 at 15:41
1
$\begingroup$

Usually, when I'm estimating a solution of a system of linear equations, I save the approximation $\vec x^{n-1}$ and use it to compute $x_{err}=max|x_i^n-x_i^{n-1}|$ over each component i. This is faster than calculating the residual, since it doesn't require a matrix-vector product. Although, if your system is ill-conditioned, it still may not give you the correct solution for a given tolerance $\epsilon$.

Improve this answer
$\endgroup$
6
  • 1
    $\begingroup$ Quite interesting. I mean, it's like a Cauchy sequence, but the problem, what I have here is that I cannot bound the error. I only know, that the alogorithm is still converging, if I am right. $\endgroup$
    – vanCompute
    Jul 25, 2012 at 20:35
  • $\begingroup$ This is bad! Now you cannot distinguish between stagnation and convergence. Most iterative methods either embed an explicit residual evaluation or have one available through a recurrence relation, therefore you can use much better methods at no extra cost. In any case, it is never okay to simply declare convergence when your method is converging slowly. It is easy to construct problems for which that approach is arbitrarily bad. $\endgroup$
    – Jed Brown
    Jul 25, 2012 at 20:37
  • 1
    $\begingroup$ @JedBrown: True, but the same can be said of residual evaluation: one can construct problems which residual fails to distinguish between the two. $\endgroup$
    – Paul
    Jul 25, 2012 at 20:43
  • $\begingroup$ Yes, if I am looking at the differences between residuals, I could run into "slow convergence", too. Is that, what you mean, Paul? $\endgroup$
    – vanCompute
    Jul 25, 2012 at 20:50
  • $\begingroup$ @JedBrown: What would be a method with no extra cost? $\endgroup$
    – vanCompute
    Jul 25, 2012 at 20:52

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.