Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
dof_level.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2003 - 2019 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #include <deal.II/base/memory_consumption.h>
17 
18 #include <deal.II/hp/dof_level.h>
19 #include <deal.II/hp/fe_collection.h>
20 
21 #include <iostream>
22 
23 DEAL_II_NAMESPACE_OPEN
24 
25 namespace internal
26 {
27  namespace hp
28  {
30  static_cast<DoFLevel::active_fe_index_type>(-1);
31 
32 
33 
34  template <int dim, int spacedim>
35  void
37  const ::hp::FECollection<dim, spacedim> &fe_collection)
38  {
39  (void)fe_collection;
40 
41  if (dof_offsets.size() == 0 || dof_indices.size() == 0)
42  return;
43 
44  // in a first run through, count how many new slots we need in the
45  // dof_indices array after compression. note that the 'cell'
46  // counter is incremented inside the loop
47  unsigned int new_size = 0;
48  for (unsigned int cell = 0; cell < dof_offsets.size();)
49  // see if this cell is active on the current level
50  if (dof_offsets[cell] != static_cast<offset_type>(-1))
51  {
52  // find the next cell active on this level
53  unsigned int next_cell = cell + 1;
54  while ((next_cell < dof_offsets.size()) &&
55  (dof_offsets[next_cell] == static_cast<offset_type>(-1)))
56  ++next_cell;
57 
58  const unsigned int next_offset =
59  (next_cell < dof_offsets.size() ? dof_offsets[next_cell] :
60  dof_indices.size());
61 
62  Assert(next_offset - dof_offsets[cell] ==
63  fe_collection[active_fe_indices[cell]]
64  .template n_dofs_per_object<dim>(),
66 
67  // see if the range of dofs for this cell can be compressed and if
68  // so how many slots we have to store for them
69  if (next_offset > dof_offsets[cell])
70  {
71  bool compressible = true;
72  for (unsigned int j = dof_offsets[cell] + 1; j < next_offset;
73  ++j)
74  if (dof_indices[j] != dof_indices[j - 1] + 1)
75  {
76  compressible = false;
77  break;
78  }
79  if (compressible == true)
80  new_size += 1;
81  else
82  new_size += (next_offset - dof_offsets[cell]);
83  }
84 
85  // then move on to the next cell
86  cell = next_cell;
87  }
88  else
89  ++cell;
90 
91  // now allocate the new array and copy into it whatever we need
92  std::vector<types::global_dof_index> new_dof_indices;
93  new_dof_indices.reserve(new_size);
94  std::vector<offset_type> new_dof_offsets(dof_offsets.size(),
95  static_cast<offset_type>(-1));
96  for (unsigned int cell = 0; cell < dof_offsets.size();)
97  // see if this cell is active on the current level
98  if (dof_offsets[cell] != static_cast<offset_type>(-1))
99  {
100  // find the next cell active on this level
101  unsigned int next_cell = cell + 1;
102  while ((next_cell < dof_offsets.size()) &&
103  (dof_offsets[next_cell] == static_cast<offset_type>(-1)))
104  ++next_cell;
105 
106  const unsigned int next_offset =
107  (next_cell < dof_offsets.size() ? dof_offsets[next_cell] :
108  dof_indices.size());
109 
110  Assert(next_offset - dof_offsets[cell] ==
111  fe_collection[active_fe_indices[cell]]
112  .template n_dofs_per_object<dim>(),
113  ExcInternalError());
114 
115  new_dof_offsets[cell] = new_dof_indices.size();
116 
117  // see if the range of dofs for this cell can be compressed and if
118  // so how many slots we have to store for them
119  if (next_offset > dof_offsets[cell])
120  {
121  bool compressible = true;
122  for (unsigned int j = dof_offsets[cell] + 1; j < next_offset;
123  ++j)
124  if (dof_indices[j] != dof_indices[j - 1] + 1)
125  {
126  compressible = false;
127  break;
128  }
129 
130  // if this cell is compressible, then copy the first index and
131  // mark this in the dof_offsets array
132  if (compressible == true)
133  {
134  new_dof_indices.push_back(dof_indices[dof_offsets[cell]]);
135 
136  // make sure that the current active_fe_index indicates
137  // that this entry hasn't been compressed yet
139  false,
140  ExcInternalError());
141 
142  // then mark the compression
143  active_fe_indices[cell] =
145  }
146  else
147  for (unsigned int i = dof_offsets[cell]; i < next_offset; ++i)
148  new_dof_indices.push_back(dof_indices[i]);
149  }
150 
151  // then move on to the next cell
152  cell = next_cell;
153  }
154  else
155  ++cell;
156 
157  // finally swap old and new content
158  Assert(new_dof_indices.size() == new_size, ExcInternalError());
159  dof_indices.swap(new_dof_indices);
160  dof_offsets.swap(new_dof_offsets);
161  }
162 
163 
164 
165  template <int dim, int spacedim>
166  void
168  const ::hp::FECollection<dim, spacedim> &fe_collection)
169  {
170  if (dof_offsets.size() == 0 || dof_indices.size() == 0)
171  return;
172 
173  // in a first run through, count how many new slots we need in the
174  // dof_indices array after uncompression.
175  unsigned int new_size = 0;
176  for (unsigned int cell = 0; cell < dof_offsets.size(); ++cell)
177  if (dof_offsets[cell] != static_cast<offset_type>(-1))
178  {
179  // we know now that the slot for this cell is used. extract the
180  // active_fe_index for it and see how many entries we need
181  new_size += fe_collection[active_fe_index(cell)]
182  .template n_dofs_per_object<dim>();
183  }
184 
185  // now allocate the new array and copy into it whatever we need
186  std::vector<types::global_dof_index> new_dof_indices;
187  new_dof_indices.reserve(new_size);
188  std::vector<offset_type> new_dof_offsets(dof_offsets.size(),
189  static_cast<offset_type>(-1));
190  for (unsigned int cell = 0; cell < dof_offsets.size();)
191  // see if this cell is active on the current level
192  if (dof_offsets[cell] != static_cast<offset_type>(-1))
193  {
194  // find the next cell active on this level
195  unsigned int next_cell = cell + 1;
196  while ((next_cell < dof_offsets.size()) &&
197  (dof_offsets[next_cell] == static_cast<offset_type>(-1)))
198  ++next_cell;
199 
200  const unsigned int next_offset =
201  (next_cell < dof_offsets.size() ? dof_offsets[next_cell] :
202  dof_indices.size());
203 
204  // set offset for this cell
205  new_dof_offsets[cell] = new_dof_indices.size();
206 
207  // see if we need to uncompress this set of dofs
208  if (is_compressed_entry(active_fe_indices[cell]) == false)
209  {
210  // apparently not. simply copy them
211  Assert(next_offset - dof_offsets[cell] ==
212  fe_collection[active_fe_indices[cell]]
213  .template n_dofs_per_object<dim>(),
214  ExcInternalError());
215  for (unsigned int i = dof_offsets[cell]; i < next_offset; ++i)
216  new_dof_indices.push_back(dof_indices[i]);
217  }
218  else
219  {
220  // apparently so. uncompress
221  Assert(next_offset - dof_offsets[cell] == 1,
222  ExcInternalError());
223  const unsigned int dofs_per_object =
224  fe_collection[get_toggled_compression_state(
225  active_fe_indices[cell])]
226  .template n_dofs_per_object<dim>();
227  for (unsigned int i = 0; i < dofs_per_object; ++i)
228  new_dof_indices.push_back(dof_indices[dof_offsets[cell]] + i);
229 
230  // then mark the uncompression
231  active_fe_indices[cell] =
233  }
234 
235  // then move on to the next cell
236  cell = next_cell;
237  }
238  else
239  ++cell;
240 
241  // verify correct size, then swap arrays
242  Assert(new_dof_indices.size() == new_size, ExcInternalError());
243  dof_indices.swap(new_dof_indices);
244  dof_offsets.swap(new_dof_offsets);
245  }
246 
247 
248 
249  std::size_t
251  {
258  }
259 
260 
261 
262  void
264  {
265  for (auto &active_fe_index : active_fe_indices)
268  }
269 
270 
271 
272  // explicit instantiations
273  template void
274  DoFLevel::compress_data(const ::hp::FECollection<1, 1> &);
275  template void
276  DoFLevel::compress_data(const ::hp::FECollection<1, 2> &);
277  template void
278  DoFLevel::compress_data(const ::hp::FECollection<1, 3> &);
279  template void
280  DoFLevel::compress_data(const ::hp::FECollection<2, 2> &);
281  template void
282  DoFLevel::compress_data(const ::hp::FECollection<2, 3> &);
283  template void
284  DoFLevel::compress_data(const ::hp::FECollection<3, 3> &);
285 
286  template void
287  DoFLevel::uncompress_data(const ::hp::FECollection<1, 1> &);
288  template void
289  DoFLevel::uncompress_data(const ::hp::FECollection<1, 2> &);
290  template void
291  DoFLevel::uncompress_data(const ::hp::FECollection<1, 3> &);
292  template void
293  DoFLevel::uncompress_data(const ::hp::FECollection<2, 2> &);
294  template void
295  DoFLevel::uncompress_data(const ::hp::FECollection<2, 3> &);
296  template void
297  DoFLevel::uncompress_data(const ::hp::FECollection<3, 3> &);
298  } // namespace hp
299 } // namespace internal
300 
301 DEAL_II_NAMESPACE_CLOSE
void uncompress_data(const ::hp::FECollection< dim, spacedim > &fe_collection)
Definition: dof_level.cc:167
std::size_t memory_consumption() const
Definition: dof_level.cc:250
static active_fe_index_type get_toggled_compression_state(const active_fe_index_type active_fe_index)
Definition: dof_level.h:393
unsigned int offset_type
Definition: dof_level.h:112
std::vector< offset_type > cell_cache_offsets
Definition: dof_level.h:198
void compress_data(const ::hp::FECollection< dim, spacedim > &fe_collection)
Definition: dof_level.cc:36
unsigned short int active_fe_index_type
Definition: dof_level.h:117
#define Assert(cond, exc)
Definition: exceptions.h:1407
Definition: hp.h:117
std::vector< offset_type > dof_offsets
Definition: dof_level.h:186
static bool is_compressed_entry(const active_fe_index_type active_fe_index)
Definition: dof_level.h:385
void normalize_active_fe_indices()
Definition: dof_level.cc:263
std::vector< types::global_dof_index > dof_indices
Definition: dof_level.h:193
std::vector< types::global_dof_index > cell_dof_indices_cache
Definition: dof_level.h:205
unsigned int active_fe_index(const unsigned int obj_index) const
Definition: dof_level.h:465
static const active_fe_index_type invalid_active_fe_index
Definition: dof_level.h:129
std::vector< active_fe_index_type > future_fe_indices
Definition: dof_level.h:172
std::vector< active_fe_index_type > active_fe_indices
Definition: dof_level.h:158
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
static ::ExceptionBase & ExcInternalError()