Actual source code: lsqr_converged.c
1: #include "private/kspimpl.h"
2: #include ../src/ksp/ksp/impls/lsqr/lsqr.h
5: PetscErrorCode KSPConvergedLSQR(KSP solksp,PetscInt iter,PetscReal rnorm,KSPConvergedReason *reason,void *ctx)
6: {
7: PetscReal atol; /* Absolute convergence criterium */
8: PetscReal dtol; /* Divergence criterium */
9: PetscReal rtol; /* Relative convergence criterium */
10: PetscReal stop1; /* Stop if: |r| < rtol*|b| + atol*|A|*|x| */
11: PetscReal stop2; /* Stop if: |A^t.r|/(|A|*|r|) < atol */
12: Vec x_sol; /* Current solution vector */
14: PetscReal arnorm, anorm, bnorm, xnorm; /* Norms of A*residual; matrix A; rhs; solution */
16: PetscInt mxiter; /* Maximum # of iterations */
20: *reason = KSP_CONVERGED_ITERATING;
21: if (iter == 0) {return(0);};
23: if (PetscIsInfOrNanScalar(rnorm)) {
24: *reason = KSP_DIVERGED_NAN;
25: return(0);
26: }
28: KSPGetTolerances( solksp, &rtol, &atol, &dtol, &mxiter );
29: if ( iter > mxiter ) {
30: *reason = KSP_DIVERGED_ITS;
31: return(0);
32: }
34: KSPGetSolution( solksp, &x_sol );
35: VecNorm(x_sol, NORM_2 , &xnorm);
37: KSPLSQRGetArnorm( solksp, &arnorm, &bnorm, &anorm);
38: if (bnorm > 0.0) {
39: stop1 = rnorm / bnorm;
40: rtol = rtol + atol * anorm*xnorm/bnorm;
41: } else {
42: stop1 = 0.0;
43: rtol = 0.0;
44: }
45: stop2 = 0.0;
46: if (rnorm > 0.0) stop2 = arnorm / (anorm * rnorm );
48: /* Test for tolerances set by the user */
49: if ( stop1 <= rtol ) *reason = KSP_CONVERGED_RTOL;
50: if ( stop2 <= atol ) *reason = KSP_CONVERGED_ATOL;
52: /* Test for machine precision */
53: if (bnorm > 0) {
54: stop1 = stop1 / (1.0 + anorm*xnorm/bnorm);
55: } else {
56: stop1 = 0.0;
57: }
58: stop1 = 1.0 + stop1;
59: stop2 = 1.0 + stop2;
60: if ( stop1 <= 1.0 ) *reason = KSP_CONVERGED_RTOL;
61: if ( stop2 <= 1.0 ) *reason = KSP_CONVERGED_ATOL;
62: return(0);
63: }