From 5e346af00bee50bcdcf2911de8fea225bc3f14b6 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Mon, 6 Apr 1998 17:31:31 +0000 Subject: [PATCH] Add new function to distribute cell to dof data. git-svn-id: https://svn.dealii.org/trunk@147 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/dofs/dof_handler.h | 27 ++++++++++++- deal.II/deal.II/source/dofs/dof_handler.cc | 47 ++++++++++++++++++++++ 2 files changed, 73 insertions(+), 1 deletion(-) diff --git a/deal.II/deal.II/include/dofs/dof_handler.h b/deal.II/deal.II/include/dofs/dof_handler.h index 1d39fb5de8..efc298267d 100644 --- a/deal.II/deal.II/include/dofs/dof_handler.h +++ b/deal.II/deal.II/include/dofs/dof_handler.h @@ -27,6 +27,7 @@ template class TriaActiveIterator; template class Triangulation; +class dVector; class dSMatrix; class dSMatrixStruct; class ConstraintMatrix; @@ -601,7 +602,24 @@ class DoFHandler : public DoFDimensionInfo { * the transfer process. */ unsigned int max_transfer_entries (const unsigned int max_level_diff) const; - + + /** + * Take a vector of values which live on + * cells (e.g. an error per cell) and + * distribute it to the dofs in such a + * way that a finite element field results, + * which can then be further processed, + * e.g. for output. + * + * It is assumed that the number of + * elements in #cell_data# equals the + * number of active cells. The size of + * #dof_data# is adjusted to the right + * size. + */ + void distribute_cell_to_dof_vector (const vector &cell_data, + dVector &dof_data) const; + /** * @name Cell iterator functions */ @@ -1034,6 +1052,13 @@ class DoFHandler : public DoFDimensionInfo { * Exception */ DeclException0 (ExcFunctionNotUseful); + /** + * Exception + */ + DeclException2 (ExcWrongSize, + int, int, + << "The dimension " << arg1 << " of the vector is wrong. " + << "It should be " << arg2); protected: /** diff --git a/deal.II/deal.II/source/dofs/dof_handler.cc b/deal.II/deal.II/source/dofs/dof_handler.cc index 9bcedc4f06..fa51a436f3 100644 --- a/deal.II/deal.II/source/dofs/dof_handler.cc +++ b/deal.II/deal.II/source/dofs/dof_handler.cc @@ -1465,6 +1465,53 @@ unsigned int DoFHandler<2>::max_transfer_entries (const unsigned int max_level_d +template +void DoFHandler::distribute_cell_to_dof_vector (const vector &cell_data, + dVector &dof_data) const { + Assert (cell_data.size()==tria->n_active_cells(), + ExcWrongSize (cell_data.size(), tria->n_active_cells())); + + // assign the right size to the output vector + dof_data.reinit (n_dofs()); + // count how often we have added a value + // in the sum for each dof + vector touch_count (n_dofs(), 0); + + active_cell_iterator cell = begin_active(), + endc = end(); + unsigned int present_cell = 0; + vector dof_indices; + const unsigned int total_dofs = selected_fe->total_dofs; + + for (; cell!=endc; ++cell, ++present_cell) + { + dof_indices.resize (0); + cell->get_dof_indices (dof_indices); + for (unsigned int i=0; i::reserve_space (const FiniteElement<1> &fe) { Assert (tria->n_levels() > 0, ExcInvalidTriangulation()); // delete all levels and set them up -- 2.39.5