]> https://gitweb.dealii.org/ - dealii.git/commitdiff
improved version of DoFTools::count_dofs_per_component for non primitives
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 17 Jun 2005 10:37:50 +0000 (10:37 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 17 Jun 2005 10:37:50 +0000 (10:37 +0000)
tests for Raviart-Thomas-Elements added

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

deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/source/dofs/dof_tools.cc
deal.II/deal.II/source/multigrid/mg_transfer_block.cc
tests/bits/dof_tools_00a.cc
tests/bits/dof_tools_common.cc
tests/multigrid/transfer.cc

index 4bb5181299b1bfe25215a9cae15ab87e835116a4..1d678a28492182a84cb2098e1180443aa562e701 100644 (file)
@@ -906,6 +906,17 @@ class DoFTools
                                      * than the total number of
                                      * degrees of freedom.
                                      *
+                                     * This behavior can be switched
+                                     * off by the optional parameter
+                                     * <tt>vector_valued_once</tt>. If
+                                     * this is <tt>true</tt>, the
+                                     * number of components of a
+                                     * nonprimitive vector valued
+                                     * element is collected only in
+                                     * the first component. All other
+                                     * components will have a count
+                                     * of zero.
+                                     *
                                      * The additional optional
                                      * argument @p target_component
                                      * allows for a re-sorting and
@@ -934,6 +945,7 @@ class DoFTools
     static void
     count_dofs_per_component (const DoFHandler<dim>&     dof_handler,
                              std::vector<unsigned int>& dofs_per_component,
+                             const bool vector_valued_once = false,
                              std::vector<unsigned int>  target_component
                              =std::vector<unsigned int>());
     
index 9728a3d8f5c02b7c0184b0b0bc547b2e505136ac..0ecf810e9417628be3165c62ae89e2929cec3a80 100644 (file)
@@ -1723,9 +1723,11 @@ void
 DoFTools::count_dofs_per_component (
   const DoFHandler<dim>&     dof_handler,
   std::vector<unsigned int>& dofs_per_component,
+  bool only_once,
   std::vector<unsigned int>  target_component)
 {
-  const unsigned int n_components = dof_handler.get_fe().n_components();
+  const FiniteElement<dim>& fe = dof_handler.get_fe();
+  const unsigned int n_components = fe.n_components();
   dofs_per_component.resize (n_components);
   std::fill (dofs_per_component.begin(), dofs_per_component.end(), 0U);
   
@@ -1750,8 +1752,9 @@ DoFTools::count_dofs_per_component (
     {
       dofs_per_component[0] = dof_handler.n_dofs();
       return;
-    };
-  
+    }
+
+      
                                   // otherwise determine the number
                                   // of dofs in each component
                                   // separately. do so in parallel
@@ -1775,12 +1778,27 @@ DoFTools::count_dofs_per_component (
   threads.join_all ();
 
                                   // next count what we got
-  for (unsigned int i=0; i<n_components; ++i)
-    dofs_per_component[target_component[i]]
-      += std::count(dofs_in_component[i].begin(),
-                   dofs_in_component[i].end(),
-                   true);
-
+  unsigned int component = 0;
+  for (unsigned int b=0;b<fe.n_base_elements();++b)
+    {
+      const FiniteElement<dim>& base = fe.base_element(b);
+                                      // Dimension of base element
+      unsigned int d = base.n_components();
+      
+      for (unsigned int m=0;m<fe.element_multiplicity(b);++m)
+       {
+         for (unsigned int dd=0;dd<d;++dd)
+           {
+             if (base.is_primitive() || (!only_once || dd==0))
+             dofs_per_component[target_component[component]]
+               += std::count(dofs_in_component[component].begin(),
+                             dofs_in_component[component].end(),
+                             true);
+             ++component;
+           }
+       }
+    }
+  
                                   // finally sanity check. this is
                                   // only valid if the finite element
                                   // is actually primitive, so
@@ -3118,7 +3136,7 @@ template
 void
 DoFTools::count_dofs_per_component<deal_II_dimension> (
   const DoFHandler<deal_II_dimension>&,
-  std::vector<unsigned int>&, std::vector<unsigned int>);
+  std::vector<unsigned int>&, bool, std::vector<unsigned int>);
 
 template
 void
index 5a74928df1d251bf1e802bd0b82b4b5ffdda1469..c7d9aac3f412080a3451733235b1f25b0f674cd6 100644 (file)
@@ -113,7 +113,7 @@ void MGTransferBlockBase::build_matrices (
   component_start.resize(target_component.size());
   DoFTools::
     count_dofs_per_component (static_cast<const DoFHandler<dim>&>(mg_dof),
-                              component_start, target_component);
+                              component_start, true, target_component);
 
   unsigned int k=0;
   for (unsigned int i=0;i<component_start.size();++i)
index 426c4ba6590bea260d041790c9a68ad37736a49a..f661e89ebcaba54be9ee8687d870c65aa40f0509 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$ 
 //
-//    Copyright (C) 2003, 2004 by the deal.II authors
+//    Copyright (C) 2003, 2004, 2005 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -38,12 +38,28 @@ check_this (const DoFHandler<dim> &dof_handler)
     deallog << ' ' << dofs_per_component[i];
   deallog << std::endl;
 
+  DoFTools::count_dofs_per_component (dof_handler,
+                                     dofs_per_component, true);
+
+  for (unsigned int i=0;i<n_components;++i)
+    deallog << ' ' << dofs_per_component[i];
+  deallog << std::endl;
+
   if (n_components>1)
     {
       std::vector<unsigned int> target_component(n_components,0U);
 
       DoFTools::count_dofs_per_component (dof_handler,
                                          dofs_per_component,
+                                         false,
+                                         target_component);
+      for (unsigned int i=0;i<n_components;++i)
+       deallog << ' ' << dofs_per_component[i];
+      deallog << std::endl;
+
+      DoFTools::count_dofs_per_component (dof_handler,
+                                         dofs_per_component,
+                                         true,
                                          target_component);
       for (unsigned int i=0;i<n_components;++i)
        deallog << ' ' << dofs_per_component[i];
@@ -54,6 +70,14 @@ check_this (const DoFHandler<dim> &dof_handler)
       
       DoFTools::count_dofs_per_component (dof_handler,
                                          dofs_per_component,
+                                         false,
+                                         target_component);
+      for (unsigned int i=0;i<n_components;++i)
+       deallog << ' ' << dofs_per_component[i];
+      deallog << std::endl;
+      DoFTools::count_dofs_per_component (dof_handler,
+                                         dofs_per_component,
+                                         true,
                                          target_component);
       for (unsigned int i=0;i<n_components;++i)
        deallog << ' ' << dofs_per_component[i];
index 93837fbb608f5f5f0faf777d94539bc90400c280..9729eb0cde0620b5eb4b2bc22c7fdf09062e2d3b 100644 (file)
@@ -26,6 +26,7 @@
 #include <fe/fe_dgq.h>
 #include <fe/fe_dgp.h>
 #include <fe/fe_nedelec.h>
+#include <fe/fe_raviart_thomas.h>
 #include <fe/fe_system.h>
 
 #include <fstream>
@@ -158,8 +159,18 @@ main()
 
       CHECK(Nedelec, 1, 2);
       CHECK(Nedelec, 1, 3);
-  
 
+      CHECK(RaviartThomas, 0, 2);
+      CHECK(RaviartThomas, 1, 2);
+      CHECK(RaviartThomas, 2, 2);
+
+      CHECK(RaviartThomasNodal, 0, 2);
+      CHECK(RaviartThomasNodal, 1, 2);
+      CHECK(RaviartThomasNodal, 2, 2);
+      CHECK(RaviartThomasNodal, 0, 3);
+      CHECK(RaviartThomasNodal, 1, 3);
+      CHECK(RaviartThomasNodal, 2, 3);
+      
       CHECK_SYS1(FE_Q<1>(1),  3,1);
       CHECK_SYS1(FE_DGQ<1>(2),2,1);
       CHECK_SYS1(FE_DGP<1>(3),1,1);
index f33bb3279084f852397670bae4a70cb06e938753..b82ca5370efdaf3ab95c7eddf346d62adc181587 100644 (file)
@@ -90,7 +90,7 @@ void check_select(const FiniteElement<dim>& fe,
   mgdof.distribute_dofs(fe);
   DoFRenumbering::component_wise(mgdof, target_component);
   vector<unsigned int> ndofs(fe.n_components());
-  DoFTools::count_dofs_per_component(mgdof, ndofs, target_component);
+  DoFTools::count_dofs_per_component(mgdof, ndofs, true, target_component);
   
   for (unsigned int l=0;l<tr.n_levels();++l)
     DoFRenumbering::component_wise(mgdof, l, mg_target_component);

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.