]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add Assert's to some quick tests 3231/head
authorDenis Davydov <davydden@gmail.com>
Thu, 13 Oct 2016 17:49:53 +0000 (19:49 +0200)
committerDenis Davydov <davydden@gmail.com>
Thu, 13 Oct 2016 17:49:53 +0000 (19:49 +0200)
tests/quick_tests/arpack.cc
tests/quick_tests/step-petsc.cc
tests/quick_tests/step-slepc.cc
tests/quick_tests/step-trilinos.cc

index 087d061d275ed95a7bdb4847bf552b36ab528838..4eda317b6b0527999010d0aa7533e43567f4bd16 100644 (file)
@@ -93,6 +93,30 @@ void test ()
                      eigenvalues,
                      eigenvectors,
                      eigenvalues.size());
+
+  {
+    const double precision = 1e-7;
+    Vector<double> Ax(eigenvectors[0]), Bx(eigenvectors[0]);
+    for (unsigned int i=0; i < eigenvectors.size(); ++i)
+      {
+        B.vmult(Bx,eigenvectors[i]);
+
+        for (unsigned int j=0; j < eigenvectors.size(); j++)
+          Assert( std::abs( eigenvectors[j] * Bx - (i==j))< precision,
+                  ExcMessage("Eigenvectors " +
+                             Utilities::int_to_string(i) +
+                             " and " +
+                             Utilities::int_to_string(j) +
+                             " are not orthonormal!"));
+
+        A.vmult(Ax,eigenvectors[i]);
+        Ax.add(-1.0*std::real(eigenvalues[i]),Bx);
+        Assert (Ax.l2_norm() < precision,
+                ExcMessage("Returned vector " +
+                           Utilities::int_to_string(i) +
+                           " is not an eigenvector!"));
+      }
+  }
 }
 
 
index 61470331a82be2ab5757e57ca9295a7be1d20e78..77ee6929a77d2eb0f550a2c30ffc3d3dd7fd710b 100644 (file)
@@ -155,8 +155,10 @@ void LaplaceProblem::solve ()
   PETScWrappers::PreconditionBlockJacobi preconditioner (A);\r
   cg_solver.solve (A, x, b, preconditioner);\r
 \r
-  // some output\r
-  // ?\r
+  PETScWrappers::Vector res(x);\r
+  A.residual(res,x,b);\r
+  AssertThrow(res.l2_norm()<1e-3,\r
+              ExcInternalError());\r
 }\r
 \r
 void LaplaceProblem::run ()\r
index 42f9d792430cc17c569b4e512aa0e4c407f02054..c1469519f5936c1506b6e2720cfdac546eb39ca3 100644 (file)
@@ -154,11 +154,37 @@ void LaplaceEigenspectrumProblem::assemble_system ()
 \r
 void LaplaceEigenspectrumProblem::solve ()\r
 {\r
-  SolverControl solver_control (1000, 1e-03);\r
+  SolverControl solver_control (1000, 1e-10);\r
   SLEPcWrappers::SolverArnoldi eigensolver (solver_control);\r
   eigensolver.set_which_eigenpairs (EPS_SMALLEST_REAL);\r
   eigensolver.solve (A, B, lambda, x, x.size());\r
 \r
+  {\r
+    const double precision = 1e-7;\r
+    PETScWrappers::Vector Ax(x[0]), Bx(x[0]);\r
+    for (unsigned int i=0; i < x.size(); ++i)\r
+      {\r
+        B.vmult(Bx,x[i]);\r
+\r
+        for (unsigned int j=0; j < x.size(); j++)\r
+          if (j!=i)\r
+            Assert( std::abs( x[j] * Bx )< precision,\r
+                    ExcMessage("Eigenvectors " +\r
+                               Utilities::int_to_string(i) +\r
+                               " and " +\r
+                               Utilities::int_to_string(j) +\r
+                               " are not orthogonal!"));\r
+\r
+        A.vmult(Ax,x[i]);\r
+        Ax.add(-1.0*lambda[i],Bx);\r
+        Assert (Ax.l2_norm() < precision,\r
+                ExcMessage("Returned vector " +\r
+                           Utilities::int_to_string(i) +\r
+                           " is not an eigenvector!"));\r
+      }\r
+  }\r
+\r
+\r
   // some output\r
   output_table.add_value ("lambda", lambda[0]);\r
   output_table.add_value ("error", std::fabs(2.-lambda[0]));\r
index df87c68ddb6f91cde6fb62c37133c1e37bba26e6..b145721504f9e0bd55a641e26da236ede3273749 100644 (file)
@@ -159,8 +159,10 @@ void LaplaceProblem::solve ()
   preconditioner.initialize (A);\r
   cg_solver.solve (A, x, b, preconditioner);\r
 \r
-  // some output\r
-  // ?\r
+  TrilinosWrappers::Vector res(x);\r
+  A.residual(res,x,b);\r
+  AssertThrow(res.l2_norm()<1e-3,\r
+              ExcInternalError());\r
 }\r
 \r
 void LaplaceProblem::run ()\r

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.