]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix an issue with the Chebyshev preconditioner.
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 29 Oct 2013 11:03:34 +0000 (11:03 +0000)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 29 Oct 2013 11:03:34 +0000 (11:03 +0000)
git-svn-id: https://svn.dealii.org/trunk@31479 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/lac/precondition.h
tests/lac/precondition_chebyshev_01.cc [new file with mode: 0644]
tests/lac/precondition_chebyshev_01.output [new file with mode: 0644]

index d2e059fd366e92a2a6151d4e316d8c57e492df40..598152b2f37c0b31c7b371c87282b190cc21dcdb 100644 (file)
@@ -1818,11 +1818,16 @@ PreconditionChebyshev<MATRIX,VECTOR>::initialize (const MATRIX &matrix,
           cg_message.erase(0, pos+5);
           std::string first = cg_message;
 
-          first.erase(cg_message.find_first_of(" "), std::string::npos);
+          if (cg_message.find_first_of(" ") != std::string::npos)
+            first.erase(cg_message.find_first_of(" "), std::string::npos);
           std::istringstream(first)      >> min_eigenvalue;
 
-          cg_message.erase(0, cg_message.find_last_of(" ")+1);
-          std::istringstream(cg_message) >> max_eigenvalue;
+          if (cg_message.find_last_of(" ") != std::string::npos)
+            {
+              cg_message.erase(0, cg_message.find_last_of(" ")+1);
+              std::istringstream(cg_message) >> max_eigenvalue;
+            }
+          else max_eigenvalue = min_eigenvalue;
         }
       else
         min_eigenvalue = max_eigenvalue = 1;
@@ -1910,7 +1915,7 @@ PreconditionChebyshev<MATRIX,VECTOR>::Tvmult (VECTOR &dst,
     (src, data.matrix_diagonal_inverse, true, 0., 1./theta, update1,
      update2, dst);
 
-  for (unsigned int k=0; k<data.degree-1; ++k)
+  for (unsigned int k=0; k<data.degree; ++k)
     {
       matrix_ptr->Tvmult (update2, dst);
       const double rhokp = 1./(2.*sigma-rhok);
diff --git a/tests/lac/precondition_chebyshev_01.cc b/tests/lac/precondition_chebyshev_01.cc
new file mode 100644 (file)
index 0000000..215456a
--- /dev/null
@@ -0,0 +1,111 @@
+// ---------------------------------------------------------------------
+// $Id$
+//
+// Copyright (C) 2005 - 2013 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+// Tests PreconditionChebyshev::vmult and PreconditionChebyshev::Tvmult
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/lac/precondition.h>
+#include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/vector.h>
+
+#include <fstream>
+#include <iomanip>
+#include <iomanip>
+#include <cmath>
+
+
+class FullMatrixModified : public FullMatrix<double>
+{
+public:
+  FullMatrixModified (unsigned int size1, unsigned int size2)
+    :
+    FullMatrix<double>(size1,size2)
+  {}
+
+  double el(unsigned int i, unsigned int j) const
+  {
+    return this->operator()(i,j);
+  }
+};
+
+
+void
+check()
+{
+  const unsigned int size = 10;
+  FullMatrixModified m(size,size);
+  for (unsigned int i=0; i<size; ++i)
+    m(i,i) = i+1;
+
+  Vector<double> in (size), out(size);
+  for (unsigned int i=0; i<size; ++i)
+    in(i) = (double)rand()/RAND_MAX;
+
+  PreconditionChebyshev<FullMatrixModified,Vector<double> > prec;
+  typename PreconditionChebyshev<FullMatrixModified,Vector<double> >::AdditionalData
+    data;
+  data.smoothing_range = 2 * size;
+  data.degree = 3;
+  prec.initialize(m, data);
+
+  deallog << "Check  vmult orig: ";
+  prec.vmult(out, in);
+  for (unsigned int i=0; i<size; ++i)
+    deallog << out(i) << " ";
+  deallog << std::endl;
+
+  deallog << "Check Tvmult orig: ";
+  prec.Tvmult(out, in);
+  for (unsigned int i=0; i<size; ++i)
+    deallog << out(i) << " ";
+  deallog << std::endl;
+
+  Vector<double> matrix_diagonal(size);
+  matrix_diagonal = 1;
+  data.matrix_diagonal_inverse = matrix_diagonal;
+  prec.initialize(m, data);
+
+  deallog << "Check  vmult diag: ";
+  prec.vmult(out, in);
+  for (unsigned int i=0; i<size; ++i)
+    deallog << out(i) << " ";
+  deallog << std::endl;
+
+  deallog << "Check Tvmult diag: ";
+  prec.Tvmult(out, in);
+  for (unsigned int i=0; i<size; ++i)
+    deallog << out(i) << " ";
+  deallog << std::endl;
+
+}
+
+
+int main()
+{
+  std::ofstream logfile("output");
+  deallog << std::fixed;
+  deallog << std::setprecision(0);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  check();
+
+  return 0;
+}
diff --git a/tests/lac/precondition_chebyshev_01.output b/tests/lac/precondition_chebyshev_01.output
new file mode 100644 (file)
index 0000000..3fcd646
--- /dev/null
@@ -0,0 +1,5 @@
+
+DEAL::Check  vmult orig: 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 
+DEAL::Check Tvmult orig: 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 
+DEAL::Check  vmult diag: 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 
+DEAL::Check Tvmult diag: 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.