<h3>lac</h3>
<ol>
+ <li> <p>Fixed: The <code class="class">SolverBicgstab</code> class
+ did not handle hitting on the solution early very
+ gracefully (it threw an exception). This is now fixed.
+ <br>
+ (Roger Young 2007/03/07)
+ </p>
+
<li> <p>Fixed: The <code class="member">SparseDirectMA27</code> class allows
to run the sparse direct solver as a separate program (the
<code>detached_ma27</code> program) rather than as part of the current
// $Id$
// Version: $Name$
//
-// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
+// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
return true;
s.equ(1., r, -alpha, v);
+
+ // check for early success, see
+ // the lac/bicgstab_early
+ // testcase as to why this is
+ // necessary
+ if (this->control().check(step, s.l2_norm()/Vb->l2_norm())
+ == SolverControl::success)
+ {
+ Vx->add(alpha, y);
+ print_vectors(step, *Vx, r, y);
+ return false;
+ }
+
precondition.vmult(z,s);
A.vmult(t,z);
rhobar = t*s;
solver eigen \
sparse_ilu sparse_ilu_t lapack_fill block_minres \
identity_matrix* trace \
- print_formatted_*
+ print_formatted_* \
+ bicgstab_early
# from above list of regular expressions, generate the real set of
--- /dev/null
+//---------------------------- bicgstab_early.cc ---------------------------
+// $Id: bicgstab_early.cc 11749 2005-11-09 19:11:20Z wolf $
+// Version: $Name$
+//
+// Copyright (C) 2007 by the deal.II authors
+//
+// This file is subject to QPL and may not be distributed
+// without copyright and license information. Please refer
+// to the file deal.II/doc/license.html for the text and
+// further information on this license.
+//
+//---------------------------- bicgstab_early.cc ---------------------------
+
+// adapted from a testcase by Roger Young, sent to the mailing list
+// 2007-03-02, that illustrates that bicgstab can't handle early
+// success
+
+#include "../tests.h"
+#include <cmath>
+#include <fstream>
+#include <iostream>
+#include <iomanip>
+#include "testmatrix.h"
+#include <base/logstream.h>
+#include <lac/sparse_matrix.h>
+#include <lac/sparse_ilu.h>
+#include <lac/solver_bicgstab.h>
+#include <lac/vector.h>
+#include <lac/precondition.h>
+
+
+int main()
+{
+ std::ofstream logfile("bicgstab_early/output");
+ logfile.precision(4);
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+ deallog.threshold_double(1.e-10);
+
+ GrowingVectorMemory<> mem;
+ SolverControl control(100, 1.e-3);
+ SolverBicgstab<> bicgstab(control, mem);
+
+ SparsityPattern sparsity_pattern(4,4,4);
+ for (unsigned int i=0; i<4; ++i)
+ for (unsigned int j=0; j<4; ++j)
+ sparsity_pattern.add(i,j);
+ sparsity_pattern.compress();
+
+ SparseMatrix<double> M(sparsity_pattern);
+ M.set(0,0,21.1);
+ M.set(0,1,0);
+ M.set(0,2,0);
+ M.set(0,3,0);
+ M.set(1,1,7.033333333);
+ M.set(1,0,0);
+ M.set(1,2,0);
+ M.set(1,3,3.516666667);
+ M.set(2,2,21.1);
+ M.set(2,0,0);
+ M.set(2,1,0);
+ M.set(2,3,0);
+ M.set(3,3,7.033333333);
+ M.set(3,0,0);
+ M.set(3,1,3.516666667);
+ M.set(3,2,0);
+
+ Vector<double> rhs(4);
+ rhs(0) = rhs(2) = 0;
+ rhs(1) = rhs(3) = 0.0975;
+
+
+ SparseILU<double> ilu (sparsity_pattern);
+ ilu.decompose (M);
+
+ Vector<double> solution (4);
+ // this would fail with elements of
+ // the solution vector being set to
+ // NaN before Roger's suggested
+ // change
+ bicgstab.solve (M, solution, rhs, ilu);
+
+ for (unsigned int i=0; i<4; ++i)
+ deallog << solution(i) << std::endl;
+
+ Vector<double> res (4);
+ M.residual (res, solution, rhs);
+ deallog << "residual=" << res.l2_norm()
+ << std::endl;
+}
+
--- /dev/null
+
+DEAL::OK
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.0001279
+DEAL:psor:Bicgstab::Convergence step 4 value 0.0002656
DEAL:psor:GMRES::Starting value 1.467
DEAL:psor:GMRES::Convergence step 5 value 0.0006649
DEAL:psor:GMRES::Starting value 1.467
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::Failure step 10 value 0.001961
-DEAL:no-fail::Iterative method reported convergence failure in step 10 with residual 0.00196051
+DEAL:no-fail:Bicgstab::Convergence step 10 value 0.0002572
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
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 11 value 0.0003990
+DEAL:no:Bicgstab::Convergence step 10 value 0.0002572
DEAL:no:GMRES::Starting value 11.00
DEAL:no:GMRES::Convergence step 43 value 0.0009788
DEAL:no:GMRES::Starting value 11.00
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 11 value 0.0003990
+DEAL:rich:Bicgstab::Convergence step 10 value 0.0002572
DEAL:rich:GMRES::Starting value 6.600
DEAL:rich:GMRES::Convergence step 42 value 0.0007803
DEAL:rich:GMRES::Starting value 6.600
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 5 value 9.315e-05
+DEAL:ssor:Bicgstab::Convergence step 4 value 0.0003834
DEAL:ssor:GMRES::Starting value 11.05
DEAL:ssor:GMRES::Convergence step 7 value 0.0006838
DEAL:ssor:GMRES::Starting value 11.05
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 14 value 0.0009987
+DEAL:sor:Bicgstab::Convergence step 12 value 0.0007246
DEAL:sor:GMRES::Starting value 7.322
DEAL:sor:GMRES::Convergence step 23 value 0.0006981
DEAL:sor:GMRES::Starting value 7.322
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 11 value 0.0001588
+DEAL:psor:Bicgstab::Convergence step 10 value 0.0003167
DEAL:psor:GMRES::Starting value 7.345
DEAL:psor:GMRES::Convergence step 20 value 0.0008491
DEAL:psor:GMRES::Starting value 7.345