]> https://gitweb.dealii.org/ - dealii.git/commitdiff
tweak PETSc tests 15767/head
authorStefano Zampini <stefano.zampini@gmail.com>
Wed, 19 Jul 2023 17:53:21 +0000 (19:53 +0200)
committerStefano Zampini <stefano.zampini@gmail.com>
Wed, 19 Jul 2023 18:03:48 +0000 (20:03 +0200)
tests/petsc/assemble_matrix_parallel_01.cc
tests/petsc/step-77-snes.with_p4est=true.with_petsc_with_mumps=true.mpirun=2.output [moved from tests/petsc/step-77-snes.with_p4est=true.mpirun=2.output with 100% similarity]

index 8ac2b5fbffa9c4efeb56919125c4dd6bda9440bf..714d4f67e7af94ab79b3ea9d0ff6307faf6abca9 100644 (file)
 // same as deal.II/assemble_matrix_parallel_01, but for PETSc matrices
 // and vectors
 //
-// This test requires PETSc to be configured with the option
-// "--with-threadsafety" in case of debug builds of PETSc
-// For optimized builds, the above option is needed only in case
-// users want PETSc to produce useful logs with "-log_view" runtime
-// option
 #include <deal.II/base/function.h>
 #include <deal.II/base/graph_coloring.h>
 #include <deal.II/base/quadrature_lib.h>
@@ -53,6 +48,8 @@
 #include <deal.II/numerics/matrix_tools.h>
 #include <deal.II/numerics/vector_tools.h>
 
+#include <petscconf.h>
+
 #include <complex>
 #include <iostream>
 
@@ -384,6 +381,12 @@ template <int dim>
 void
 LaplaceProblem<dim>::assemble_test()
 {
+  // This test requires PETSc to be configured with the option
+  // "--with-threadsafety" in case of debug builds of PETSc.
+  // For optimized builds, the above option is needed only in case
+  // users want PETSc to produce useful logs with "-log_view" runtime
+  // option
+#if !defined(PETSC_USE_DEBUG) || defined(PETSC_HAVE_THREADSAFETY)
   test_matrix = 0;
   test_rhs    = 0;
 
@@ -405,11 +408,18 @@ LaplaceProblem<dim>::assemble_test()
   test_rhs.compress(VectorOperation::add);
 
   test_matrix.add(-1, reference_matrix);
-
-  // there should not even be roundoff difference between matrices
-  deallog << "error in matrix: " << test_matrix.frobenius_norm() << std::endl;
   test_rhs.add(-1., reference_rhs);
-  deallog << "error in vector: " << test_rhs.l2_norm() << std::endl;
+
+  auto mnorm = test_matrix.frobenius_norm();
+  auto vnorm = test_rhs.l2_norm();
+#else
+  auto mnorm = 0.0;
+  auto vnorm = 0.0;
+#endif
+
+  // there should not even be roundoff difference
+  deallog << "error in matrix: " << mnorm << std::endl;
+  deallog << "error in vector: " << vnorm << std::endl;
 }
 
 

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.