From: bangerth Date: Tue, 28 Aug 2007 21:36:37 +0000 (+0000) Subject: Fixed: the function X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=94033384f62cd49adbfb5d03a7082b3ee5bfebf8;p=dealii-svn.git Fixed: the function DataOut::build_patches had a quadratic algorithm when generatic cell-data (as opposed to DoF data). This algorithm became a bottleneck when generating output on meshes with large number of cells. This is now fixed. git-svn-id: https://svn.dealii.org/trunk@15079 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/source/numerics/data_out.cc b/deal.II/deal.II/source/numerics/data_out.cc index 181117637b..56ab3c83a5 100644 --- a/deal.II/deal.II/source/numerics/data_out.cc +++ b/deal.II/deal.II/source/numerics/data_out.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors +// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -409,19 +409,41 @@ void DataOut::build_some_patches (Data &data) const unsigned int n_q_points = patch_points.n_quadrature_points; - typename std::vector< dealii::DataOutBase::Patch >::iterator patch = this->patches.begin(); + typename std::vector< dealii::DataOutBase::Patch >::iterator + patch = this->patches.begin(); cell_iterator cell=first_cell(); + + // keep track of the index of the + // current cell so we can + // efficiently evaluate cell-based + // data (as opposed to DoF-based + // data). we do so only if + // this->cell_data.size() != 0 + unsigned int cell_index = (this->cell_data.size() != 0 + ? + std::distance (this->dofs->begin_active(), + active_cell_iterator (cell)) + : + deal_II_numbers::invalid_unsigned_int); + // get first cell in this thread for (unsigned int i=0; (idofs->end()); ++i) { ++patch; + + const cell_iterator old_cell = cell; + cell = next_cell(cell); + + if (this->cell_data.size() != 0) + cell_index += std::distance (active_cell_iterator(old_cell), + active_cell_iterator(cell)); } // now loop over all cells and // actually create the patches - for (;cell != this->dofs->end();) + for (; cell != this->dofs->end();) { Assert (patch != this->patches.end(), ExcInternalError()); @@ -483,8 +505,8 @@ void DataOut::build_some_patches (Data &data) for (unsigned int q=0; qdata(dataset*data.n_components+component,q) = data.patch_values_system[q](component); - }; - }; + } + } // then do the cell data. only // compute the number of a cell if @@ -495,15 +517,11 @@ void DataOut::build_some_patches (Data &data) if (this->cell_data.size() != 0) { Assert (!cell->has_children(), ExcNotImplemented()); -//TODO: This introduces a quadratic component into the algorithm. We may be doing better if we always kept track of the cell_number when 'cell' is increased, but we have to make sure that we only do this if cell_number is actually used, since we don't want to enforce the activeness constraints otherwise - const unsigned int cell_number - = std::distance (this->dofs->begin_active(), - active_cell_iterator (cell)); for (unsigned int dataset=0; datasetcell_data.size(); ++dataset) { const double value - = this->cell_data[dataset]->get_cell_data_value (cell_number); + = this->cell_data[dataset]->get_cell_data_value (cell_index); for (unsigned int q=0; qdata(dataset+this->dof_data.size()*data.n_components,q) = value; @@ -564,7 +582,7 @@ void DataOut::build_some_patches (Data &data) // for the neighbor index patch->neighbors[f] = this->patches[(*data.cell_to_patch_index_map) [neighbor->level()][neighbor->index()]].patch_index; - }; + } // next cell (patch) in this // thread @@ -572,7 +590,14 @@ void DataOut::build_some_patches (Data &data) (idofs->end()); ++i) { ++patch; + + const cell_iterator old_cell = cell; + cell = next_cell(cell); + + if (this->cell_data.size() != 0) + cell_index += std::distance (active_cell_iterator(old_cell), + active_cell_iterator(cell)); } } } diff --git a/deal.II/doc/news/changes.html b/deal.II/doc/news/changes.html index 4fbfee7886..7c15ac8058 100644 --- a/deal.II/doc/news/changes.html +++ b/deal.II/doc/news/changes.html @@ -1115,6 +1115,16 @@ inconvenience this causes.

deal.II

    +
  1. Fixed: the function + DataOut::build_patches + had a quadratic algorithm when generatic cell-data (as opposed + to DoF data). This algorithm became a bottleneck when + generating output on meshes with large number of cells. This is + now fixed. +
    + (WB 2007/08/28) +

    +
  2. New: the function DoFTools::get_active_fe_indices extracts for each cell the active finite element index used on it.