]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Change the order of indices in the {sub,}face_interpolation matrices to make it confo...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 16 Sep 2006 12:14:20 +0000 (12:14 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 16 Sep 2006 12:14:20 +0000 (12:14 +0000)
git-svn-id: https://svn.dealii.org/trunk@13912 0785d39b-7218-0410-832d-ea1e28bc413d

12 files changed:
deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_dgp.h
deal.II/deal.II/include/fe/fe_dgp_nonparametric.h
deal.II/deal.II/include/fe/fe_dgq.h
deal.II/deal.II/include/fe/fe_q.h
deal.II/deal.II/include/fe/fe_system.h
deal.II/deal.II/source/dofs/dof_tools.cc
deal.II/deal.II/source/fe/fe_dgp.cc
deal.II/deal.II/source/fe/fe_dgp_monomial.cc
deal.II/deal.II/source/fe/fe_dgp_nonparametric.cc
deal.II/deal.II/source/fe/fe_q.cc
deal.II/deal.II/source/fe/fe_system.cc

index cdb37a4239f973319b176f138de75203701e9e67..aa3efe0c05e6fe5bd7a63b2b5861ab91390c8a4d 100644 (file)
@@ -952,8 +952,8 @@ class FiniteElement : public Subscriptor,
                                      * of one element to the face of
                                      * the neighboring element. 
                                      * The size of the matrix is
-                                     * then @p dofs_per_face times
-                                     * <tt>source.dofs_per_face</tt>.
+                                     * then <tt>source.dofs_per_face</tt> times
+                                     * <tt>this->dofs_per_face</tt>.
                                      *
                                      * Derived elements will have to
                                      * implement this function. They
@@ -978,8 +978,8 @@ class FiniteElement : public Subscriptor,
                                      * of one element to the subface of
                                      * the neighboring element. 
                                      * The size of the matrix is
-                                     * then @p dofs_per_face times
-                                     * <tt>source.dofs_per_face</tt>.
+                                     * then <tt>source.dofs_per_face</tt> times
+                                     * <tt>this->dofs_per_face</tt>.
                                      *
                                      * Derived elements will have to
                                      * implement this function. They
index 874010a29bf5670cd21a93156bad2323d30f9dde..a8ee8f2c83b5d6c04a40b9c370401a20a883d278 100644 (file)
@@ -218,8 +218,8 @@ class FE_DGP : public FE_Poly<PolynomialSpace<dim>,dim>
                                      * of one element to the face of
                                      * the neighboring element. 
                                      * The size of the matrix is
-                                     * then @p dofs_per_face times
-                                     * <tt>source.dofs_per_face</tt>.
+                                     * then <tt>source.dofs_per_face</tt> times
+                                     * <tt>this->dofs_per_face</tt>.
                                      *
                                      * Derived elements will have to
                                      * implement this function. They
@@ -243,8 +243,8 @@ class FE_DGP : public FE_Poly<PolynomialSpace<dim>,dim>
                                      * of one element to the face of
                                      * the neighboring element. 
                                      * The size of the matrix is
-                                     * then @p dofs_per_face times
-                                     * <tt>source.dofs_per_face</tt>.
+                                     * then <tt>source.dofs_per_face</tt> times
+                                     * <tt>this->dofs_per_face</tt>.
                                      *
                                      * Derived elements will have to
                                      * implement this function. They
index 1f043d1844a8fb964fd0741e24581e7f1bda7db9..b1510bf9ba214047f170abba146fb602fbe07f55 100644 (file)
@@ -211,8 +211,8 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                      * of one element to the face of
                                      * the neighboring element. 
                                      * The size of the matrix is
-                                     * then @p dofs_per_face times
-                                     * <tt>source.dofs_per_face</tt>.
+                                     * then <tt>source.dofs_per_face</tt> times
+                                     * <tt>this->dofs_per_face</tt>.
                                      *
                                      * Derived elements will have to
                                      * implement this function. They
@@ -236,8 +236,8 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                      * of one element to the face of
                                      * the neighboring element. 
                                      * The size of the matrix is
-                                     * then @p dofs_per_face times
-                                     * <tt>source.dofs_per_face</tt>.
+                                     * then <tt>source.dofs_per_face</tt> times
+                                     * <tt>this->dofs_per_face</tt>.
                                      *
                                      * Derived elements will have to
                                      * implement this function. They
index b73f97ba45624dea5c75097a393896061ad74c02..b4498b4c3386893156f80af456dfb865fa2e759c 100644 (file)
@@ -128,8 +128,8 @@ class FE_DGQ : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                      * of one element to the face of
                                      * the neighboring element. 
                                      * The size of the matrix is
-                                     * then @p dofs_per_face times
-                                     * <tt>source.dofs_per_face</tt>.
+                                     * then <tt>source.dofs_per_face</tt> times
+                                     * <tt>this->dofs_per_face</tt>.
                                      *
                                      * Derived elements will have to
                                      * implement this function. They
@@ -153,8 +153,8 @@ class FE_DGQ : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                      * of one element to the face of
                                      * the neighboring element. 
                                      * The size of the matrix is
-                                     * then @p dofs_per_face times
-                                     * <tt>source.dofs_per_face</tt>.
+                                     * then <tt>source.dofs_per_face</tt> times
+                                     * <tt>this->dofs_per_face</tt>.
                                      *
                                      * Derived elements will have to
                                      * implement this function. They
index 31733b0a537ba11a4401fa6b6564b6f0d171f265..257749bbe0ff9d4887793f36d93ea9064287cc3b 100644 (file)
@@ -280,8 +280,8 @@ class FE_Q : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                      * of one element to the face of
                                      * the neighboring element. 
                                      * The size of the matrix is
-                                     * then @p dofs_per_face times
-                                     * <tt>source.dofs_per_face</tt>.
+                                     * then <tt>source.dofs_per_face</tt> times
+                                     * <tt>this->dofs_per_face</tt>.
                                      *
                                      * Derived elements will have to
                                      * implement this function. They
@@ -305,8 +305,8 @@ class FE_Q : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                      * of one element to the face of
                                      * the neighboring element. 
                                      * The size of the matrix is
-                                     * then @p dofs_per_face times
-                                     * <tt>source.dofs_per_face</tt>.
+                                     * then <tt>source.dofs_per_face</tt> times
+                                     * <tt>this->dofs_per_face</tt>.
                                      *
                                      * Derived elements will have to
                                      * implement this function. They
index b9e764201b3cbbf9633767ab52e30c638ea4baf8..77673d4352e29a3ff07adda6cfba7c81ae2c9aa9 100644 (file)
@@ -415,8 +415,8 @@ class FESystem : public FiniteElement<dim>
                                      * of one element to the face of
                                      * the neighboring element. 
                                      * The size of the matrix is
-                                     * then @p dofs_per_face times
-                                     * <tt>source.dofs_per_face</tt>.
+                                     * then <tt>source.dofs_per_face</tt> times
+                                     * <tt>this->dofs_per_face</tt>.
                                      *
                                      * Base elements of this element will
                                      * have to implement this function. They
@@ -442,8 +442,8 @@ class FESystem : public FiniteElement<dim>
                                      * of one element to the subface of
                                      * the neighboring element. 
                                      * The size of the matrix is
-                                     * then @p dofs_per_face times
-                                     * <tt>source.dofs_per_face</tt>.
+                                     * then <tt>source.dofs_per_face</tt> times
+                                     * <tt>this->dofs_per_face</tt>.
                                      *
                                      * Base elements of this element will
                                      * have to implement this function. They
index 4b81357ae4c9027ec936abee18ec8acc65dc09c4..6953c58638eb492fde5bece02d9eeb6586685f62 100644 (file)
@@ -1575,10 +1575,10 @@ namespace internal
                                         // previous steps of the hp hanging node
                                         // procedure.
 
-       Assert (face_constraints.m () == dofs_mother.size (),
+       Assert (face_constraints.n () == dofs_mother.size (),
                ExcDimensionMismatch(dofs_mother.size (),
                                     face_constraints.n()));
-       Assert (face_constraints.n () == dofs_child.size (),
+       Assert (face_constraints.m () == dofs_child.size (),
                ExcDimensionMismatch(dofs_child.size (),
                                     face_constraints.m()));
 
@@ -1588,7 +1588,7 @@ namespace internal
        Assert (n_dofs_mother <= n_dofs_child,
                ExcInternalError ());
 
-       for (unsigned int col=0; col!=n_dofs_child; ++col
+       for (unsigned int row=0; row!=n_dofs_child; ++row
          {
            bool constraint_already_satisfied = false;
 
@@ -1598,8 +1598,8 @@ namespace internal
                                             // the corresponding global dof
                                             // indices
            for (unsigned int i=0; i<n_dofs_mother; ++i)
-             if (face_constraints (i,col) == 1.0)
-               if (dofs_mother[i] == dofs_child[col])
+             if (face_constraints (row,i) == 1.0)
+               if (dofs_mother[i] == dofs_child[row])
                  {
                    constraint_already_satisfied = true;
                    break;
@@ -1607,19 +1607,19 @@ namespace internal
 
            if (constraint_already_satisfied == false)
              { 
-               constraints.add_line (dofs_child[col]);
+               constraints.add_line (dofs_child[row]);
                for (unsigned int i=0; i!=n_dofs_mother; ++i)
-                  if (face_constraints(i,col) != 0)
+                  if (face_constraints(row,i) != 0)
                    {
 #ifdef WOLFGANG
-                     std::cout << "   " << dofs_child[col]
-                               << " -> " << face_constraints (i,col) << " * "
+                     std::cout << "   " << dofs_child[row]
+                               << " -> " << face_constraints (row,i) << " * "
                                << dofs_mother[i]
                                << std::endl;
 #endif
-                     constraints.add_entry (dofs_child[col],
+                     constraints.add_entry (dofs_child[row],
                                             dofs_mother[i],
-                                            face_constraints (i,col));
+                                            face_constraints (row,i));
                    }
              }  
          }      
@@ -2241,8 +2241,8 @@ namespace internal
                                                         // properties of a
                                                         // finite element onto
                                                         // that mesh
-                       face_constraints.reinit (n_dofs_on_mother,
-                                                n_dofs_on_children);
+                       face_constraints.reinit (n_dofs_on_children,
+                                                n_dofs_on_mother);
                        cell->get_fe()
                          .get_subface_interpolation_matrix (subface->get_fe(subface_fe_index),
                                                             c, face_constraints);
@@ -2473,13 +2473,13 @@ namespace internal
                    dofs_on_mother.resize (n_dofs_on_mother);
                    cell->face(face)->get_dof_indices (dofs_on_mother, cell->active_fe_index ());
 
-                   FullMatrix<double> fc_mother_ipol (n_dofs_on_children,
-                                                      n_dofs_on_mother);
-                   FullMatrix<double> fc_mother_sface (n_dofs_on_children,
-                                                       n_dofs_on_mother);
+                   FullMatrix<double> fc_mother_ipol (n_dofs_on_mother,
+                                                      n_dofs_on_children);
+                   FullMatrix<double> fc_mother_sface (n_dofs_on_mother,
+                                                       n_dofs_on_children);
                    dominating_fe.get_face_interpolation_matrix (cell->get_fe(),
                                                                 fc_mother_ipol);
-                   fc_ipol_sface.mmult (fc_mother_sface, fc_mother_ipol);
+                   fc_mother_ipol.mmult (fc_mother_sface, fc_ipol_sface);
 
                                                     // Add constraints to global constraint
                                                     // matrix.
@@ -2522,15 +2522,15 @@ namespace internal
                                              
                                                           // Now create the element
                                                           // constraint for this subface.
-                         FullMatrix<double> fc_child_sface_ipol (n_dofs_on_children,
-                                                                 n_dofs_on_mother);
-                         FullMatrix<double> fc_child_sface_sface (n_dofs_on_children,
-                                                                  n_dofs_on_mother);
+                         FullMatrix<double> fc_child_sface_ipol (n_dofs_on_mother,
+                                                                 n_dofs_on_children);
+                         FullMatrix<double> fc_child_sface_sface (n_dofs_on_mother,
+                                                                  n_dofs_on_children);
                          dominating_fe.get_subface_interpolation_matrix 
                            (subface->get_fe(subface_fe_index),
                             c, fc_child_sface_ipol);
 
-                         fc_ipol_sface.mmult (fc_child_sface_sface, fc_child_sface_ipol);
+                         fc_child_sface_ipol.mmult (fc_child_sface_sface, fc_ipol_sface);
 
                                                           // Add constraints to global constraint
                                                           // matrix.
@@ -2606,8 +2606,8 @@ namespace internal
                                                  
                                                        // Now create the element
                                                        // constraint for this subface.
-                      FullMatrix<double> face_constraints (dofs_on_mother.size(),
-                                                           dofs_on_children.size());
+                      FullMatrix<double> face_constraints (dofs_on_children.size(),
+                                                           dofs_on_mother.size());
                       cell->get_fe().get_face_interpolation_matrix (neighbor->get_fe (),
                                                                     face_constraints);
 
index 88e9491620ec825a619316c0989480e416457efa..9f21d1801b3e775badfdb38a916b6209359a463a 100644 (file)
@@ -111,7 +111,7 @@ get_face_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
          ExcDimensionMismatch (interpolation_matrix.m(),
                                0));
   Assert (interpolation_matrix.n() == 0,
-         ExcDimensionMismatch (interpolation_matrix.m(),
+         ExcDimensionMismatch (interpolation_matrix.n(),
                                0));
 }
 
@@ -140,7 +140,7 @@ get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
          ExcDimensionMismatch (interpolation_matrix.m(),
                                0));
   Assert (interpolation_matrix.n() == 0,
-         ExcDimensionMismatch (interpolation_matrix.m(),
+         ExcDimensionMismatch (interpolation_matrix.n(),
                                0));
 }
 
index 9495615b87e708b3e3da6d7a52ed121ca648bf23..fb42e547cb0c7579032d0cf1ea98809bf859bf8f 100644 (file)
@@ -272,7 +272,7 @@ get_face_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
          ExcDimensionMismatch (interpolation_matrix.m(),
                                0));
   Assert (interpolation_matrix.n() == 0,
-         ExcDimensionMismatch (interpolation_matrix.m(),
+         ExcDimensionMismatch (interpolation_matrix.n(),
                                0));
 }
 
@@ -301,7 +301,7 @@ get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
          ExcDimensionMismatch (interpolation_matrix.m(),
                                0));
   Assert (interpolation_matrix.n() == 0,
-         ExcDimensionMismatch (interpolation_matrix.m(),
+         ExcDimensionMismatch (interpolation_matrix.n(),
                                0));
 }
 
index 3934f5f75f65eabf8f6eaf3418e796b27450f83f..402e9ca2a7c69c9a6c10e380195abaa547457176 100644 (file)
@@ -441,7 +441,7 @@ get_face_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
          ExcDimensionMismatch (interpolation_matrix.m(),
                                0));
   Assert (interpolation_matrix.n() == 0,
-         ExcDimensionMismatch (interpolation_matrix.m(),
+         ExcDimensionMismatch (interpolation_matrix.n(),
                                0));
 }
 
@@ -470,7 +470,7 @@ get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
          ExcDimensionMismatch (interpolation_matrix.m(),
                                0));
   Assert (interpolation_matrix.n() == 0,
-         ExcDimensionMismatch (interpolation_matrix.m(),
+         ExcDimensionMismatch (interpolation_matrix.n(),
                                0));
 }
 
index 0f4d12c95e36cd2bebd47eaea0d37c188c7d1421..2338673e4f5fcdb0df753636df912a10832c5df3 100644 (file)
@@ -364,10 +364,10 @@ get_face_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
                typename FiniteElement<dim>::
                ExcInterpolationNotImplemented());
   
-  Assert (interpolation_matrix.m() == this->dofs_per_face,
-         ExcDimensionMismatch (interpolation_matrix.m(),
+  Assert (interpolation_matrix.n() == this->dofs_per_face,
+         ExcDimensionMismatch (interpolation_matrix.n(),
                                this->dofs_per_face));
-  Assert (interpolation_matrix.n() == x_source_fe.dofs_per_face,
+  Assert (interpolation_matrix.m() == x_source_fe.dofs_per_face,
          ExcDimensionMismatch (interpolation_matrix.m(),
                                x_source_fe.dofs_per_face));
 
@@ -431,11 +431,11 @@ get_face_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
          if (fabs (matrix_entry) < eps)
            matrix_entry = 0.0;
 
-         interpolation_matrix(j,i) = matrix_entry;
+         interpolation_matrix(i,j) = matrix_entry;
        }  
     }
 
-                                  // make sure that the column sum of
+                                  // make sure that the row sum of
                                   // each of the matrices is 1 at
                                   // this point. this must be so
                                   // since the shape functions sum up
@@ -445,7 +445,7 @@ get_face_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
       double sum = 0.;
 
       for (unsigned int i=0; i<this->dofs_per_face; ++i)
-        sum += interpolation_matrix(i,j);
+        sum += interpolation_matrix(j,i);
 
       Assert (std::fabs(sum-1) < 2e-14*this->degree*(dim-1),
               ExcInternalError());
@@ -469,10 +469,10 @@ get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
                typename FiniteElement<dim>::
                ExcInterpolationNotImplemented());
   
-  Assert (interpolation_matrix.m() == this->dofs_per_face,
-         ExcDimensionMismatch (interpolation_matrix.m(),
+  Assert (interpolation_matrix.n() == this->dofs_per_face,
+         ExcDimensionMismatch (interpolation_matrix.n(),
                                this->dofs_per_face));
-  Assert (interpolation_matrix.n() == x_source_fe.dofs_per_face,
+  Assert (interpolation_matrix.m() == x_source_fe.dofs_per_face,
          ExcDimensionMismatch (interpolation_matrix.m(),
                                x_source_fe.dofs_per_face));
 
@@ -533,11 +533,11 @@ get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
          if (fabs (matrix_entry) < eps)
            matrix_entry = 0.0;
 
-         interpolation_matrix(j,i) = matrix_entry;
+         interpolation_matrix(i,j) = matrix_entry;
        }
     }
 
-                                  // make sure that the column sum of
+                                  // make sure that the row sum of
                                   // each of the matrices is 1 at
                                   // this point. this must be so
                                   // since the shape functions sum up
@@ -547,7 +547,7 @@ get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
       double sum = 0.;
 
       for (unsigned int i=0; i<this->dofs_per_face; ++i)
-        sum += interpolation_matrix(i,j);
+        sum += interpolation_matrix(j,i);
 
       Assert (std::fabs(sum-1) < 2e-14*this->degree*(dim-1),
               ExcInternalError());
index 126de218c70ca0d9d09b79c88736ba2e44d4cb9b..0f17d4981e3987918c46eb3ce0476fe88a25cc7b 100644 (file)
@@ -1872,10 +1872,10 @@ get_face_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
               typename FiniteElement<dim>::
               ExcInterpolationNotImplemented());
 
-  Assert (interpolation_matrix.m() == this->dofs_per_face,
-         ExcDimensionMismatch (interpolation_matrix.m(),
+  Assert (interpolation_matrix.n() == this->dofs_per_face,
+         ExcDimensionMismatch (interpolation_matrix.n(),
                                this->dofs_per_face));
-  Assert (interpolation_matrix.n() == x_source_fe.dofs_per_face,
+  Assert (interpolation_matrix.m() == x_source_fe.dofs_per_face,
          ExcDimensionMismatch (interpolation_matrix.m(),
                                x_source_fe.dofs_per_face));
 
@@ -1922,8 +1922,8 @@ get_face_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
              ExcNotImplemented());
 
                                       // get the interpolation from the bases
-      base_to_base_interpolation.reinit (base.dofs_per_face,
-                                        base_other.dofs_per_face);
+      base_to_base_interpolation.reinit (base_other.dofs_per_face,
+                                        base.dofs_per_face);
       base.get_face_interpolation_matrix (base_other,
                                          base_to_base_interpolation);
       
@@ -1940,9 +1940,9 @@ get_face_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
            if (fe_other_system->face_system_to_base_index(j).first
                ==
                std::make_pair (base_index_other, multiplicity_other))
-             interpolation_matrix(i, j)
-               = base_to_base_interpolation(this->face_system_to_base_index(i).second,
-                                            fe_other_system->face_system_to_base_index(j).second);
+             interpolation_matrix(j, i)
+               = base_to_base_interpolation(fe_other_system->face_system_to_base_index(j).second,
+                                            this->face_system_to_base_index(i).second);
          
                                       // advance to the next base element
                                       // for this and the other
@@ -1999,10 +1999,10 @@ get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
               typename FiniteElement<dim>::
               ExcInterpolationNotImplemented());
 
-  Assert (interpolation_matrix.m() == this->dofs_per_face,
-         ExcDimensionMismatch (interpolation_matrix.m(),
+  Assert (interpolation_matrix.n() == this->dofs_per_face,
+         ExcDimensionMismatch (interpolation_matrix.n(),
                                this->dofs_per_face));
-  Assert (interpolation_matrix.n() == x_source_fe.dofs_per_face,
+  Assert (interpolation_matrix.m() == x_source_fe.dofs_per_face,
          ExcDimensionMismatch (interpolation_matrix.m(),
                                x_source_fe.dofs_per_face));
 
@@ -2049,8 +2049,8 @@ get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
              ExcNotImplemented());
 
                                       // get the interpolation from the bases
-      base_to_base_interpolation.reinit (base.dofs_per_face,
-                                        base_other.dofs_per_face);
+      base_to_base_interpolation.reinit (base_other.dofs_per_face,
+                                        base.dofs_per_face);
       base.get_subface_interpolation_matrix (base_other,
                                             subface,
                                             base_to_base_interpolation);
@@ -2068,9 +2068,9 @@ get_subface_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
            if (fe_other_system->face_system_to_base_index(j).first
                ==
                std::make_pair (base_index_other, multiplicity_other))
-             interpolation_matrix(i, j)
-               = base_to_base_interpolation(this->face_system_to_base_index(i).second,
-                                            fe_other_system->face_system_to_base_index(j).second);
+             interpolation_matrix(j, i)
+               = base_to_base_interpolation(fe_other_system->face_system_to_base_index(j).second,
+                                            this->face_system_to_base_index(i).second);
          
                                       // advance to the next base element
                                       // for this and the other

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.