]> https://gitweb.dealii.org/ - dealii.git/commitdiff
And again the fix for the project_bv_curl_conf_04 test 9000/head
authorWinnifried Wollner <wollner@mathematik.tu-darmstadt.de>
Mon, 9 Mar 2020 13:18:33 +0000 (14:18 +0100)
committerWolfgang Bangerth <bangerth@colostate.edu>
Sun, 29 Mar 2020 16:03:35 +0000 (10:03 -0600)
doc/news/changes/minor/20200228Wollner
include/deal.II/numerics/vector_tools.templates.h
tests/numerics/project_bv_curl_conf_04.cc
tests/numerics/project_bv_curl_conf_04.output

index 3c1eda51c803f9e207f00e137bbed80dd4133ed3..989cd54a2583fcbf35e168f87543d0d09c84b535 100644 (file)
@@ -1,6 +1,4 @@
 Changed: Fixed VectorTools::project_boundary_values_curl_conforming_l2 to
 work with FESystems containing an FE_Nedelec.
 <br>
-
 (Winnifried Wollner, 2020/02/28)
-
index 5e5ff781d36d684dd0c8eec5ed35e7ef648be70c..73b53f30cbb59ad91a66f20bd930a913efc33e5f 100644 (file)
@@ -5442,7 +5442,7 @@ namespace VectorTools
       //
       // We call the map associated_edge_dof_to_face_dof
       std::vector<unsigned int> associated_edge_dof_to_face_dof(
-        degree + 1, numbers::invalid_dof_index);
+        degree + 1, numbers::invalid_unsigned_int);
 
       // Lowest DoF in the base element allowed for this edge:
       const unsigned int lower_bound =
@@ -5477,7 +5477,7 @@ namespace VectorTools
           // Note, assuming that the edge orientations are "standard"
           //       i.e. cell->line_orientation(line) = true.
           Assert(cell->line_orientation(line),
-                 ExcMessage("Edge orientation doesnot meet expectation."));
+                 ExcMessage("Edge orientation does not meet expectation."));
           // Next, translate from face to cell. Note, this might be assuming
           // that the edge orientations are "standard" (not sure any more at
           // this time), i.e.
@@ -5994,7 +5994,7 @@ namespace VectorTools
                       AssertIndexRange(associated_face_dof_index,
                                        associated_face_dof_to_face_dof.size());
                       associated_face_dof_to_face_dof
-                        [associated_face_dof_index] = quad_dof_idx;
+                        [associated_face_dof_index] = face_idx;
                       ++associated_face_dof_index;
                     }
                 }
index 8f93314a55c5ba9dc9550616a414a17e726a2034..559ce11db9b40cadf0954f66f025e81923d94694 100644 (file)
@@ -1,3 +1,23 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2014 - 2019 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// Checks that project_boundary_values_curl_conforming
+// works for FESystems with FE_Nedelec and FE_Q mixed
+
+
 #include <deal.II/dofs/dof_handler.h>
 
 #include <deal.II/fe/fe_nedelec.h>
@@ -11,7 +31,9 @@
 #include <deal.II/lac/affine_constraints.h>
 
 #include <deal.II/numerics/vector_tools.h>
-using namespace dealii;
+
+#include "../tests.h"
+
 
 const unsigned int dim = 3;
 
@@ -24,7 +46,7 @@ apply_boundary_values(DoFHandler<dim> &          dof_handler,
                       Vector<double> &           dst)
 {
   constraints.clear();
-  for (unsigned int i = 0; i < 6; i++)
+  for (unsigned int i = 0; i < 6; ++i)
     VectorTools::project_boundary_values_curl_conforming_l2(
       dof_handler,
       start_comp,
@@ -37,6 +59,8 @@ apply_boundary_values(DoFHandler<dim> &          dof_handler,
   constraints.distribute(dst);
 }
 
+
+
 bool
 test_boundary_values(DoFHandler<dim> &   dof_handler,
                      Mapping<dim> &      mapping,
@@ -47,24 +71,24 @@ test_boundary_values(DoFHandler<dim> &   dof_handler,
 {
   double ret = true;
   // Initialize
-  QGaussLobatto<dim - 1>                         quadrature(3);
-  FEFaceValues<dim>                              fe_values(mapping,
+  QGaussLobatto<dim - 1> quadrature(3);
+  FEFaceValues<dim>      fe_values(mapping,
                               fe,
                               quadrature,
                               update_values | update_normal_vectors |
                                 update_quadrature_points | update_JxW_values);
-  typename DoFHandler<dim>::active_cell_iterator cell =
-                                                   dof_handler.begin_active(),
-                                                 endc = dof_handler.end();
-  unsigned                    num_cells               = 0;
+
+  unsigned                    num_cells = 0;
   std::vector<Vector<double>> local_averages(
     6, Vector<double>(dof_handler.n_dofs()));
   const FEValuesExtractors::Vector nedelec(start_comp);
 
   // For testing we take only one cell!
-  assert(dof_handler.n_dofs() == cell->get_fe().dofs_per_cell);
+  Assert(dof_handler.n_dofs() ==
+           dof_handler.begin_active()->get_fe().dofs_per_cell,
+         ExcInternalError());
 
-  for (; cell != endc; ++cell)
+  for (const auto &cell : dof_handler.active_cell_iterators())
     {
       // For testing we take only one cell! Otherwise, local to global is
       // missing
@@ -114,15 +138,15 @@ test_boundary_values(DoFHandler<dim> &   dof_handler,
       num_cells++;
     }
   // Test if ok
-  for (unsigned int bc = 0; bc < 6; bc++)
+  for (unsigned int bc = 0; bc < 6; ++bc)
     {
       for (unsigned int basis = 0; basis < dof_handler.n_dofs(); basis++)
         {
           if (fabs(local_averages[bc](basis)) > 1.e-10)
             {
-              std::cout << "Wrong values on boundary " << bc
-                        << " in basis function " << basis << " of size "
-                        << fabs(local_averages[bc](basis)) << std::endl;
+              deallog << "Wrong values on boundary " << bc
+                      << " in basis function " << basis << " of size "
+                      << fabs(local_averages[bc](basis)) << std::endl;
               ret = false;
             }
         }
@@ -130,9 +154,12 @@ test_boundary_values(DoFHandler<dim> &   dof_handler,
   return ret;
 }
 
+
 int
 main(int /*argc*/, char ** /*argv*/)
 {
+  initlog();
+
   Triangulation<dim> triangulation;
   GridGenerator::hyper_cube(triangulation, 0, 1, true);
   MappingQ<dim>             mapping(1);
@@ -142,7 +169,7 @@ main(int /*argc*/, char ** /*argv*/)
 
   {
     // Test 1: Only FE_Nedelec
-    std::cout << "Testing FE_Nedelec(0): ";
+    deallog << "Testing FE_Nedelec(0): ";
     FE_Nedelec<dim> fe(0);
     DoFHandler<dim> dof_handler(triangulation);
     dof_handler.distribute_dofs(fe);
@@ -150,17 +177,17 @@ main(int /*argc*/, char ** /*argv*/)
     apply_boundary_values(dof_handler, constraints, 3, 0, mapping, test_vec2);
     if (test_boundary_values(dof_handler, mapping, fe, 3, 0, test_vec2))
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
   }
   {
     // Test 2: FESystem(FE_Nedelec(0))
-    std::cout << "Testing FESystem(FE_Nedelec(0)): ";
+    deallog << "Testing FESystem(FE_Nedelec(0)): ";
     FESystem<dim>   fe_system(FE_Nedelec<dim>(0), 1);
     DoFHandler<dim> dof_handler(triangulation);
     dof_handler.distribute_dofs(fe_system);
@@ -168,16 +195,16 @@ main(int /*argc*/, char ** /*argv*/)
     apply_boundary_values(dof_handler, constraints, 3, 0, mapping, test_vec);
     if (test_boundary_values(dof_handler, mapping, fe_system, 3, 0, test_vec))
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
     // Now, in the case of only one element both initializations should give
     // identical vectors
-    std::cout << "Checking Consistency: ";
+    deallog << "Checking Consistency: ";
     bool equal = (test_vec.size() == test_vec2.size());
     if (equal)
       {
@@ -188,17 +215,17 @@ main(int /*argc*/, char ** /*argv*/)
       }
     if (equal)
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
   }
   {
     // Test 1b: Only FE_Nedelec
-    std::cout << "Testing FE_Nedelec(1): ";
+    deallog << "Testing FE_Nedelec(1): ";
     FE_Nedelec<dim> fe(1);
     DoFHandler<dim> dof_handler(triangulation);
     dof_handler.distribute_dofs(fe);
@@ -206,17 +233,17 @@ main(int /*argc*/, char ** /*argv*/)
     apply_boundary_values(dof_handler, constraints, 3, 0, mapping, test_vec2);
     if (test_boundary_values(dof_handler, mapping, fe, 3, 0, test_vec2))
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
   }
   {
     // Test 2b: FESystem(FE_Nedelec(1))
-    std::cout << "Testing FESystem(FE_Nedelec(1)): ";
+    deallog << "Testing FESystem(FE_Nedelec(1)): ";
     FESystem<dim>   fe_system(FE_Nedelec<dim>(1), 1);
     DoFHandler<dim> dof_handler(triangulation);
     dof_handler.distribute_dofs(fe_system);
@@ -224,16 +251,16 @@ main(int /*argc*/, char ** /*argv*/)
     apply_boundary_values(dof_handler, constraints, 3, 0, mapping, test_vec);
     if (test_boundary_values(dof_handler, mapping, fe_system, 3, 0, test_vec))
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
     // Now, in the case of only one element both initializations should give
     // identical vectors
-    std::cout << "Checking Consistency: ";
+    deallog << "Checking Consistency: ";
     bool equal = (test_vec.size() == test_vec2.size());
     if (equal)
       {
@@ -244,17 +271,17 @@ main(int /*argc*/, char ** /*argv*/)
       }
     if (equal)
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
   }
   {
     // Test 3: FESystem(FE_Nedelec(0), FE_Nedelec(0))
-    std::cout << "Testing FESystem(FE_Nedelec(0), FE_Nedelec(0)): ";
+    deallog << "Testing FESystem(FE_Nedelec(0), FE_Nedelec(0)): ";
     FESystem<dim>   fe_system(FE_Nedelec<dim>(0), 1, FE_Nedelec<dim>(0), 1);
     DoFHandler<dim> dof_handler(triangulation);
     dof_handler.distribute_dofs(fe_system);
@@ -262,17 +289,17 @@ main(int /*argc*/, char ** /*argv*/)
     apply_boundary_values(dof_handler, constraints, 6, 3, mapping, test_vec);
     if (test_boundary_values(dof_handler, mapping, fe_system, 6, 3, test_vec))
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
   }
   {
     // Test 4: FESystem(FE_Nedelec(1), FE_Nedelec(0))
-    std::cout << "Testing FESystem(FE_Nedelec(1), FE_Nedelec(0)): ";
+    deallog << "Testing FESystem(FE_Nedelec(1), FE_Nedelec(0)): ";
     FESystem<dim>   fe_system(FE_Nedelec<dim>(1), 1, FE_Nedelec<dim>(0), 1);
     DoFHandler<dim> dof_handler(triangulation);
     dof_handler.distribute_dofs(fe_system);
@@ -280,17 +307,17 @@ main(int /*argc*/, char ** /*argv*/)
     apply_boundary_values(dof_handler, constraints, 6, 3, mapping, test_vec);
     if (test_boundary_values(dof_handler, mapping, fe_system, 6, 3, test_vec))
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
   }
   {
     // Test 5: FESystem(FE_Nedelec(0), FE_Nedelec(1))
-    std::cout << "Testing FESystem(FE_Nedelec(0), FE_Nedelec(1)): ";
+    deallog << "Testing FESystem(FE_Nedelec(0), FE_Nedelec(1)): ";
     FESystem<dim>   fe_system(FE_Nedelec<dim>(0), 1, FE_Nedelec<dim>(1), 1);
     DoFHandler<dim> dof_handler(triangulation);
     dof_handler.distribute_dofs(fe_system);
@@ -298,17 +325,17 @@ main(int /*argc*/, char ** /*argv*/)
     apply_boundary_values(dof_handler, constraints, 6, 3, mapping, test_vec);
     if (test_boundary_values(dof_handler, mapping, fe_system, 6, 3, test_vec))
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
   }
   {
     // Test 6: FESystem(FE_Nedelec(0), FE_Q(1))
-    std::cout << "Testing FESystem(FE_Nedelec(0), FE_Q(1)): ";
+    deallog << "Testing FESystem(FE_Nedelec(0), FE_Q(1)): ";
     FESystem<dim>   fe_system(FE_Nedelec<dim>(0), 1, FE_Q<dim>(1), 1);
     DoFHandler<dim> dof_handler(triangulation);
     dof_handler.distribute_dofs(fe_system);
@@ -316,17 +343,17 @@ main(int /*argc*/, char ** /*argv*/)
     apply_boundary_values(dof_handler, constraints, 4, 0, mapping, test_vec);
     if (test_boundary_values(dof_handler, mapping, fe_system, 4, 0, test_vec))
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
   }
   {
     // Test 7: FESystem(FE_Q(1), FE_Nedelec(0))
-    std::cout << "Testing FESystem(FE_Q(1), FE_Nedelec(0)): ";
+    deallog << "Testing FESystem(FE_Q(1), FE_Nedelec(0)): ";
     FESystem<dim>   fe_system(FE_Q<dim>(1), 1, FE_Nedelec<dim>(0), 1);
     DoFHandler<dim> dof_handler(triangulation);
     dof_handler.distribute_dofs(fe_system);
@@ -334,17 +361,17 @@ main(int /*argc*/, char ** /*argv*/)
     apply_boundary_values(dof_handler, constraints, 4, 1, mapping, test_vec);
     if (test_boundary_values(dof_handler, mapping, fe_system, 4, 1, test_vec))
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
   }
   {
     // Test 8: FESystem(FE_Nedelec(0), FE_Q(2))
-    std::cout << "Testing FESystem(FE_Nedelec(0), FE_Q(2)): ";
+    deallog << "Testing FESystem(FE_Nedelec(0), FE_Q(2)): ";
     FESystem<dim>   fe_system(FE_Nedelec<dim>(0), 1, FE_Q<dim>(2), 1);
     DoFHandler<dim> dof_handler(triangulation);
     dof_handler.distribute_dofs(fe_system);
@@ -352,17 +379,17 @@ main(int /*argc*/, char ** /*argv*/)
     apply_boundary_values(dof_handler, constraints, 4, 0, mapping, test_vec);
     if (test_boundary_values(dof_handler, mapping, fe_system, 4, 0, test_vec))
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
   }
   {
     // Test 9: FESystem(FE_Q(2), FE_Nedelec(0))
-    std::cout << "Testing FESystem(FE_Q(2), FE_Nedelec(0)): ";
+    deallog << "Testing FESystem(FE_Q(2), FE_Nedelec(0)): ";
     FESystem<dim>   fe_system(FE_Q<dim>(2), 1, FE_Nedelec<dim>(0), 1);
     DoFHandler<dim> dof_handler(triangulation);
     dof_handler.distribute_dofs(fe_system);
@@ -370,17 +397,17 @@ main(int /*argc*/, char ** /*argv*/)
     apply_boundary_values(dof_handler, constraints, 4, 1, mapping, test_vec);
     if (test_boundary_values(dof_handler, mapping, fe_system, 4, 1, test_vec))
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
   }
   {
     // Test 10: FESystem(FE_Nedelec(1), FE_Q(1))
-    std::cout << "Testing FESystem(FE_Nedelec(1), FE_Q(1)): ";
+    deallog << "Testing FESystem(FE_Nedelec(1), FE_Q(1)): ";
     FESystem<dim>   fe_system(FE_Nedelec<dim>(1), 1, FE_Q<dim>(1), 1);
     DoFHandler<dim> dof_handler(triangulation);
     dof_handler.distribute_dofs(fe_system);
@@ -388,17 +415,17 @@ main(int /*argc*/, char ** /*argv*/)
     apply_boundary_values(dof_handler, constraints, 4, 0, mapping, test_vec);
     if (test_boundary_values(dof_handler, mapping, fe_system, 4, 0, test_vec))
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
   }
   {
     // Test 11: FESystem(FE_Q(1), FE_Nedelec(1))
-    std::cout << "Testing FESystem(FE_Q(1), FE_Nedelec(1)): ";
+    deallog << "Testing FESystem(FE_Q(1), FE_Nedelec(1)): ";
     FESystem<dim>   fe_system(FE_Q<dim>(1), 1, FE_Nedelec<dim>(1), 1);
     DoFHandler<dim> dof_handler(triangulation);
     dof_handler.distribute_dofs(fe_system);
@@ -406,11 +433,11 @@ main(int /*argc*/, char ** /*argv*/)
     apply_boundary_values(dof_handler, constraints, 4, 1, mapping, test_vec);
     if (test_boundary_values(dof_handler, mapping, fe_system, 4, 1, test_vec))
       {
-        std::cout << "OK" << std::endl;
+        deallog << "OK" << std::endl;
       }
     else
       {
-        std::cout << "Failed" << std::endl;
+        deallog << "Failed" << std::endl;
         abort();
       }
   }
index b8d4974d031f4cbd4033388cf73ceb85d4e05011..245c9111eb11f9e97837b1c94f0402326bf1a8fe 100644 (file)
@@ -1,15 +1,16 @@
-Testing FE_Nedelec(0): OK
-Testing FESystem(FE_Nedelec(0)): OK
-Checking Consistency: OK
-Testing FE_Nedelec(1): OK
-Testing FESystem(FE_Nedelec(1)): OK
-Checking Consistency: OK
-Testing FESystem(FE_Nedelec(0), FE_Nedelec(0)): OK
-Testing FESystem(FE_Nedelec(1), FE_Nedelec(0)): OK
-Testing FESystem(FE_Nedelec(0), FE_Nedelec(1)): OK
-Testing FESystem(FE_Nedelec(0), FE_Q(1)): OK
-Testing FESystem(FE_Q(1), FE_Nedelec(0)): OK
-Testing FESystem(FE_Nedelec(0), FE_Q(2)): OK
-Testing FESystem(FE_Q(2), FE_Nedelec(0)): OK
-Testing FESystem(FE_Nedelec(1), FE_Q(1)): OK
-Testing FESystem(FE_Q(1), FE_Nedelec(1)): OK
+
+DEAL::Testing FE_Nedelec(0): OK
+DEAL::Testing FESystem(FE_Nedelec(0)): OK
+DEAL::Checking Consistency: OK
+DEAL::Testing FE_Nedelec(1): OK
+DEAL::Testing FESystem(FE_Nedelec(1)): OK
+DEAL::Checking Consistency: OK
+DEAL::Testing FESystem(FE_Nedelec(0), FE_Nedelec(0)): OK
+DEAL::Testing FESystem(FE_Nedelec(1), FE_Nedelec(0)): OK
+DEAL::Testing FESystem(FE_Nedelec(0), FE_Nedelec(1)): OK
+DEAL::Testing FESystem(FE_Nedelec(0), FE_Q(1)): OK
+DEAL::Testing FESystem(FE_Q(1), FE_Nedelec(0)): OK
+DEAL::Testing FESystem(FE_Nedelec(0), FE_Q(2)): OK
+DEAL::Testing FESystem(FE_Q(2), FE_Nedelec(0)): OK
+DEAL::Testing FESystem(FE_Nedelec(1), FE_Q(1)): OK
+DEAL::Testing FESystem(FE_Q(1), FE_Nedelec(1)): OK

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.