]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Periodic bc: Relocate first_vector_components parameter
authorMatthias Maier <tamiko@43-1.org>
Fri, 21 Aug 2015 21:24:21 +0000 (16:24 -0500)
committerMatthias Maier <tamiko@43-1.org>
Sun, 23 Aug 2015 01:52:00 +0000 (20:52 -0500)
This commit removes the first_vector_components parameter from
GridTools::collect_periodic_faces(), as well as, the
first_vector_components data field from the struct
GridTools::PeriodicFacePair. The parameter is no added to all varianst of
DoFTools::make_periodicity_constraints instead.

This is an incompatible change made for consistency: A PeriodicFacePair
should not store first_vector_components and with that information about
the underlying finite element space.

include/deal.II/dofs/dof_tools.h
include/deal.II/grid/grid_tools.h
source/dofs/dof_tools_constraints.cc
source/dofs/dof_tools_constraints.inst.in
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in

index b627eeed71f1fce6d6c88a284a06f1febc9c6248..b1ae1229679c65e11d06e40cc96564e03228975a 100644 (file)
@@ -896,22 +896,20 @@ namespace DoFTools
    * and any combination of that...
    * @endcode
    *
-   * Optionally a matrix @p matrix along with an std::vector @p
-   * first_vector_components can be specified that describes how DoFs on @p
-   * face_1 should be modified prior to constraining to the DoFs of @p face_2.
-   * Here, two declarations are possible: If the std::vector @p
-   * first_vector_components is non empty the matrix is interpreted as a @p
-   * dim $\times$ @p dim rotation matrix that is applied to all vector valued
-   * blocks listed in @p first_vector_components of the FESystem. If @p
-   * first_vector_components is empty the matrix is interpreted as an
+   * Optionally a matrix @p matrix along with an std::vector
+   * @p first_vector_components can be specified that describes how DoFs on
+   * @p face_1 should be modified prior to constraining to the DoFs of
+   * @p face_2. Here, two declarations are possible: If the std::vector
+   * @p first_vector_components is non empty the matrix is interpreted as a
+   * @p dim $\times$ @p dim rotation matrix that is applied to all vector
+   * valued blocks listed in @p first_vector_components of the FESystem. If
+   * @p first_vector_components is empty the matrix is interpreted as an
    * interpolation matrix with size no_face_dofs $\times$ no_face_dofs.
    *
    * Detailed information can be found in the see
    * @ref GlossPeriodicConstraints "Glossary entry on periodic boundary conditions".
    *
-   * @todo: Reference to soon be written example step and glossary article.
-   *
-   * @author Matthias Maier, 2012, 2014
+   * @author Matthias Maier, 2012 - 2015
    */
   template<typename FaceIterator>
   void
@@ -947,15 +945,16 @@ namespace DoFTools
    * @ref GlossPeriodicConstraints "Glossary entry on periodic boundary conditions"
    * for further information.
    *
-   * @author Daniel Arndt, Matthias Maier, 2013, 2014
+   * @author Daniel Arndt, Matthias Maier, 2013 - 2015
    */
   template<typename DH>
   void
   make_periodicity_constraints
   (const std::vector<GridTools::PeriodicFacePair<typename DH::cell_iterator> >
-   &periodic_faces,
-   dealii::ConstraintMatrix &constraint_matrix,
-   const ComponentMask      &component_mask = ComponentMask());
+                                   &periodic_faces,
+   dealii::ConstraintMatrix        &constraint_matrix,
+   const ComponentMask             &component_mask = ComponentMask(),
+   const std::vector<unsigned int> &first_vector_components = std::vector<unsigned int>());
 
 
 
index ca335e247386d032a2a12b220836a09585451fbe..3239513e252c7971702879601fcce69b3899fffd 100644 (file)
@@ -936,24 +936,18 @@ namespace GridTools
     std::bitset<3> orientation;
 
     /**
-     * A matrix that describes how vector valued DoFs of the first face should
-     * be modified prior to constraining to the DoFs of the second face. If
-     * the std::vector first_vector_components is non empty the matrix is
-     * interpreted as a @p dim $\times$ @p dim rotation matrix that is applied
-     * to all vector valued blocks listed in @p first_vector_components of the
-     * FESystem. If @p first_vector_components is empty the matrix is
-     * interpreted as an interpolation matrix with size no_face_dofs $\times$
-     * no_face_dofs. For more details see make_periodicity_constraints() and
-     * the glossary
-     * @ref GlossPeriodicConstraints "glossary entry on periodic conditions".
+     * A matrix that describes how vector valued DoFs of the first face
+     * should be modified prior to constraining to the DoFs of the second
+     * face. If the std::vector first_vector_components (supplied as a
+     * parameter to DofTools::make_periodicity_constraints()) is non empty
+     * the matrix is interpreted as a @p dim $\times$ @p dim rotation
+     * matrix that is applied to all vector valued blocks listed in @p
+     * first_vector_components of the FESystem. Alternatively, if @p
+     * first_vector_components is empty the matrix is interpreted as an
+     * interpolation matrix with size no_face_dofs $\times$ no_face_dofs.
+     * For more details see DoFTools::make_periodicity_constraints().
      */
     FullMatrix<double> matrix;
-
-    /**
-     * A vector of unsigned ints pointing to the first components of all
-     * vector valued blocks that should be rotated by matrix.
-     */
-    std::vector<unsigned int> first_vector_components;
   };
 
 
@@ -1047,9 +1041,10 @@ namespace GridTools
 
 
   /**
-   * This function will collect periodic face pairs on the coarsest mesh level
-   * of the given @p container (a Triangulation or DoFHandler) and add them to
-   * the vector @p matched_pairs leaving the original contents intact.
+   * This function will collect periodic face pairs on the coarsest mesh
+   * level of the given @p container (a Triangulation or DoFHandler) and
+   * add them to the vector @p matched_pairs leaving the original contents
+   * intact.
    *
    * Define a 'first' boundary as all boundary faces having boundary_id @p
    * b_id1 and a 'second' boundary consisting of all faces belonging to @p
@@ -1060,25 +1055,28 @@ namespace GridTools
    * orthogonal_equality().
    *
    * The bitset that is returned inside of PeriodicFacePair encodes the
-   * _relative_ orientation of the first face with respect to the second face,
-   * see the documentation of orthogonal_equality() for further details.
+   * _relative_ orientation of the first face with respect to the second
+   * face, see the documentation of orthogonal_equality() for further
+   * details.
    *
    * The @p direction refers to the space direction in which periodicity is
    * enforced.
    *
    * The @p offset is a vector tangential to the faces that is added to the
    * location of vertices of the 'first' boundary when attempting to match
-   * them to the corresponding vertices of the 'second' boundary. This can be
-   * used to implement conditions such as $u(0,y)=u(1,y+1)$.
-   *
-   * Optionally a (dim x dim) rotation matrix @p matrix along with a vector @p
-   * first_vector_components can be specified that describes how vector valued
-   * DoFs of the first face should be modified prior to constraining to the
-   * DoFs of the second face. If @p first_vector_components is non empty the
-   * matrix is interpreted as a rotation matrix that is applied to all vector
-   * valued blocks listed in @p first_vector_components of the FESystem. For
-   * more details see make_periodicity_constraints() and the glossary
-   * @ref GlossPeriodicConstraints "glossary entry on periodic conditions".
+   * them to the corresponding vertices of the 'second' boundary. This can
+   * be used to implement conditions such as $u(0,y)=u(1,y+1)$.
+   *
+   * Optionally a @p matrix can be specified that describes how vector
+   * valued DoFs of the first face should be modified prior to constraining
+   * to the DoFs of the second face. If the matrix has size
+   * $n\_face\_dofs\times n\_face\_dofs$, the periodicity constraints are
+   * applied as $dofs\_2 = matrix\cdot dofs\_1$. If the matrix has size
+   * $dim\times dim$, the matrix is interpreted as a rotation matrix for
+   * vector valued components.
+   * For more details see DoFTools::make_periodicity_constraints(), the
+   * glossary @ref GlossPeriodicConstraints "glossary entry on periodic
+   * conditions" and @ref step_45 "step-45".
    *
    * @tparam Container A type that satisfies the requirements of a mesh
    * container (see
@@ -1090,11 +1088,11 @@ namespace GridTools
    * periodicity algebraically.
    *
    * @note Because elements will be added to @p matched_pairs (and existing
-   * entries will be preserved), it is possible to call this function several
-   * times with different boundary ids to generate a vector with all periodic
-   * pairs.
+   * entries will be preserved), it is possible to call this function
+   * several times with different boundary ids to generate a vector with
+   * all periodic pairs.
    *
-   * @author Daniel Arndt, Matthias Maier, 2013, 2014
+   * @author Daniel Arndt, Matthias Maier, 2013 - 2015
    */
   template <typename CONTAINER>
   void
@@ -1105,8 +1103,7 @@ namespace GridTools
    const int                                                          direction,
    std::vector<PeriodicFacePair<typename CONTAINER::cell_iterator> > &matched_pairs,
    const Tensor<1,CONTAINER::space_dimension>                        &offset = dealii::Tensor<1,CONTAINER::space_dimension>(),
-   const FullMatrix<double>                                          &matrix = FullMatrix<double>(),
-   const std::vector<unsigned int>                                   &first_vector_components = std::vector<unsigned int>());
+   const FullMatrix<double>                                          &matrix = FullMatrix<double>());
 
 
   /**
@@ -1150,8 +1147,7 @@ namespace GridTools
    const int                                                          direction,
    std::vector<PeriodicFacePair<typename CONTAINER::cell_iterator> > &matched_pairs,
    const dealii::Tensor<1,CONTAINER::space_dimension>                &offset = dealii::Tensor<1,CONTAINER::space_dimension>(),
-   const FullMatrix<double>                                          &matrix = FullMatrix<double>(),
-   const std::vector<unsigned int>                                   &first_vector_components = std::vector<unsigned int>());
+   const FullMatrix<double>                                          &matrix = FullMatrix<double>());
 
   /*@}*/
   /**
index 3384d27aa9b7c61747f633d631634517570ce1e1..195caade26de4a527f6288e36b92403b9528793c 100644 (file)
@@ -2267,9 +2267,10 @@ namespace DoFTools
   void
   make_periodicity_constraints
   (const std::vector<GridTools::PeriodicFacePair<typename DH::cell_iterator> >
-   &periodic_faces,
-   dealii::ConstraintMatrix &constraint_matrix,
-   const ComponentMask      &component_mask)
+                                   &periodic_faces,
+   dealii::ConstraintMatrix        &constraint_matrix,
+   const ComponentMask             &component_mask,
+   const std::vector<unsigned int> &first_vector_components)
   {
     typedef std::vector<GridTools::PeriodicFacePair<typename DH::cell_iterator> >
     FaceVector;
@@ -2300,7 +2301,7 @@ namespace DoFTools
                                      it->orientation[1],
                                      it->orientation[2],
                                      it->matrix,
-                                     it->first_vector_components);
+                                     first_vector_components);
       }
   }
 
@@ -2310,12 +2311,12 @@ namespace DoFTools
 
   template<typename DH>
   void
-  make_periodicity_constraints (const DH                       &dof_handler,
-                                const types::boundary_id       b_id1,
-                                const types::boundary_id       b_id2,
-                                const int                      direction,
-                                dealii::ConstraintMatrix       &constraint_matrix,
-                                const ComponentMask            &component_mask)
+  make_periodicity_constraints (const DH                        &dof_handler,
+                                const types::boundary_id        b_id1,
+                                const types::boundary_id        b_id2,
+                                const int                       direction,
+                                dealii::ConstraintMatrix        &constraint_matrix,
+                                const ComponentMask             &component_mask)
   {
     static const int space_dim = DH::space_dimension;
     (void)space_dim;
@@ -2341,11 +2342,11 @@ namespace DoFTools
 
   template<typename DH>
   void
-  make_periodicity_constraints (const DH                       &dof_handler,
-                                const types::boundary_id       b_id,
-                                const int                      direction,
-                                dealii::ConstraintMatrix       &constraint_matrix,
-                                const ComponentMask            &component_mask)
+  make_periodicity_constraints (const DH                        &dof_handler,
+                                const types::boundary_id         b_id,
+                                const int                        direction,
+                                dealii::ConstraintMatrix        &constraint_matrix,
+                                const ComponentMask             &component_mask)
   {
     static const int dim = DH::dimension;
     static const int space_dim = DH::space_dimension;
index 420155d5d211e27123e2f099d316495c520b8db4..b652ebd906fe1c33da8deb430dd79eb2b46ec6db 100644 (file)
@@ -39,7 +39,8 @@ for (DH : DOFHANDLERS; deal_II_dimension : DIMENSIONS)
   DoFTools::make_periodicity_constraints<DH>
   (const std::vector<GridTools::PeriodicFacePair<DH::cell_iterator> > &,
    dealii::ConstraintMatrix &,
-   const ComponentMask &);
+   const ComponentMask &,
+   const std::vector<unsigned int> &);
 
 
   template
index 9d289634d0c6567228f395100f8ea65cd57f68ba..7a04c594b9015358092f3d7436d6929bce01f299 100644 (file)
@@ -2702,8 +2702,7 @@ next_cell:
    const int                                        direction,
    std::vector<PeriodicFacePair<CellIterator> >     &matched_pairs,
    const dealii::Tensor<1,CellIterator::AccessorType::space_dimension> &offset,
-   const FullMatrix<double>                         &matrix,
-   const std::vector<unsigned int>                  &first_vector_components)
+   const FullMatrix<double>                         &matrix)
   {
     static const int space_dim = CellIterator::AccessorType::space_dimension;
     (void)space_dim;
@@ -2741,8 +2740,7 @@ next_cell:
                   {cell1, cell2},
                   {face_idx1, face_idx2},
                   orientation,
-                  matrix,
-                  first_vector_components
+                  matrix
                 };
                 matched_pairs.push_back(matched_face);
                 pairs2.erase(it2);
@@ -2768,8 +2766,7 @@ next_cell:
    const int                                  direction,
    std::vector<PeriodicFacePair<typename CONTAINER::cell_iterator> > &matched_pairs,
    const Tensor<1,CONTAINER::space_dimension> &offset,
-   const FullMatrix<double>                   &matrix,
-   const std::vector<unsigned int>            &first_vector_components)
+   const FullMatrix<double>                   &matrix)
   {
     static const int dim = CONTAINER::dimension;
     static const int space_dim = CONTAINER::space_dimension;
@@ -2810,8 +2807,8 @@ next_cell:
             ExcMessage ("Unmatched faces on periodic boundaries"));
 
     // and call match_periodic_face_pairs that does the actual matching:
-    match_periodic_face_pairs(pairs1, pairs2, direction, matched_pairs,
-                              offset, matrix, first_vector_components);
+    match_periodic_face_pairs(pairs1, pairs2, direction, matched_pairs, offset,
+                              matrix);
   }
 
 
@@ -2824,8 +2821,7 @@ next_cell:
    const int                                  direction,
    std::vector<PeriodicFacePair<typename CONTAINER::cell_iterator> > &matched_pairs,
    const Tensor<1,CONTAINER::space_dimension> &offset,
-   const FullMatrix<double>                   &matrix,
-   const std::vector<unsigned int>            &first_vector_components)
+   const FullMatrix<double>                   &matrix)
   {
     static const int dim = CONTAINER::dimension;
     static const int space_dim = CONTAINER::space_dimension;
@@ -2873,8 +2869,8 @@ next_cell:
 #endif
 
     // and call match_periodic_face_pairs that does the actual matching:
-    match_periodic_face_pairs(pairs1, pairs2, direction, matched_pairs,
-                              offset, matrix, first_vector_components);
+    match_periodic_face_pairs(pairs1, pairs2, direction, matched_pairs, offset,
+                              matrix);
 
 #ifdef DEBUG
     //check for standard orientation
index 8acb27b95c1b26e9f381d61c17eadbf9cc50a1ee..c1c5fe585c539ecada7c3360565abfe87e702532 100644 (file)
@@ -24,7 +24,7 @@ for (X : TRIANGULATION_AND_DOFHANDLERS; deal_II_dimension : DIMENSIONS ; deal_II
   template
     unsigned int
     find_closest_vertex (const X &,
-                        const Point<deal_II_space_dimension> &);
+                         const Point<deal_II_space_dimension> &);
 
   template
     std::vector<dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type>
@@ -37,19 +37,19 @@ for (X : TRIANGULATION_AND_DOFHANDLERS; deal_II_dimension : DIMENSIONS ; deal_II
   template
     std::pair<dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type, Point<deal_II_dimension> >
     find_active_cell_around_point (const Mapping<deal_II_dimension, deal_II_space_dimension> &,
-                                  const X &,
-                                  const Point<deal_II_space_dimension> &);
+                                   const X &,
+                                   const Point<deal_II_space_dimension> &);
 
   template
     std::list<std::pair<X::cell_iterator, X::cell_iterator> >
     get_finest_common_cells (const X &mesh_1,
-                            const X &mesh_2);
+                               const X &mesh_2);
 
 
   template
     bool
     have_same_coarse_mesh (const X &mesh_1,
-                          const X &mesh_2);
+                           const X &mesh_2);
 
   \}
 
@@ -67,33 +67,33 @@ for (deal_II_dimension : DIMENSIONS ; deal_II_space_dimension : SPACE_DIMENSIONS
   template
     unsigned int
     find_closest_vertex (const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &,
-                        const Point<deal_II_space_dimension> &);
+                         const Point<deal_II_space_dimension> &);
 
   template
     std::vector<dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension> >::type>
     find_cells_adjacent_to_vertex(const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &,
-                                 const unsigned int);
+                                  const unsigned int);
   template
     dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension> >::type
     find_active_cell_around_point (const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &,
-                                  const Point<deal_II_space_dimension> &p);
+                                   const Point<deal_II_space_dimension> &p);
 
   template
     std::pair<dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension> >::type, Point<deal_II_dimension> >
     find_active_cell_around_point (const Mapping<deal_II_dimension, deal_II_space_dimension> &,
-                                  const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &,
-                                  const Point<deal_II_space_dimension> &);
+                                   const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &,
+                                   const Point<deal_II_space_dimension> &);
 
   template
     std::list<std::pair<parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension>::cell_iterator, parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension>::cell_iterator> >
     get_finest_common_cells (const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &mesh_1,
-                            const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &mesh_2);
+                               const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &mesh_2);
 
 
   template
     bool
     have_same_coarse_mesh (const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &mesh_1,
-                          const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &mesh_2);
+                           const parallel::distributed::Triangulation<deal_II_dimension,deal_II_space_dimension> &mesh_2);
 
   \}
 
@@ -107,13 +107,13 @@ for (deal_II_space_dimension : SPACE_DIMENSIONS)
 
     dealii::internal::ActiveCellIterator<deal_II_space_dimension, deal_II_space_dimension, parallel::distributed::Triangulation<deal_II_space_dimension, deal_II_space_dimension> >::type
     find_active_cell_around_point (const parallel::distributed::Triangulation<deal_II_space_dimension> &,
-                                  const Point<deal_II_space_dimension> &p);
+                                   const Point<deal_II_space_dimension> &p);
 
 
     std::pair<dealii::internal::ActiveCellIterator<deal_II_space_dimension, deal_II_space_dimension, parallel::distributed::Triangulation<deal_II_space_dimension, deal_II_space_dimension> >::type, Point<deal_II_space_dimension> >
     find_active_cell_around_point (const Mapping<deal_II_space_dimension> &,
-                                  const parallel::distributed::Triangulation<deal_II_space_dimension> &,
-                                  const Point<deal_II_space_dimension> &);
+                                   const parallel::distributed::Triangulation<deal_II_space_dimension> &,
+                                   const Point<deal_II_space_dimension> &);
 }
 
 
@@ -135,28 +135,28 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
 
     template
       void delete_unused_vertices (std::vector<Point<deal_II_space_dimension> > &,
-                                  std::vector<CellData<deal_II_dimension> > &,
-                                  SubCellData &);
+                                   std::vector<CellData<deal_II_dimension> > &,
+                                   SubCellData &);
 
     template
       void delete_duplicated_vertices (std::vector<Point<deal_II_space_dimension> > &,
-                                      std::vector<CellData<deal_II_dimension> > &,
-                                      SubCellData &,
-                                      std::vector<unsigned int> &,
-                                      double);
+                                       std::vector<CellData<deal_II_dimension> > &,
+                                       SubCellData &,
+                                       std::vector<unsigned int> &,
+                                       double);
 
     template
       void shift<deal_II_dimension> (const Tensor<1,deal_II_space_dimension> &,
-                                               Triangulation<deal_II_dimension, deal_II_space_dimension> &);
+                                                Triangulation<deal_II_dimension, deal_II_space_dimension> &);
 
     template
       void scale<deal_II_dimension> (const double,
-                                    Triangulation<deal_II_dimension, deal_II_space_dimension> &);
+                                     Triangulation<deal_II_dimension, deal_II_space_dimension> &);
 
     template
       void distort_random<deal_II_dimension> (const double,
-                                    Triangulation<deal_II_dimension, deal_II_space_dimension> &,
-                                    const bool);
+                                     Triangulation<deal_II_dimension, deal_II_space_dimension> &,
+                                     const bool);
 
     template
       void get_face_connectivity_of_cells
@@ -175,16 +175,16 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
 
     template
       void partition_triangulation (const unsigned int,
-                              Triangulation<deal_II_dimension, deal_II_space_dimension> &);
+                               Triangulation<deal_II_dimension, deal_II_space_dimension> &);
 
     template
       void partition_triangulation (const unsigned int,
-                                   const SparsityPattern &,
-                                   Triangulation<deal_II_dimension, deal_II_space_dimension> &);
+                                    const SparsityPattern &,
+                                    Triangulation<deal_II_dimension, deal_II_space_dimension> &);
 
     template
       std::pair<hp::DoFHandler<deal_II_dimension, deal_II_space_dimension>::active_cell_iterator,
-               Point<deal_II_dimension> >
+                Point<deal_II_dimension> >
       find_active_cell_around_point
       (const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension> &,
        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
@@ -194,7 +194,7 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
 
     template
       void get_subdomain_association (const Triangulation<deal_II_dimension, deal_II_space_dimension>  &,
-                                     std::vector<types::subdomain_id> &);
+                                      std::vector<types::subdomain_id> &);
 
     template
     unsigned int count_cells_with_subdomain_association(
@@ -303,8 +303,7 @@ for (X : TRIANGULATION_AND_DOFHANDLERS; deal_II_dimension : DIMENSIONS ; deal_II
                                       const int,
                                       std::vector<PeriodicFacePair<X::cell_iterator> > &,
                                       const Tensor<1,X::space_dimension> &,
-                                      const FullMatrix<double> &,
-                                      const std::vector<unsigned int> &);
+                                      const FullMatrix<double> &);
 
       template
       void collect_periodic_faces<X> (const X &,
@@ -312,8 +311,7 @@ for (X : TRIANGULATION_AND_DOFHANDLERS; deal_II_dimension : DIMENSIONS ; deal_II
                                       const int,
                                       std::vector<PeriodicFacePair<X::cell_iterator> > &,
                                       const Tensor<1,X::space_dimension> &,
-                                      const FullMatrix<double> &,
-                                      const std::vector<unsigned int> &);
+                                      const FullMatrix<double> &);
 
     #endif
 
@@ -336,8 +334,7 @@ for (deal_II_dimension : DIMENSIONS ; deal_II_space_dimension : SPACE_DIMENSIONS
                                   const int,
                                   std::vector<PeriodicFacePair<parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension>::cell_iterator> > &,
                                   const Tensor<1,parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension>::space_dimension> &,
-                                  const FullMatrix<double> &,
-                                  const std::vector<unsigned int> &);
+                                  const FullMatrix<double> &);
 
       template
       void
@@ -347,8 +344,7 @@ for (deal_II_dimension : DIMENSIONS ; deal_II_space_dimension : SPACE_DIMENSIONS
                                   const int,
                                   std::vector<PeriodicFacePair<parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension>::cell_iterator> > &,
                                   const Tensor<1,parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension>::space_dimension> &,
-                                  const FullMatrix<double> &,
-                                  const std::vector<unsigned int> &);
+                                  const FullMatrix<double> &);
      \}
    #endif
 #endif

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.