16 #include <deal.II/base/memory_consumption.h> 18 #include <deal.II/hp/dof_level.h> 19 #include <deal.II/hp/fe_collection.h> 23 DEAL_II_NAMESPACE_OPEN
29 template <
int dim,
int spacedim>
32 const ::hp::FECollection<dim, spacedim> &fe_collection)
42 unsigned int new_size = 0;
43 for (
unsigned int cell = 0; cell <
dof_offsets.size();)
48 unsigned int next_cell = cell + 1;
53 const unsigned int next_offset =
59 .
template n_dofs_per_object<dim>(),
66 bool compressible =
true;
67 for (
unsigned int j =
dof_offsets[cell] + 1; j < next_offset;
74 if (compressible ==
true)
87 std::vector<types::global_dof_index> new_dof_indices;
88 new_dof_indices.reserve(new_size);
89 std::vector<offset_type> new_dof_offsets(
dof_offsets.size(),
91 for (
unsigned int cell = 0; cell <
dof_offsets.size();)
96 unsigned int next_cell = cell + 1;
101 const unsigned int next_offset =
107 .
template n_dofs_per_object<dim>(),
110 new_dof_offsets[cell] = new_dof_indices.size();
116 bool compressible =
true;
117 for (
unsigned int j =
dof_offsets[cell] + 1; j < next_offset;
121 compressible =
false;
127 if (compressible ==
true)
142 for (
unsigned int i =
dof_offsets[cell]; i < next_offset; ++i)
160 template <
int dim,
int spacedim>
163 const ::hp::FECollection<dim, spacedim> &fe_collection)
170 unsigned int new_size = 0;
171 for (
unsigned int cell = 0; cell <
dof_offsets.size(); ++cell)
177 .template n_dofs_per_object<dim>();
181 std::vector<types::global_dof_index> new_dof_indices;
182 new_dof_indices.reserve(new_size);
183 std::vector<offset_type> new_dof_offsets(
dof_offsets.size(),
185 for (
unsigned int cell = 0; cell <
dof_offsets.size();)
190 unsigned int next_cell = cell + 1;
195 const unsigned int next_offset =
200 new_dof_offsets[cell] = new_dof_indices.size();
208 .
template n_dofs_per_object<dim>(),
210 for (
unsigned int i =
dof_offsets[cell]; i < next_offset; ++i)
218 const unsigned int dofs_per_object =
221 .template n_dofs_per_object<dim>();
222 for (
unsigned int i = 0; i < dofs_per_object; ++i)
296 DEAL_II_NAMESPACE_CLOSE
void uncompress_data(const ::hp::FECollection< dim, spacedim > &fe_collection)
static active_fe_index_type get_toggled_compression_state(const active_fe_index_type active_fe_index)
std::vector< offset_type > cell_cache_offsets
void compress_data(const ::hp::FECollection< dim, spacedim > &fe_collection)
#define Assert(cond, exc)
std::vector< offset_type > dof_offsets
static bool is_compressed_entry(const active_fe_index_type active_fe_index)
void normalize_active_fe_indices()
std::vector< types::global_dof_index > dof_indices
std::size_t memory_consumption() const
std::vector< types::global_dof_index > cell_dof_indices_cache
std::vector< active_fe_index_type > active_fe_indices
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
static::ExceptionBase & ExcInternalError()
unsigned int active_fe_index(const unsigned int obj_index) const