Fixed a bug in the Nedelec elements

git-svn-id: https://svn.dealii.org/trunk@22820 0785d39b-7218-0410-832d-ea1e28bc413d

diff --git a/deal.II/include/deal.II/numerics/vectors.templates.h b/deal.II/include/deal.II/numerics/vectors.templates.h
index f1340a2929b946528d70d64122de9190a63ed315..8bd2b0d6126b719c83960210d0b2fb2b9f5bda2f 100644 (file)
--- a/deal.II/include/deal.II/numerics/vectors.templates.h
+++ b/deal.II/include/deal.II/numerics/vectors.templates.h
@@ -2882,19 +2882,11 @@ namespace internals {
               global_face_coordinate_directions[GeometryInfo<3>::faces_per_cell][2]
               = { { 1, 2 },
                   { 1, 2 },
-                  { 0, 2 },
-                  { 0, 2 },
+                  { 2, 0 },
+                  { 2, 0 },
                   { 0, 1 },
                   { 0, 1 } };
-            const unsigned int
-              local_face_coordinate_directions[GeometryInfo<3>::faces_per_cell][2]
-              = { { 1, 0 },
-                  { 1, 0 },
-                  { 0, 1 },
-                  { 0, 1 },
-                  { 1, 0 },
-                  { 1, 0 } };
-
+            
                                                             // The projection is
                                                             // divided into two steps.
                                                             // In the first step we
@@ -2933,10 +2925,10 @@ namespace internals {
 
                 for (unsigned int i = 0; i < 2; ++i)
                   for (unsigned int j = 0; j <= degree; ++j)
-                    tmp -= dof_values[(i + 2 * local_face_coordinate_directions[face][0]) * superdegree + j]
+                    tmp -= dof_values[(i + 2) * superdegree + j]
                            * fe_values[vec].value (cell->get_fe ().face_to_cell_index
-                                                   ((i + 2 * local_face_coordinate_directions[face][0])
-                                                    * superdegree + j, face), q_point);
+                                                   ((i + 2) * superdegree + j,
+                                                    face), q_point);
 
                 const double JxW
                   = std::sqrt (fe_values.JxW (q_point)
@@ -3034,10 +3026,9 @@ namespace internals {
                  for (unsigned int i = 0; i < 2; ++i)
                    for (unsigned int j = 0; j <= degree; ++j)
                      tmp
-                       -= dof_values[(i + 2 * local_face_coordinate_directions[face][1]) * superdegree + j]
+                       -= dof_values[i * superdegree + j]
                        * fe_values[vec].value (cell->get_fe ().face_to_cell_index
-                                               ((i + 2 * local_face_coordinate_directions[face][1])
-                                                * superdegree + j, face), q_point);
+                                               (i * superdegree + j, face), q_point);
 
                  const double JxW
                    = std::sqrt (fe_values.JxW (q_point)
@@ -3556,10 +3547,12 @@ project_boundary_values_curl_conforming (const hp::DoFHandler<dim>& dof_handler,
               constraints.add_line (dof);
               constraints.set_inhomogeneity (dof, computed_constraints[dof]);
             }
+        
+        break;
       }
 
       default:
-            Assert (false, ExcNotImplemented ());
+        Assert (false, ExcNotImplemented ());
     }
 }
 

diff --git a/deal.II/source/base/polynomials_nedelec.cc b/deal.II/source/base/polynomials_nedelec.cc
index 15baff767bc652434fa07b8c2956383b0ae91bf6..d8a9db0d67863869f94de8d84fed62cf3f68519d 100644 (file)
--- a/deal.II/source/base/polynomials_nedelec.cc
+++ b/deal.II/source/base/polynomials_nedelec.cc
@@ -432,12 +432,12 @@ const
                                  * my_degree + j
                                  + GeometryInfo<dim>::lines_per_cell][2]
                             = p2_values[i + (j + k * (my_degree + 2) + 2)
-                                        * (my_degree + 1)];
+                                          * (my_degree + 1)]; 
                           values[i + (j + (2 * k + 5) * my_degree
                                       + GeometryInfo<dim>::lines_per_cell)
                                  * (my_degree + 1)][0]
-                            = unit_point_values[i + ((j + 2) * (my_degree + 2)
-                                                     + k) * (my_degree + 1)];
+                            = unit_point_values[i + ((j + 2) * (my_degree + 2) + k)
+                                                  * (my_degree + 1)];
                           values[(i + 2 * (k + 4) * (my_degree + 1)
                                   + GeometryInfo<dim>::lines_per_cell)
                                  * my_degree + j
@@ -1085,7 +1085,7 @@ const
                                                     + GeometryInfo<dim>::lines_per_cell)
                                                    * my_degree + j
                                                    + GeometryInfo<dim>::lines_per_cell]
-                                                  [o + n][l][m] = 0.0;
+                                                  [o + k][l][m] = 0.0;
                                       }
 
                                     grad_grads[(i + 2 * k * (my_degree + 1)

diff --git a/deal.II/source/fe/fe_nedelec.cc b/deal.II/source/fe/fe_nedelec.cc
index f744500769673bb0cf5a2f659339e7319057e200..323ead2db4b626cc779717e2a000232aef0234ab 100644 (file)
--- a/deal.II/source/fe/fe_nedelec.cc
+++ b/deal.II/source/fe/fe_nedelec.cc
@@ -2064,16 +2064,16 @@ FE_Nedelec<dim>::has_support_on_face (const unsigned int shape_index,
                       >= (GeometryInfo<3>::lines_per_cell + 2 * deg)
                          * this->degree) &&
                    (shape_index
-                      < (GeometryInfo<3>::lines_per_cell + 5 * deg)
+                      < (GeometryInfo<3>::lines_per_cell + 4 * deg)
                         * this->degree)) ||
                   ((shape_index
-                      >= (GeometryInfo<3>::lines_per_cell + 6 * deg)
+                      >= (GeometryInfo<3>::lines_per_cell + 5 * deg)
                          * this->degree) &&
                    (shape_index
-                      < (GeometryInfo<3>::lines_per_cell + 7 * deg)
+                      < (GeometryInfo<3>::lines_per_cell + 6 * deg)
                         * this->degree)) ||
                   ((shape_index
-                      >= (GeometryInfo<3>::lines_per_cell + 8 * deg)
+                      >= (GeometryInfo<3>::lines_per_cell + 7 * deg)
                          * this->degree) &&
                    (shape_index
                       < (GeometryInfo<3>::lines_per_cell + 9 * deg)
@@ -2102,19 +2102,13 @@ FE_Nedelec<dim>::has_support_on_face (const unsigned int shape_index,
                       >= (GeometryInfo<3>::lines_per_cell + 2 * deg)
                          * this->degree) &&
                    (shape_index
-                      < (GeometryInfo<3>::lines_per_cell + 4 * deg)
+                      < (GeometryInfo<3>::lines_per_cell + 5 * deg)
                         * this->degree)) ||
                    ((shape_index
-                       >= (GeometryInfo<3>::lines_per_cell + 5 * deg)
-                          * this->degree) &&
-                    (shape_index
-                       < (GeometryInfo<3>::lines_per_cell + 6 * deg)
-                         * this->degree)) ||
-                   ((shape_index
-                       >= (GeometryInfo<3>::lines_per_cell + 7 * deg)
+                       >= (GeometryInfo<3>::lines_per_cell + 6 * deg)
                           * this->degree) &&
                     (shape_index
-                       < (GeometryInfo<3>::lines_per_cell + 8 * deg)
+                       < (GeometryInfo<3>::lines_per_cell + 7 * deg)
                          * this->degree)) ||
                    ((shape_index
                        >= (GeometryInfo<3>::lines_per_cell + 9 * deg)
@@ -2219,19 +2213,13 @@ FE_Nedelec<dim>::has_support_on_face (const unsigned int shape_index,
                       < (GeometryInfo<3>::lines_per_cell + 3 * deg)
                         * this->degree)) ||
                   ((shape_index
-                      >= (GeometryInfo<3>::lines_per_cell + 4 * deg)
+                      >= (GeometryInfo<3>::lines_per_cell + 5 * deg)
                          * this->degree) &&
                    (shape_index
-                      < (GeometryInfo<3>::lines_per_cell + 5 * deg)
-                        * this->degree)) ||
-                  ((shape_index
-                      >= (GeometryInfo<3>::lines_per_cell + 6 * deg)
-                         * this->degree) &&
-                   (shape_index
-                      < (GeometryInfo<3>::lines_per_cell + 7 * deg)
+                      < (GeometryInfo<3>::lines_per_cell + 6 * deg)
                         * this->degree)) ||
                   ((shape_index
-                      >= (GeometryInfo<3>::lines_per_cell + 8 * deg)
+                      >= (GeometryInfo<3>::lines_per_cell + 7 * deg)
                          * this->degree) &&
                    (shape_index
                       < (GeometryInfo<3>::lines_per_cell + 10 * deg)
@@ -2253,16 +2241,16 @@ FE_Nedelec<dim>::has_support_on_face (const unsigned int shape_index,
                       < (GeometryInfo<3>::lines_per_cell + 3 * deg)
                         * this->degree)) ||
                   ((shape_index
-                      >= (GeometryInfo<3>::lines_per_cell + 4 * deg)
+                      >= (GeometryInfo<3>::lines_per_cell + 5 * deg)
                          * this->degree) &&
                    (shape_index
-                      < (GeometryInfo<3>::lines_per_cell + 5 * deg)
+                      < (GeometryInfo<3>::lines_per_cell + 6 * deg)
                         * this->degree)) ||
                   ((shape_index
-                      >= (GeometryInfo<3>::lines_per_cell + 6 * deg)
+                      >= (GeometryInfo<3>::lines_per_cell + 7 * deg)
                          * this->degree) &&
                    (shape_index
-                      < (GeometryInfo<3>::lines_per_cell + 7 * deg)
+                      < (GeometryInfo<3>::lines_per_cell + 8 * deg)
                         * this->degree)) ||
                   ((shape_index
                       >= (GeometryInfo<3>::lines_per_cell + 10 * deg)