<h3>General</h3>
<ol>
+<li> New: The DoFTools::count_dofs_per_block function now also works
+for objects of type hp::DoFHandler.
+<br>
+(Jason Sheldon, 2011/11/13)
+
<li> Fixed: The Intel compiler complains that it can't copy Trilinos vector
reference objects, preventing the compiling of step-32. This is now fixed.
<br>
const unsigned int level,
const std::vector<bool>& selected_dofs = std::vector<bool>(),
unsigned int offset = 0);
-
+
/**
* Create a sparsity pattern which
* in each row lists the degrees
const bool boundary_patches = false,
const bool level_boundary_patches = false,
const bool single_cell_patches = false);
-
+
/**
* Create a sparsity pattern which
* in each row lists the degrees of
* conditions.
*
* Since the patches are defined
- * through refinement, th
+ * through refinement, th
*
* @arg <tt>block_list</tt>: the
* SparsityPattern into which the
const DH& dof_handler,
const unsigned int level,
const bool interior_dofs_only = false);
-
+
//@}
/**
* Extract a vector that represents the
* in the target_blocks argument
* if given.
*/
- template <int dim, int spacedim>
+ template <class DH>
void
- count_dofs_per_block (const DoFHandler<dim,spacedim>& dof_handler,
- std::vector<unsigned int>& dofs_per_block,
- std::vector<unsigned int> target_blocks
- = std::vector<unsigned int>());
+ count_dofs_per_block (const DH &dof,
+ std::vector<unsigned int> &dofs_per_block,
+ std::vector<unsigned int> target_block
+ = std::vector<unsigned int>());
/**
* @deprecated See the previous
- template <int dim, int spacedim>
+ template <class DH>
void
- count_dofs_per_block (const DoFHandler<dim,spacedim>& dof_handler,
- std::vector<unsigned int> &dofs_per_block,
- std::vector<unsigned int> target_block)
+ count_dofs_per_block (const DH &dof_handler,
+ std::vector<unsigned int> &dofs_per_block,
+ std::vector<unsigned int> target_block)
{
- const FiniteElement<dim,spacedim>& fe = dof_handler.get_fe();
-
- std::fill (dofs_per_block.begin(), dofs_per_block.end(), 0U);
+ const dealii::hp::FECollection<DH::dimension,DH::space_dimension>
+ fe_collection (dof_handler.get_fe());
+ Assert (fe_collection.size() < 256, ExcNotImplemented());
- // If the empty vector was given as
- // default argument, set up this
- // vector as identity.
- if (target_block.size()==0)
+ for (unsigned int this_fe=0; this_fe<fe_collection.size(); ++this_fe)
{
- target_block.resize(fe.n_blocks());
- for (unsigned int i=0; i<fe.n_blocks(); ++i)
- target_block[i] = i;
- }
- else
- Assert (target_block.size()==fe.n_blocks(),
- ExcDimensionMismatch(target_block.size(),
- fe.n_blocks()));
+ const FiniteElement<DH::dimension,DH::space_dimension>& fe = fe_collection[this_fe];
+ std::fill (dofs_per_block.begin(), dofs_per_block.end(), 0U);
+ // If the empty vector was given as
+ // default argument, set up this
+ // vector as identity.
+ if (target_block.size()==0)
+ {
+ target_block.resize(fe.n_blocks());
+ for (unsigned int i=0; i<fe.n_blocks(); ++i)
+ target_block[i] = i;
+ }
+ else
+ Assert (target_block.size()==fe.n_blocks(),
+ ExcDimensionMismatch(target_block.size(),
+ fe.n_blocks()));
- const unsigned int max_block
- = *std::max_element (target_block.begin(),
- target_block.end());
- const unsigned int n_target_blocks = max_block + 1;
- const unsigned int n_blocks = fe.n_blocks();
- AssertDimension (dofs_per_block.size(), n_target_blocks);
+ const unsigned int max_block
+ = *std::max_element (target_block.begin(),
+ target_block.end());
+ const unsigned int n_target_blocks = max_block + 1;
+ const unsigned int n_blocks = fe.n_blocks();
- // special case for only one
- // block. treat this first
- // since it does not require any
- // computations
- if (n_blocks == 1)
- {
- dofs_per_block[0] = dof_handler.n_dofs();
- return;
- }
- // otherwise determine the number
- // of dofs in each block
- // separately.
- std::vector<unsigned char> dofs_by_block (dof_handler.n_locally_owned_dofs());
- internal::extract_dofs_by_component (dof_handler, std::vector<bool>(),
- true, dofs_by_block);
+ AssertDimension (dofs_per_block.size(), n_target_blocks);
- // next count what we got
- for (unsigned int block=0; block<fe.n_blocks(); ++block)
- dofs_per_block[target_block[block]]
- += std::count(dofs_by_block.begin(), dofs_by_block.end(),
- block);
+ // special case for only one
+ // block. treat this first
+ // since it does not require any
+ // computations
+ if (n_blocks == 1)
+ {
+ dofs_per_block[0] = dof_handler.n_dofs();
+ return;
+ }
+ // otherwise determine the number
+ // of dofs in each block
+ // separately.
+ std::vector<unsigned char> dofs_by_block (dof_handler.n_locally_owned_dofs());
+ internal::extract_dofs_by_component (dof_handler, std::vector<bool>(),
+ true, dofs_by_block);
+
+ // next count what we got
+ for (unsigned int block=0; block<fe.n_blocks(); ++block)
+ dofs_per_block[target_block[block]]
+ += std::count(dofs_by_block.begin(), dofs_by_block.end(),
+ block);
#ifdef DEAL_II_USE_P4EST
#if DEAL_II_COMPILER_SUPPORTS_MPI
- // if we are working on a parallel
- // mesh, we now need to collect
- // this information from all
- // processors
- if (const parallel::distributed::Triangulation<dim> * tria
- = (dynamic_cast<const parallel::distributed::Triangulation<dim>*>
- (&dof_handler.get_tria())))
- {
- std::vector<unsigned int> local_dof_count = dofs_per_block;
- MPI_Allreduce ( &local_dof_count[0], &dofs_per_block[0], n_target_blocks,
- MPI_UNSIGNED, MPI_SUM, tria->get_communicator());
- }
+ // if we are working on a parallel
+ // mesh, we now need to collect
+ // this information from all
+ // processors
+ if (const parallel::distributed::Triangulation<DH::dimension,DH::space_dimension> * tria
+ = (dynamic_cast<const parallel::distributed::Triangulation<DH::dimension,DH::space_dimension>*>
+ (&dof_handler.get_tria())))
+ {
+ std::vector<unsigned int> local_dof_count = dofs_per_block;
+ MPI_Allreduce ( &local_dof_count[0], &dofs_per_block[0], n_target_blocks,
+ MPI_UNSIGNED, MPI_SUM, tria->get_communicator());
+ }
#endif
#endif
+ }
}
typename DH::cell_iterator cell;
typename DH::cell_iterator endc = dof_handler.end(level);
std::vector<unsigned int> indices;
-
+
unsigned int i=0;
for (cell=dof_handler.begin(level); cell != endc; ++i, ++cell)
{
}
}
}
-
-
+
+
template <class DH>
void make_single_patch(
SparsityPattern& block_list,
block_list.reinit(1, dof_handler.n_dofs(level), dof_handler.n_dofs(level));
typename DH::cell_iterator cell;
typename DH::cell_iterator endc = dof_handler.end(level);
-
+
std::vector<unsigned int> indices;
std::vector<bool> exclude;
-
+
for (cell=dof_handler.begin(level); cell != endc;++cell)
{
indices.resize(cell->get_fe().dofs_per_cell);
cell->get_mg_dof_indices(indices);
-
+
if (interior_only)
{
// Exclude degrees of
exclude.resize(fe.dofs_per_cell);
std::fill(exclude.begin(), exclude.end(), false);
const unsigned int dpf = fe.dofs_per_face;
-
+
for (unsigned int face=0;face<GeometryInfo<DH::dimension>::faces_per_cell;++face)
if (cell->at_boundary(face) || cell->neighbor(face)->level() != cell->level())
for (unsigned int i=0;i<dpf;++i)
}
}
-
+
template <class DH>
void make_child_patches(
SparsityPattern& block_list,
typename DH::cell_iterator cell;
typename DH::cell_iterator pcell = dof_handler.begin(level-1);
typename DH::cell_iterator endc = dof_handler.end(level-1);
-
+
std::vector<unsigned int> indices;
std::vector<bool> exclude;
-
+
for (unsigned int block = 0;pcell != endc;++pcell)
{
if (!pcell->has_children()) continue;
// which are on faces
// of the parent
const unsigned int dpf = fe.dofs_per_face;
-
+
for (unsigned int d=0;d<DH::dimension;++d)
{
const unsigned int face = GeometryInfo<DH::dimension>::vertex_to_face[child][d];
for (unsigned int i=0;i<dpf;++i)
exclude[fe.face_to_cell_index(i,face)] = false;
}
-
+
for (unsigned int i=0;i<n_dofs;++i)
if (!exclude[i])
block_list.add(block, indices[i]);
++block;
}
}
-
-
-
+
+
+
template <class DH>
void make_vertex_patches(
SparsityPattern& block_list,
// Estimate for the number of
// dofs at this point
std::vector<unsigned int> vertex_dof_count(dof_handler.get_tria().n_vertices(), 0);
-
+
// Identify all vertices active
// on this level and remember
// some data about them
// From now on, only vertices
// with positive dof count are
// "in".
-
+
// Remove vertices at boundaries
// or in corners
for (unsigned int vg=0;vg<vertex_dof_count.size();++vg)
for (unsigned int vg=0;vg<vertex_mapping.size();++vg)
if (vertex_dof_count[vg] != 0)
vertex_dof_count[vertex_mapping[vg]] = vertex_dof_count[vg];
-
+
// Now that we have all the data,
// we reduce it to the part we
// actually want
vertex_dof_count.resize(i);
-
+
// At this point, the list of
// patches is ready. Now we enter
// the dofs into the sparsity
for (unsigned int i=0;i<vertex_dof_count.size();++i)
deallog << ' ' << vertex_dof_count[i];
deallog << std::endl;
-
+
std::vector<unsigned int> indices;
std::vector<bool> exclude;
-
+
for (cell=dof_handler.begin(level); cell != endc; ++cell)
{
const FiniteElement<DH::dimension>& fe = cell->get_fe();
indices.resize(fe.dofs_per_cell);
cell->get_mg_dof_indices(indices);
-
+
for (unsigned int v=0;v<GeometryInfo<DH::dimension>::vertices_per_cell;++v)
{
const unsigned int vg = cell->vertex_index(v);
const unsigned int block = vertex_mapping[vg];
if (block == numbers::invalid_unsigned_int)
continue;
-
+
if (interior_only)
{
// Exclude degrees of
exclude.resize(fe.dofs_per_cell);
std::fill(exclude.begin(), exclude.end(), false);
const unsigned int dpf = fe.dofs_per_face;
-
+
for (unsigned int d=0;d<DH::dimension;++d)
{
const unsigned int a_face = GeometryInfo<DH::dimension>::vertex_to_face[v][d];
template
void
-DoFTools::count_dofs_per_block<deal_II_dimension> (
+DoFTools::count_dofs_per_block<DoFHandler<deal_II_dimension> > (
const DoFHandler<deal_II_dimension>&,
std::vector<unsigned int>&, std::vector<unsigned int>);
+template
+void
+DoFTools::count_dofs_per_block<hp::DoFHandler<deal_II_dimension> > (
+ const hp::DoFHandler<deal_II_dimension>&,
+ std::vector<unsigned int>&, std::vector<unsigned int>);
+
+template
+void
+DoFTools::count_dofs_per_block<MGDoFHandler<deal_II_dimension> > (
+ const MGDoFHandler<deal_II_dimension>&,
+ std::vector<unsigned int>&, std::vector<unsigned int>);
template
void
template
void
-DoFTools::count_dofs_per_block<deal_II_dimension,deal_II_dimension+1> (
+DoFTools::count_dofs_per_block<DoFHandler<deal_II_dimension,deal_II_dimension+1> > (
const DoFHandler<deal_II_dimension,deal_II_dimension+1>&,
std::vector<unsigned int>&, std::vector<unsigned int>);
+template
+void
+DoFTools::count_dofs_per_block<hp::DoFHandler<deal_II_dimension,deal_II_dimension+1> > (
+ const hp::DoFHandler<deal_II_dimension,deal_II_dimension+1>&,
+ std::vector<unsigned int>&, std::vector<unsigned int>);
+
+template
+void
+DoFTools::count_dofs_per_block<MGDoFHandler<deal_II_dimension,deal_II_dimension+1> > (
+ const MGDoFHandler<deal_II_dimension,deal_II_dimension+1>&,
+ std::vector<unsigned int>&, std::vector<unsigned int>);
+
#endif
template
--- /dev/null
+//---------------------------- count_dofs_per_block_01.cc ---------------------------
+// $Id$
+// Version: $Name$
+//
+// Copyright (C) 2009, 2011 by the deal.II authors
+//
+// This file is subject to QPL and may not be distributed
+// without copyright and license information. Please refer
+// to the file deal.II/doc/license.html for the text and
+// further information on this license.
+//
+//---------------------------- count_dofs_per_block_01.cc ---------------------------
+
+// like bits/count_dofs_per_block_01 but for hp::DoFHandler
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/hp/dof_handler.h>
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_dgp.h>
+#include <deal.II/fe/fe_nedelec.h>
+#include <deal.II/fe/fe_raviart_thomas.h>
+#include <deal.II/fe/fe_system.h>
+
+#include <fstream>
+#include <iomanip>
+#include <iomanip>
+#include <string>
+
+
+std::string output_file_name = "count_dofs_per_block_01/output";
+
+
+
+void print (const std::vector<unsigned int> &v)
+{
+ deallog << v.size();
+ for (unsigned int i=0; i<v.size(); ++i)
+ deallog << ' ' << v[i];
+ deallog << std::endl;
+}
+
+
+
+template <int dim>
+void
+check ()
+{
+ // create tria and dofhandler
+ // objects. set different boundary
+ // and sub-domain ids
+ Triangulation<dim> tria;
+ GridGenerator::hyper_cube(tria, 0., 1.);
+ tria.refine_global (1);
+ for (int i=0; i<2; ++i)
+ {
+ tria.begin_active()->set_refine_flag();
+ tria.execute_coarsening_and_refinement ();
+ }
+
+ // produce two copies of the Taylor-Hood
+ // element
+ hp::FECollection<dim> fe;
+ fe.push_back (FESystem<dim> (FE_Q<dim>(1), dim, FE_DGQ<dim>(0), 1));
+ fe.push_back (FESystem<dim> (FE_Q<dim>(1), dim, FE_DGQ<dim>(0), 1));
+
+ hp::DoFHandler<dim> dof_handler (tria);
+
+ typename hp::DoFHandler<dim>::active_cell_iterator
+ cell = dof_handler.begin_active(),
+ endc = dof_handler.end();
+ for (; cell!=endc; ++cell)
+ cell->set_active_fe_index (std::rand() % fe.size());
+
+ dof_handler.distribute_dofs (fe);
+
+ // no grouping
+ {
+ std::vector<unsigned int> dpc(dim+1);
+ DoFTools::count_dofs_per_component (dof_handler, dpc);
+ print (dpc);
+ }
+
+ {
+ std::vector<unsigned int> dpc(dim+1);
+ DoFTools::count_dofs_per_block (dof_handler, dpc);
+ print (dpc);
+ }
+
+
+ // grouping into less groups than
+ // components
+ {
+ std::vector<unsigned int> group(dim+1, 0);
+ group[dim] = 1;
+ std::vector<unsigned int> dpc(2);
+ DoFTools::count_dofs_per_component (dof_handler, dpc, false, group);
+ Assert (dpc.size() == 2, ExcInternalError());
+ print (dpc);
+ }
+
+ {
+ std::vector<unsigned int> group(dim+1, 0);
+ group[dim] = 1;
+ std::vector<unsigned int> dpc(2);
+ DoFTools::count_dofs_per_block (dof_handler, dpc, group);
+ Assert (dpc.size() == 2, ExcInternalError());
+ print (dpc);
+ }
+
+ // grouping into more groups than
+ // components
+ {
+ std::vector<unsigned int> group(dim+1, 2*dim);
+ group[dim] = 0;
+ std::vector<unsigned int> dpc(2*dim+1);
+ DoFTools::count_dofs_per_component (dof_handler, dpc, false, group);
+ Assert (dpc.size() == 2*dim+1, ExcInternalError());
+ print (dpc);
+ }
+
+ {
+ std::vector<unsigned int> group(dim+1, 2*dim);
+ group[dim] = 0;
+ std::vector<unsigned int> dpc(2*dim+1);
+ DoFTools::count_dofs_per_block (dof_handler, dpc, group);
+ Assert (dpc.size() == 2*dim+1, ExcInternalError());
+ print (dpc);
+ }
+}
+
+
+
+int main()
+{
+ std::ofstream logfile(output_file_name.c_str());
+ logfile << std::setprecision (2);
+ deallog << std::setprecision (2);
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+ deallog.threshold_double(1.e-10);
+
+ check<1> ();
+ check<2> ();
+ check<3> ();
+}
--- /dev/null
+
+DEAL::2 5 4
+DEAL::2 5 4
+DEAL::2 5 4
+DEAL::2 5 4
+DEAL::3 4 0 5
+DEAL::3 4 0 5
+DEAL::3 18 18 10
+DEAL::3 18 18 10
+DEAL::2 36 10
+DEAL::2 36 10
+DEAL::5 10 0 0 0 36
+DEAL::5 10 0 0 0 36
+DEAL::4 60 60 60 22
+DEAL::4 60 60 60 22
+DEAL::2 180 22
+DEAL::2 180 22
+DEAL::7 22 0 0 0 0 0 180
+DEAL::7 22 0 0 0 0 0 180
--- /dev/null
+//---------------------------- count_dofs_per_block_02.cc ---------------------------
+// $Id$
+// Version: $Name$
+//
+// Copyright (C) 2009, 2011 by the deal.II authors
+//
+// This file is subject to QPL and may not be distributed
+// without copyright and license information. Please refer
+// to the file deal.II/doc/license.html for the text and
+// further information on this license.
+//
+//---------------------------- count_dofs_per_block_02.cc ---------------------------
+
+// like _01 but with different elements
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/hp/dof_handler.h>
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_dgp.h>
+#include <deal.II/fe/fe_nedelec.h>
+#include <deal.II/fe/fe_raviart_thomas.h>
+#include <deal.II/fe/fe_system.h>
+
+#include <fstream>
+#include <iomanip>
+#include <iomanip>
+#include <string>
+
+
+std::string output_file_name = "count_dofs_per_block_02/output";
+
+
+
+void print (const std::vector<unsigned int> &v)
+{
+ deallog << v.size();
+ for (unsigned int i=0; i<v.size(); ++i)
+ deallog << ' ' << v[i];
+ deallog << std::endl;
+}
+
+
+
+template <int dim>
+void
+check ()
+{
+ // create tria and dofhandler
+ // objects. set different boundary
+ // and sub-domain ids
+ Triangulation<dim> tria;
+ GridGenerator::hyper_cube(tria, 0., 1.);
+ tria.refine_global (1);
+ for (int i=0; i<2; ++i)
+ {
+ tria.begin_active()->set_refine_flag();
+ tria.execute_coarsening_and_refinement ();
+ }
+
+ // produce two copies of the Taylor-Hood
+ // element
+ hp::FECollection<dim> fe;
+ fe.push_back (FESystem<dim> (FE_Q<dim>(1), dim, FE_DGQ<dim>(0), 1));
+ fe.push_back (FESystem<dim> (FE_Q<dim>(2), dim, FE_DGQ<dim>(1), 1));
+
+ hp::DoFHandler<dim> dof_handler (tria);
+
+ typename hp::DoFHandler<dim>::active_cell_iterator
+ cell = dof_handler.begin_active(),
+ endc = dof_handler.end();
+ for (; cell!=endc; ++cell)
+ cell->set_active_fe_index (std::rand() % fe.size());
+
+ dof_handler.distribute_dofs (fe);
+
+ // no grouping
+ {
+ std::vector<unsigned int> dpc(dim+1);
+ DoFTools::count_dofs_per_component (dof_handler, dpc);
+ print (dpc);
+ }
+
+ {
+ std::vector<unsigned int> dpc(dim+1);
+ DoFTools::count_dofs_per_block (dof_handler, dpc);
+ print (dpc);
+ }
+
+
+ // grouping into less groups than
+ // components
+ {
+ std::vector<unsigned int> group(dim+1, 0);
+ group[dim] = 1;
+ std::vector<unsigned int> dpc(2);
+ DoFTools::count_dofs_per_component (dof_handler, dpc, false, group);
+ Assert (dpc.size() == 2, ExcInternalError());
+ print (dpc);
+ }
+
+ {
+ std::vector<unsigned int> group(dim+1, 0);
+ group[dim] = 1;
+ std::vector<unsigned int> dpc(2);
+ DoFTools::count_dofs_per_block (dof_handler, dpc, group);
+ Assert (dpc.size() == 2, ExcInternalError());
+ print (dpc);
+ }
+
+ // grouping into more groups than
+ // components
+ {
+ std::vector<unsigned int> group(dim+1, 2*dim);
+ group[dim] = 0;
+ std::vector<unsigned int> dpc(2*dim+1);
+ DoFTools::count_dofs_per_component (dof_handler, dpc, false, group);
+ Assert (dpc.size() == 2*dim+1, ExcInternalError());
+ print (dpc);
+ }
+
+ {
+ std::vector<unsigned int> group(dim+1, 2*dim);
+ group[dim] = 0;
+ std::vector<unsigned int> dpc(2*dim+1);
+ DoFTools::count_dofs_per_block (dof_handler, dpc, group);
+ Assert (dpc.size() == 2*dim+1, ExcInternalError());
+ print (dpc);
+ }
+}
+
+
+
+int main()
+{
+ std::ofstream logfile(output_file_name.c_str());
+ logfile << std::setprecision (2);
+ deallog << std::setprecision (2);
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+ deallog.threshold_double(1.e-10);
+
+ check<1> ();
+ check<2> ();
+ check<3> ();
+}
--- /dev/null
+
+DEAL::2 8 7
+DEAL::2 8 7
+DEAL::2 8 7
+DEAL::2 8 7
+DEAL::3 7 0 8
+DEAL::3 7 0 8
+DEAL::3 45 45 28
+DEAL::3 45 45 28
+DEAL::2 90 28
+DEAL::2 90 28
+DEAL::5 28 0 0 0 90
+DEAL::5 28 0 0 0 90
+DEAL::4 205 205 205 92
+DEAL::4 205 205 205 92
+DEAL::2 615 92
+DEAL::2 615 92
+DEAL::7 92 0 0 0 0 0 615
+DEAL::7 92 0 0 0 0 0 615