template <int dim, int spacedim, bool lda>
static
void
- update_cell_dof_indices_cache (const DoFCellAccessor<dealii::hp::DoFHandler<dim,spacedim>, lda> &)
+ update_cell_dof_indices_cache (const DoFCellAccessor<dealii::hp::DoFHandler<dim,spacedim>, lda> &accessor)
{
-//TODO[WB]: should implement a dof indices cache for hp as well
+ // caches are only for cells with DoFs, i.e., for active ones
+ if (accessor.has_children())
+ return;
- // not implemented, but should also
- // not be called
- Assert (false, ExcNotImplemented());
+ const unsigned int dofs_per_cell = accessor.get_fe().dofs_per_cell;
+
+ // make sure the cache is at least
+ // as big as we need it when
+ // writing to the last element of
+ // this cell
+ Assert (static_cast<unsigned int>(accessor.present_index)
+ <
+ accessor.dof_handler->levels[accessor.present_level]
+ ->cell_cache_offsets.size(),
+ ExcInternalError());
+ Assert (accessor.dof_handler->levels[accessor.present_level]
+ ->cell_cache_offsets[accessor.present_index]
+ <=
+ accessor.dof_handler->levels[accessor.present_level]
+ ->cell_dof_indices_cache.size(),
+ ExcInternalError());
+
+ std::vector<types::global_dof_index> dof_indices (dofs_per_cell);
+ accessor.dealii::DoFAccessor<dim,dealii::hp::DoFHandler<dim,spacedim>,lda>::get_dof_indices
+ (dof_indices, accessor.active_fe_index());
+
+ types::global_dof_index *next_dof_index
+ = &accessor.dof_handler->levels[accessor.present_level]
+ ->cell_dof_indices_cache[accessor.dof_handler->levels[accessor.present_level]
+ ->cell_cache_offsets[accessor.present_index]];
+ for (unsigned int i=0; i<dofs_per_cell; ++i, ++next_dof_index)
+ *next_dof_index = dof_indices[i];
}
+
+
/**
* Implement setting dof
* indices on a
get_dof_indices (const DoFCellAccessor<dealii::hp::DoFHandler<dim,spacedim>, lda> &accessor,
std::vector<types::global_dof_index> &dof_indices)
{
- // no caching for
- // hp::DoFHandler
- // implemented
- typedef
- dealii::DoFAccessor<dim,dealii::hp::DoFHandler<dim,spacedim>, lda>
- DoFAccessor;
- accessor.DoFAccessor::get_dof_indices (dof_indices,
- accessor.active_fe_index());
+ AssertDimension (dof_indices.size(), accessor.get_fe().dofs_per_cell);
+
+ Assert (!accessor.has_children(),
+ ExcMessage ("Cell must be active."));
+
+ types::global_dof_index *cache = &accessor.dof_handler->levels[accessor.level()]
+ ->cell_dof_indices_cache[accessor.dof_handler->levels[accessor.level()]->cell_cache_offsets[accessor.present_index]];
+ for (unsigned int i=0; i<accessor.get_fe().dofs_per_cell; ++i, ++cache)
+ dof_indices[i] = *cache;
}
/**
= std::vector<DoFLevel::offset_type> (
dof_handler.tria->n_raw_lines(level),
(DoFLevel::offset_type)(-1));
+ dof_handler.levels[level]->cell_cache_offsets
+ = std::vector<DoFLevel::offset_type> (
+ dof_handler.tria->n_raw_lines(level),
+ (DoFLevel::offset_type)(-1));
types::global_dof_index next_free_dof = 0;
+ types::global_dof_index cache_size = 0;
for (typename DoFHandler<dim,spacedim>::active_cell_iterator
cell=dof_handler.begin_active(level);
cell!=dof_handler.end_active(level); ++cell)
{
dof_handler.levels[level]->dof_offsets[cell->index()] = next_free_dof;
next_free_dof += cell->get_fe().dofs_per_line;
+
+ dof_handler.levels[level]->cell_cache_offsets[cell->index()] = cache_size;
+ cache_size += cell->get_fe().dofs_per_cell;
}
dof_handler.levels[level]->dof_indices
= std::vector<types::global_dof_index> (next_free_dof,
DoFHandler<dim,spacedim>::invalid_dof_index);
+ dof_handler.levels[level]->cell_dof_indices_cache
+ = std::vector<types::global_dof_index> (cache_size,
+ DoFHandler<dim,spacedim>::invalid_dof_index);
}
// safety check: make sure that
= std::vector<DoFLevel::offset_type> (
dof_handler.tria->n_raw_quads(level),
(DoFLevel::offset_type)(-1));
+ dof_handler.levels[level]->cell_cache_offsets
+ = std::vector<DoFLevel::offset_type> (
+ dof_handler.tria->n_raw_quads(level),
+ (DoFLevel::offset_type)(-1));
types::global_dof_index next_free_dof = 0;
+ types::global_dof_index cache_size = 0;
for (typename DoFHandler<dim,spacedim>::active_cell_iterator
cell=dof_handler.begin_active(level);
cell!=dof_handler.end_active(level); ++cell)
{
dof_handler.levels[level]->dof_offsets[cell->index()] = next_free_dof;
next_free_dof += cell->get_fe().dofs_per_quad;
+
+ dof_handler.levels[level]->cell_cache_offsets[cell->index()] = cache_size;
+ cache_size += cell->get_fe().dofs_per_cell;
}
dof_handler.levels[level]->dof_indices
= std::vector<types::global_dof_index> (next_free_dof,
DoFHandler<dim,spacedim>::invalid_dof_index);
+ dof_handler.levels[level]->cell_dof_indices_cache
+ = std::vector<types::global_dof_index> (cache_size,
+ DoFHandler<dim,spacedim>::invalid_dof_index);
}
// safety check: make sure that
for (unsigned int level=0; level<dof_handler.tria->n_levels(); ++level)
{
dof_handler.levels[level]->dof_offsets
- = std::vector<DoFLevel::offset_type>
- (dof_handler.tria->n_raw_hexs(level),
- (DoFLevel::offset_type)(-1));
+ = std::vector<DoFLevel::offset_type> (
+ dof_handler.tria->n_raw_hexs(level),
+ (DoFLevel::offset_type)(-1));
+ dof_handler.levels[level]->cell_cache_offsets
+ = std::vector<DoFLevel::offset_type> (
+ dof_handler.tria->n_raw_hexs(level),
+ (DoFLevel::offset_type)(-1));
types::global_dof_index next_free_dof = 0;
+ types::global_dof_index cache_size = 0;
for (typename DoFHandler<dim,spacedim>::active_cell_iterator
cell=dof_handler.begin_active(level);
cell!=dof_handler.end_active(level); ++cell)
{
dof_handler.levels[level]->dof_offsets[cell->index()] = next_free_dof;
next_free_dof += cell->get_fe().dofs_per_hex;
+
+ dof_handler.levels[level]->cell_cache_offsets[cell->index()] = cache_size;
+ cache_size += cell->get_fe().dofs_per_cell;
}
dof_handler.levels[level]->dof_indices
= std::vector<types::global_dof_index> (next_free_dof,
DoFHandler<dim,spacedim>::invalid_dof_index);
+ dof_handler.levels[level]->cell_dof_indices_cache
+ = std::vector<types::global_dof_index> (cache_size,
+ DoFHandler<dim,spacedim>::invalid_dof_index);
}
// safety check: make sure that
= std::vector<IndexSet> (1,
number_cache.locally_owned_dofs);
+ // update the cache used for cell dof indices and compress the data on the levels
+ for (active_cell_iterator cell = begin_active();
+ cell != end(); ++cell)
+ cell->update_cell_dof_indices_cache ();
+
for (unsigned int level=0; level<levels.size(); ++level)
levels[level]->compress_data (*finite_elements);
// do the renumbering
renumber_dofs_internal (new_numbers, dealii::internal::int2type<dim>());
+ // update the cache used for cell dof indices
+ for (active_cell_iterator cell = begin_active();
+ cell != end(); ++cell)
+ cell->update_cell_dof_indices_cache ();
+
// now re-compress the dof indices
for (unsigned int level=0; level<levels.size(); ++level)
levels[level]->compress_data (*finite_elements);