]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add support for nonprimitive elements
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 26 Jan 2006 11:16:29 +0000 (11:16 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 26 Jan 2006 11:16:29 +0000 (11:16 +0000)
git-svn-id: https://svn.dealii.org/trunk@12179 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/source/dofs/dof_tools.cc

index 1fffaf1d5e36c4f61b2eaf32b0a9c9eedfbbfbfa..4639fec4707dfc6c5cf5ba5159f9950621acf70b 100644 (file)
@@ -1478,6 +1478,45 @@ class DoFTools
     map_support_points_to_dofs (const Mapping<dim>    &mapping,
                                const DH<dim>         &dof_handler,
                                std::map<Point<dim>, unsigned int, Comp> &point_to_index_map);
+
+                                    /**
+                                     * Map a coupling table from the
+                                     * user friendly organization by
+                                     * components to the organization
+                                     * by blocks. Specializations of
+                                     * this function for DoFHandler
+                                     * and hp::DoFHandler are
+                                     * required due to the different
+                                     * results of their finite
+                                     * element access.
+                                     *
+                                     * The return vector will be
+                                     * initialized to the correct
+                                     * length inside this function.
+                                     */
+    template <int dim>
+    void convert_couplings_to_blocks (const hp::DoFHandler<dim>& dof_handler,
+                                     const Table<2, Coupling>& table_by_component,
+                                     std::vector<Table<2,Coupling> >& tables_by_block);
+                                    /**
+                                     * Map a coupling table from the
+                                     * user friendly organization by
+                                     * components to the organization
+                                     * by blocks. Specializations of
+                                     * this function for DoFHandler
+                                     * and hp::DoFHandler are
+                                     * required due to the different
+                                     * results of their finite
+                                     * element access.
+                                     *
+                                     * The return vector will be
+                                     * initialized to the correct
+                                     * length inside this function.
+                                     */
+    template <int dim>
+    void convert_couplings_to_blocks (const DoFHandler<dim>& dof_handler,
+                                     const Table<2, Coupling>& table_by_component,
+                                     std::vector<Table<2,Coupling> >& tables_by_block);
     
                                     /**
                                      * Exception
@@ -1594,6 +1633,36 @@ class DoFTools
 };
 
 
+/**
+ * Operator computing the maximum coupling out of two.
+ *
+ * @relates DoFTools
+ */
+DoFTools::Coupling operator |= (const DoFTools::Coupling& c1,
+                               const DoFTools::Coupling& c2)
+{
+  if (c1 == DoFTools::always || c2 == DoFTools::always)
+    return DoFTools::always;
+  if (c1 == DoFTools::nonzero || c2 == DoFTools::nonzero)
+    return DoFTools::nonzero;
+  return DoFTools::none;
+}
+
+
+/**
+ * Operator computing the maximum coupling out of two.
+ *
+ * @relates DoFTools
+ */
+DoFTools::Coupling operator | (const DoFTools::Coupling& c1,
+                              const DoFTools::Coupling& c2)
+{
+  if (c1 == DoFTools::always || c2 == DoFTools::always)
+    return DoFTools::always;
+  if (c1 == DoFTools::nonzero || c2 == DoFTools::nonzero)
+    return DoFTools::nonzero;
+  return DoFTools::none;
+}
 
 
 // ---------------------- inline and template functions --------------------
index d880f8ee5b3f13640c13cb87064dbb69e8637514..ce2fe38089a7252263436e075db3afa3c4164306 100644 (file)
@@ -3550,6 +3550,62 @@ DoFTools::compute_dof_couplings(
 }
 
 
+template<int dim>
+void
+DoFTools::convert_couplings_to_blocks (
+  const DoFHandler<dim>& dof_handler,
+  const Table<2, Coupling>& table,
+  std::vector<Table<2,Coupling> >& tables_by_block)
+{
+  const FiniteElement<dim>& fe = dof_handler.get_fe();
+  const unsigned int nb = fe.n_blocks();
+  
+  tables_by_block.resize(1);
+  tables_by_block[0].reinit(nb, nb);
+  tables_by_block[0].fill(none);
+  
+  for (unsigned int i=0;i<fe.n_components();++i)
+    {
+      const unsigned int ib = fe.component_to_block(i);
+      for (unsigned int j=0;j<fe.n_components();++j)
+       {
+         const unsigned int jb = fe.component_to_block(j);
+         tables_by_block[0](ib,jb) |= table(i,j);
+       }
+    }
+}
+
+
+template<int dim>
+void
+DoFTools::convert_couplings_to_blocks (
+  const hp::DoFHandler<dim>& dof_handler,
+  const Table<2, Coupling>& table,
+  std::vector<Table<2,Coupling> >& tables_by_block)
+{
+  const hp::FECollection<dim>& coll = dof_handler.get_fe();
+  tables_by_block.resize(coll.n_finite_elements());
+
+  for (unsigned int f=0;f<coll.n_finite_elements();++n)
+    {
+      const FiniteElement<dim>& fe = coll.ge_fe(f);
+      
+      const unsigned int nb = fe.n_blocks();  
+      tables_by_block[f].reinit(nb, nb);
+      tables_by_block[f].fill(none);
+      for (unsigned int i=0;i<fe.n_components();++i)
+       {
+         const unsigned int ib = fe.component_to_block(i);
+         for (unsigned int j=0;j<fe.n_components();++j)
+           {
+             const unsigned int jb = fe.component_to_block(j);
+             tables_by_block[f](ib,jb) |= table(i,j);
+           }
+       }
+    }
+}
+
+
 // explicit instantiations
 template void
 DoFTools::compute_row_length_vector(
@@ -4041,13 +4097,25 @@ DoFTools::map_dof_to_boundary_indices<deal_II_dimension>
 template
 void
 DoFTools::map_dofs_to_support_points<deal_II_dimension>
-(const Mapping<deal_II_dimension>       &mapping,
- const DoFHandler<deal_II_dimension>    &dof_handler,
- std::vector<Point<deal_II_dimension> > &support_points);
+(const Mapping<deal_II_dimension>&,
+ const DoFHandler<deal_II_dimension>&,
+ std::vector<Point<deal_II_dimension> >&);
 
 template
 void
 DoFTools::compute_dof_couplings<deal_II_dimension>(
-  Table<2,Coupling>&       dof_couplings,
-  const Table<2,Coupling>& component_couplings,
-  const FiniteElement<deal_II_dimension>&      fe);
+  Table<2,Coupling>&, const Table<2,Coupling>&,
+  const FiniteElement<deal_II_dimension>&);
+
+template
+void
+DoFTools::convert_couplings_to_blocks (
+  const DoFHandler<deal_II_dimension>&, const Table<2, Coupling>&,
+  std::vector<Table<2,Coupling> >&);
+
+template
+void
+DoFTools::convert_couplings_to_blocks (
+  const hp::DoFHandler<deal_II_dimension>&, const Table<2, Coupling>&,
+  std::vector<Table<2,Coupling> >&);
+

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.