]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Tests: Fix petsc_complex/solver_real_(01|02) to test for convergence
authorMatthias Maier <tamiko@43-1.org>
Tue, 12 Jul 2016 15:34:41 +0000 (10:34 -0500)
committerMatthias Maier <tamiko@43-1.org>
Tue, 12 Jul 2016 15:34:41 +0000 (10:34 -0500)
tests/petsc_complex/solver_real_01.cc
tests/petsc_complex/solver_real_01.output
tests/petsc_complex/solver_real_02.cc
tests/petsc_complex/solver_real_02.output

index 383e017417ad86814ce06c7552cbfc20a8e443bf..5d1564580a98c5e42570aa50ad0b357fef678edc 100644 (file)
 
 // Note: This is (almost) a clone of the tests/petsc/solver_01.cc
 
-
 #include "../tests.h"
 #include "../lac/testmatrix.h"
-#include <cmath>
-#include <fstream>
-#include <iostream>
-#include <iomanip>
-#include <deal.II/base/logstream.h>
+
+#include <deal.II/lac/petsc_precondition.h>
+#include <deal.II/lac/petsc_solver.h>
 #include <deal.II/lac/petsc_sparse_matrix.h>
 #include <deal.II/lac/petsc_vector.h>
-#include <deal.II/lac/petsc_solver.h>
-#include <deal.II/lac/petsc_precondition.h>
-#include <deal.II/lac/vector_memory.h>
-#include <typeinfo>
-
-template<class SOLVER, class MATRIX, class VECTOR, class PRECONDITION>
-void
-check_solve( SOLVER &solver, const MATRIX &A,
-             VECTOR &u, VECTOR &f, const PRECONDITION &P)
-{
-  deallog << "Solver type: " << typeid(solver).name() << std::endl;
-
-  u = 0.;
-  f = 1.;
-  try
-    {
-      solver.solve(A,u,f,P);
-    }
-  catch (dealii::SolverControl::NoConvergence &e)
-    {
-      deallog << "Exception: " << e.get_exc_name() << std::endl;
-    }
-
-  deallog << "Solver stopped after " << solver.control().last_step()
-          << " iterations" << std::endl;
-}
-
 
 int main(int argc, char **argv)
 {
-  std::ofstream logfile("output");
-  deallog.attach(logfile);
-  deallog << std::setprecision(4);
-  deallog.depth_console(0);
-  deallog.threshold_double(1.e-10);
+  initlog();
 
   Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
   {
-    SolverControl control(100, 1.e-3);
-
     const unsigned int size = 32;
     unsigned int dim = (size-1)*(size-1);
 
@@ -80,13 +44,24 @@ int main(int argc, char **argv)
 
     PETScWrappers::Vector  f(dim);
     PETScWrappers::Vector  u(dim);
+
     f = 1.;
+    u = 0.;
+
     A.compress (VectorOperation::insert);
     f.compress (VectorOperation::insert);
     u.compress (VectorOperation::insert);
 
+    // Richardson is a tricky smoother for the kind of FD matrix we use in
+    // this test. So, simply test that we're able to reduce the residual to
+    // a reasonably small value of 1.e-4.
+    SolverControl control(2500, 1.e-4);
+
     PETScWrappers::SolverRichardson solver(control);
     PETScWrappers::PreconditionJacobi preconditioner(A);
-    check_solve (solver, A,u,f, preconditioner);
+
+    check_solver_within_range(
+      solver.solve(A,u,f, preconditioner),
+      control.last_step(), 2295, 2300);
   }
 }
index dbab15846c21d4adc218d78c53b97ac68bb0b1e8..60a06149f0ca766a1cad3e9061a3cae514b02c2a 100644 (file)
@@ -1,7 +1,3 @@
 
 DEAL::Size 32 Unknowns 961
-DEAL::Solver type: N6dealii13PETScWrappers16SolverRichardsonE
-DEAL::Starting value 7.750
-DEAL::Failure step 100 value 4.004
-DEAL::Exception: SolverControl::NoConvergence (solver_control.last_step(), solver_control.last_value())
-DEAL::Solver stopped after 100 iterations
+DEAL::Solver stopped within 2295 - 2300 iterations
index 8536e318f16c713b55ebf1013e313c5659afa54a..d8b1040b02e6bd50627000d0e8b186b281eb93f0 100644 (file)
 
 #include "../tests.h"
 #include "../lac/testmatrix.h"
-#include <cmath>
-#include <fstream>
-#include <iostream>
-#include <iomanip>
-#include <deal.II/base/logstream.h>
+
+#include <deal.II/lac/petsc_precondition.h>
+#include <deal.II/lac/petsc_solver.h>
 #include <deal.II/lac/petsc_sparse_matrix.h>
 #include <deal.II/lac/petsc_vector.h>
-#include <deal.II/lac/petsc_solver.h>
-#include <deal.II/lac/petsc_precondition.h>
-#include <deal.II/lac/vector_memory.h>
-#include <typeinfo>
-
-template<class SOLVER, class MATRIX, class VECTOR, class PRECONDITION>
-void
-check_solve( SOLVER &solver, const MATRIX &A,
-             VECTOR &u, VECTOR &f, const PRECONDITION &P)
-{
-  deallog << "Solver type: " << typeid(solver).name() << std::endl;
-
-  u = 0.;
-  f = 1.;
-  try
-    {
-      solver.solve(A,u,f,P);
-    }
-  catch (dealii::SolverControl::NoConvergence &e)
-    {
-      deallog << "Exception: " << e.get_exc_name() << std::endl;
-    }
-
-  deallog << "Solver stopped after " << solver.control().last_step()
-          << " iterations" << std::endl;
-}
-
 
 int main(int argc, char **argv)
 {
-  std::ofstream logfile("output");
-  deallog.attach(logfile);
-  deallog << std::setprecision(4);
-  deallog.depth_console(0);
-  deallog.threshold_double(1.e-10);
+  initlog();
 
   Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
   {
-    SolverControl control(100, 1.e-3);
-
     const unsigned int size = 32;
     unsigned int dim = (size-1)*(size-1);
 
@@ -79,14 +44,24 @@ int main(int argc, char **argv)
 
     PETScWrappers::Vector  f(dim);
     PETScWrappers::Vector  u(dim);
+
     f = 1.;
+    u = 0.;
+
     A.compress (VectorOperation::insert);
     f.compress (VectorOperation::insert);
     u.compress (VectorOperation::insert);
 
+    // Chebychev is a tricky smoother for the kind of FD matrix we use in
+    // this test. So, simply test that we're able to reduce the residual to
+    // a reasonably small value of 1.e-3.
+    SolverControl control(2500, 1.e-3);
+
     PETScWrappers::SolverChebychev solver(control);
     PETScWrappers::PreconditionJacobi preconditioner(A);
-    check_solve (solver, A,u,f, preconditioner);
-  }
 
+    check_solver_within_range(
+      solver.solve(A,u,f, preconditioner),
+      control.last_step(), 1120, 1125);
+  }
 }
index 00c1e2226c669c55b2cdfe4fbf09011a659f140f..8e5ffc60890ccf502e3b24fcb8a223ec952eda33 100644 (file)
@@ -1,7 +1,3 @@
 
 DEAL::Size 32 Unknowns 961
-DEAL::Solver type: N6dealii13PETScWrappers15SolverChebychevE
-DEAL::Starting value 7.551
-DEAL::Failure step 100 value 2.937
-DEAL::Exception: SolverControl::NoConvergence (solver_control.last_step(), solver_control.last_value())
-DEAL::Solver stopped after 100 iterations
+DEAL::Solver stopped within 1120 - 1125 iterations

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.