Fix early success in Bicgstab: Should not scale residual by norm of rhs vector b...

git-svn-id: https://svn.dealii.org/trunk@29472 0785d39b-7218-0410-832d-ea1e28bc413d

diff --git a/deal.II/include/deal.II/lac/solver_bicgstab.h b/deal.II/include/deal.II/lac/solver_bicgstab.h
index 3d2bb4bbe4cbefdbebd354776c733b362a72279e..f81d517d52090c7584f28bc566993e869e20f41a 100644 (file)
--- a/deal.II/include/deal.II/lac/solver_bicgstab.h
+++ b/deal.II/include/deal.II/lac/solver_bicgstab.h
@@ -63,28 +63,22 @@ class SolverBicgstab : public Solver<VECTOR>
 {
 public:
   /**
-   * There are two possibilities to
-   * compute the residual: one is an
-   * estimate using the computed value @p
-   * tau. The other is exact computation
-   * using another matrix vector
-   * multiplication. This increases the
-   * costs of the algorithm, so it is
-   * should be set to false whenever the
-   * problem allows it.
+   * There are two possibilities to compute the residual: one is an estimate
+   * using the computed value @p tau. The other is exact computation using
+   * another matrix vector multiplication. This increases the costs of the
+   * algorithm, so it is should be set to false whenever the problem allows
+   * it.
    *
-   * Bicgstab is susceptible to breakdowns, so
-   * we need a parameter telling us, which
-   * numbers are considered zero.
+   * Bicgstab is susceptible to breakdowns, so we need a parameter telling us,
+   * which numbers are considered zero.
    */
   struct AdditionalData
   {
     /**
      * Constructor.
      *
-     * The default is to perform an
-     * exact residual computation and
-     * breakdown parameter 1e-10.
+     * The default is to perform an exact residual computation and breakdown
+     * parameter 1e-10.
      */
     AdditionalData(const bool   exact_residual = true,
                    const double breakdown      = 1.e-10) :
@@ -109,9 +103,8 @@ public:
                   const AdditionalData &data=AdditionalData());
 
   /**
-   * Constructor. Use an object of
-   * type GrowingVectorMemory as
-   * a default to allocate memory.
+   * Constructor. Use an object of type GrowingVectorMemory as a default to
+   * allocate memory.
    */
   SolverBicgstab (SolverControl        &cn,
                   const AdditionalData &data=AdditionalData());
@@ -139,13 +132,9 @@ protected:
   double criterion (const MATRIX &A, const VECTOR &x, const VECTOR &b);
 
   /**
-   * Interface for derived class.
-   * This function gets the current
-   * iteration vector, the residual
-   * and the update vector in each
-   * step. It can be used for a
-   * graphical output of the
-   * convergence history.
+   * Interface for derived class.  This function gets the current iteration
+   * vector, the residual and the update vector in each step. It can be used
+   * for a graphical output of the convergence history.
    */
   virtual void print_vectors(const unsigned int step,
                              const VECTOR &x,
@@ -329,9 +318,8 @@ SolverBicgstab<VECTOR>::iterate(const MATRIX &A,
   VECTOR &t = *Vt;
   VECTOR &v = *Vv;
 
-  v = 0;
-  p = 0;
   rbar = r;
+  bool startup = true;
 
   do
     {
@@ -340,7 +328,14 @@ SolverBicgstab<VECTOR>::iterate(const MATRIX &A,
       rhobar = r*rbar;
       beta   = rhobar * alpha / (rho * omega);
       rho    = rhobar;
-      p.sadd(beta, 1., r, -beta*omega, v);
+      if (startup == true)
+        {
+          p = r;
+          startup = false;
+        }
+      else
+        p.sadd(beta, 1., r, -beta*omega, v);
+
       precondition.vmult(y,p);
       A.vmult(v,y);
       rhobar = rbar * v;
@@ -354,12 +349,9 @@ SolverBicgstab<VECTOR>::iterate(const MATRIX &A,
 
       r.add(-alpha, v);
 
-      // check for early success, see
-      // the lac/bicgstab_early
-      // testcase as to why this is
-      // necessary
-      if (this->control().check(step, r.l2_norm()/Vb->l2_norm())
-          == SolverControl::success)
+      // check for early success, see the lac/bicgstab_early testcase as to
+      // why this is necessary
+      if (this->control().check(step, r.l2_norm()) == SolverControl::success)
         {
           Vx->add(alpha, y);
           print_vectors(step, *Vx, r, y);
@@ -437,8 +429,7 @@ SolverBicgstab<VECTOR>::solve(const MATRIX &A,
 
   deallog.pop();
 
-  // in case of failure: throw
-  // exception
+  // in case of failure: throw exception
   if (this->control().last_check() != SolverControl::success)
     throw SolverControl::NoConvergence (this->control().last_step(),
                                         this->control().last_value());

diff --git a/tests/lac/solver/cmp/generic b/tests/lac/solver/cmp/generic
index 8cef71a4b49fc976125132fe4fcdccb31c2cffeb..73f5d12b8928de46095df1cada5261ccf8f8bac7 100644 (file)
--- a/tests/lac/solver/cmp/generic
+++ b/tests/lac/solver/cmp/generic
@@ -72,7 +72,7 @@ DEAL:psor:cg::Starting value 3.000
 DEAL:psor:cg::Failure step 100 value 0.1024
 DEAL:psor::Iterative method reported convergence failure in step 100 with residual 0.102391
 DEAL:psor:Bicgstab::Starting value 3.000
-DEAL:psor:Bicgstab::Convergence step 4 value 0.0002656
+DEAL:psor:Bicgstab::Convergence step 4 value 0.0007969
 DEAL:psor:GMRES::Starting value 1.467
 DEAL:psor:GMRES::Convergence step 5 value 0.0006649
 DEAL:psor:GMRES::Starting value 3.000
@@ -83,7 +83,9 @@ DEAL:no-fail:cg::Starting value 11.00
 DEAL:no-fail:cg::Failure step 10 value 0.1496
 DEAL:no-fail::Iterative method reported convergence failure in step 10 with residual 0.149566
 DEAL:no-fail:Bicgstab::Starting value 11.00
-DEAL:no-fail:Bicgstab::Convergence step 10 value 0.0002572
+DEAL:no-fail:Bicgstab::Failure step 10 value 0.002830
+DEAL:no-fail:Bicgstab::Failure step 10 value 0.001961
+DEAL:no-fail::Iterative method reported convergence failure in step 10 with residual 0.00196051
 DEAL:no-fail:GMRES::Starting value 11.00
 DEAL:no-fail:GMRES::Failure step 10 value 0.8414
 DEAL:no-fail::Iterative method reported convergence failure in step 10 with residual 0.841354
@@ -99,7 +101,7 @@ DEAL:no::Iterative method reported convergence failure in step 100 with residual
 DEAL:no:cg::Starting value 11.00
 DEAL:no:cg::Convergence step 15 value 0.0005794
 DEAL:no:Bicgstab::Starting value 11.00
-DEAL:no:Bicgstab::Convergence step 10 value 0.0002572
+DEAL:no:Bicgstab::Convergence step 11 value 0.0006357
 DEAL:no:GMRES::Starting value 11.00
 DEAL:no:GMRES::Convergence step 43 value 0.0009788
 DEAL:no:GMRES::Starting value 11.00
@@ -112,7 +114,7 @@ DEAL:rich::Iterative method reported convergence failure in step 100 with residu
 DEAL:rich:cg::Starting value 11.00
 DEAL:rich:cg::Convergence step 15 value 0.0005794
 DEAL:rich:Bicgstab::Starting value 11.00
-DEAL:rich:Bicgstab::Convergence step 10 value 0.0002572
+DEAL:rich:Bicgstab::Convergence step 11 value 0.0006357
 DEAL:rich:GMRES::Starting value 6.600
 DEAL:rich:GMRES::Convergence step 42 value 0.0007803
 DEAL:rich:GMRES::Starting value 11.00
@@ -126,7 +128,7 @@ DEAL:ssor:Richardson::Convergence step 48 value 0.0009502
 DEAL:ssor:cg::Starting value 11.00
 DEAL:ssor:cg::Convergence step 8 value 0.0005816
 DEAL:ssor:Bicgstab::Starting value 11.00
-DEAL:ssor:Bicgstab::Convergence step 4 value 0.0003834
+DEAL:ssor:Bicgstab::Convergence step 5 value 0.0002498
 DEAL:ssor:GMRES::Starting value 13.27
 DEAL:ssor:GMRES::Convergence step 7 value 0.0008206
 DEAL:ssor:GMRES::Starting value 11.00
@@ -141,7 +143,7 @@ DEAL:sor:cg::Starting value 11.00
 DEAL:sor:cg::Failure step 100 value 5.643
 DEAL:sor::Iterative method reported convergence failure in step 100 with residual 5.64329
 DEAL:sor:Bicgstab::Starting value 11.00
-DEAL:sor:Bicgstab::Convergence step 12 value 0.0007246
+DEAL:sor:Bicgstab::Convergence step 14 value 0.0009987
 DEAL:sor:GMRES::Starting value 7.322
 DEAL:sor:GMRES::Convergence step 23 value 0.0006981
 DEAL:sor:GMRES::Starting value 11.00
@@ -154,7 +156,7 @@ DEAL:psor:cg::Starting value 11.00
 DEAL:psor:cg::Failure step 100 value 2.935
 DEAL:psor::Iterative method reported convergence failure in step 100 with residual 2.93488
 DEAL:psor:Bicgstab::Starting value 11.00
-DEAL:psor:Bicgstab::Convergence step 10 value 0.0003167
+DEAL:psor:Bicgstab::Convergence step 11 value 0.0005151
 DEAL:psor:GMRES::Starting value 7.345
 DEAL:psor:GMRES::Convergence step 20 value 0.0008491
 DEAL:psor:GMRES::Starting value 11.00

diff --git a/tests/lac/solver_selector_01/cmp/generic b/tests/lac/solver_selector_01/cmp/generic
index 200218e23028bbb68394adab4d6c82b25dd9de1b..8f2bd5c7499b0248b61ede61c2ef6e47503d2e84 100644 (file)
--- a/tests/lac/solver_selector_01/cmp/generic
+++ b/tests/lac/solver_selector_01/cmp/generic
@@ -3,7 +3,7 @@ DEAL::Size 37 Unknowns 1296
 DEAL:cg::Starting value 36.00
 DEAL:cg::Convergence step 21 value 0.003095
 DEAL:Bicgstab::Starting value 36.00
-DEAL:Bicgstab::Convergence step 10 value 0.003158
+DEAL:Bicgstab::Convergence step 16 value 0.003278
 DEAL:GMRES::Starting value 34.45
 DEAL:GMRES::Convergence step 21 value 0.001416
 DEAL:FGMRES::Starting value 36.00
@@ -11,7 +11,7 @@ DEAL:FGMRES::Convergence step 21 value 0.002535
 DEAL:cg::Starting value 36.00
 DEAL:cg::Convergence step 39 value 7.335e-08
 DEAL:Bicgstab::Starting value 36.00
-DEAL:Bicgstab::Convergence step 27 value 1.148e-08
+DEAL:Bicgstab::Convergence step 28 value 9.087e-08
 DEAL:GMRES::Starting value 34.45
 DEAL:GMRES::Convergence step 39 value 6.627e-08
 DEAL:FGMRES::Starting value 36.00

diff --git a/tests/petsc/deal_solver_02/cmp/generic b/tests/petsc/deal_solver_02/cmp/generic
index 4cbd6cc5d7e76fbbf267d8f1a61d64d278968124..d79a37c605f959f1d27f66b0854776b40e64e3c8 100644 (file)
--- a/tests/petsc/deal_solver_02/cmp/generic
+++ b/tests/petsc/deal_solver_02/cmp/generic
@@ -2,5 +2,5 @@
 DEAL::Size 32 Unknowns 961
 DEAL::Solver type: N6dealii14SolverBicgstabINS_13PETScWrappers6VectorEEE
 DEAL:Bicgstab::Starting value 31.00
-DEAL:Bicgstab::Convergence step 29 value 0.0008735
-DEAL::Solver stopped after 29 iterations
+DEAL:Bicgstab::Convergence step 35 value 0.0005391
+DEAL::Solver stopped after 35 iterations

diff --git a/tests/trilinos/deal_solver_02/cmp/generic b/tests/trilinos/deal_solver_02/cmp/generic
index 664748dd19961a5d48a4415d23990eac6f3ea915..9bead2f43a0ca1dcf7da912e345f3def0c05b111 100644 (file)
--- a/tests/trilinos/deal_solver_02/cmp/generic
+++ b/tests/trilinos/deal_solver_02/cmp/generic
@@ -2,5 +2,5 @@
 DEAL::Size 32 Unknowns 961
 DEAL::Solver type: N6dealii14SolverBicgstabINS_16TrilinosWrappers6VectorEEE
 DEAL:Bicgstab::Starting value 31.00
-DEAL:Bicgstab::Convergence step 29 value 0.0008735
-DEAL::Solver stopped after 29 iterations
+DEAL:Bicgstab::Convergence step 35 value 0.0005391
+DEAL::Solver stopped after 35 iterations