Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
block_list.h
1 // ---------------------------------------------------------------------
2 // @f$Id: block_list.h 31932 2013-12-08 02:15:54Z heister @f$
3 //
4 // Copyright (C) 2010 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__block_list_h
18 #define __deal2__block_list_h
19 
20 
21 #include <deal.II/base/subscriptor.h>
22 #include <deal.II/base/template_constraints.h>
23 #include <deal.II/lac/sparsity_pattern.h>
24 #include <deal.II/fe/fe.h>
25 
26 #include <vector>
27 
29 
52 class BlockList :
53  public Subscriptor
54 {
55 public:
60 
62  typedef std::vector<size_type> block_container;
64  typedef block_container::const_iterator const_iterator;
65 
78  void create_sparsity_pattern(SparsityPattern &sparsity, size_type n) const;
79 
86  void add(size_type block, const std::vector<size_type> &indices);
87 
96  void add(size_type block,
97  const std::vector<size_type> &indices,
98  const std::vector<bool> &selected_indices,
99  size_type offset = 0);
100 
105  void initialize(size_type n_blocks);
106 
122  template <typename ITERATOR>
123  void initialize(size_type n_blocks,
124  const ITERATOR begin,
125  const typename identity<ITERATOR>::type end);
126 
145  template <typename ITERATOR>
146  void initialize_mg(size_type n_blocks,
147  const ITERATOR begin,
148  const typename identity<ITERATOR>::type end) DEAL_II_DEPRECATED;
149 
181  template <typename ITERATOR>
182  void initialize(size_type n_blocks,
183  const ITERATOR begin,
184  const typename identity<ITERATOR>::type end,
185  const std::vector<bool> &selected_dofs,
186  size_type offset = 0) DEAL_II_DEPRECATED;
218  template <typename ITERATOR>
219  void initialize_mg(size_type n_blocks,
220  const ITERATOR begin,
221  const typename identity<ITERATOR>::type end,
222  const std::vector<bool> &selected_dofs,
223  size_type offset = 0) DEAL_II_DEPRECATED;
224 
236  template <int dim, typename ITERATOR>
237  void initialize_vertex_patches_mg(size_type n_blocks,
238  const ITERATOR begin,
239  const typename identity<ITERATOR>::type end,
240  const std::vector<bool> &selected_dofs = std::vector<bool>(),
241  size_type offset = 0) DEAL_II_DEPRECATED;
242 
250  template <int dim, typename ITERATOR>
251  unsigned int count_vertex_patches(
252  const ITERATOR begin,
253  const typename identity<ITERATOR>::type end,
254  bool same_level_only) const DEAL_II_DEPRECATED;
255 
261  template <int dim, typename ITERATOR>
262  bool cell_generates_vertex_patch(const ITERATOR cell,
263  bool same_level_only) const DEAL_II_DEPRECATED;
264 
268  size_type size() const;
269 
273  size_type block_size(size_type block) const;
274 
278  const_iterator begin(size_type block) const;
282  const_iterator end(size_type block) const;
290  size_type local_index(size_type block, size_type index) const;
291 
292 private:
296  std::vector<block_container> index_sets;
297 } DEAL_II_DEPRECATED;
298 
299 
300 inline
301 void
302 BlockList::create_sparsity_pattern(SparsityPattern &sparsity, size_type n) const
303 {
304  std::vector<unsigned int> sizes(size());
305  for (size_type b=0; b<size(); ++b)
306  sizes[b] = block_size(b);
307 
308  sparsity.reinit(size(), n, sizes);
309  for (size_type b=0; b<size(); ++b)
310  {
311  for (const_iterator i = begin(b); i != end(b); ++i)
312  sparsity.add(b,*i);
313  }
314  sparsity.compress();
315 }
316 
317 
318 inline
319 void
320 BlockList::add(const size_type block, const std::vector<size_type> &indices)
321 {
322  AssertIndexRange(block, index_sets.size());
323 
324  for (size_type i=0; i<indices.size(); ++i)
325  {
326  const size_type k = indices[i];
328  continue;
329  if (std::find(index_sets[block].begin(), index_sets[block].end(), k)
330  == index_sets[block].end())
331  index_sets[block].push_back(k);
332  }
333 }
334 
335 
336 inline
337 void
339  const size_type block,
340  const std::vector<size_type> &indices,
341  const std::vector<bool> &selected,
342  size_type offset)
343 {
344  AssertIndexRange(block, index_sets.size());
345  AssertDimension(indices.size(), selected.size());
346 
347  for (size_type i=0; i<indices.size(); ++i)
348  {
349  const size_type k = indices[i];
351  continue;
352  if (selected[i] && std::find(index_sets[block].begin(), index_sets[block].end(), k-offset)
353  == index_sets[block].end())
354  index_sets[block].push_back(k-offset);
355  }
356 }
357 
358 
359 inline
360 void
361 BlockList::initialize(size_type n_blocks)
362 {
363  index_sets.resize(n_blocks);
364 }
365 
366 
367 template <typename ITERATOR>
368 inline
369 void
370 BlockList::initialize(size_type n_blocks, const ITERATOR begin, const typename identity<ITERATOR>::type end)
371 {
372  index_sets.resize(n_blocks);
373  std::vector<size_type> indices;
374  size_type k = 0;
375  for (ITERATOR cell = begin; cell != end; ++cell, ++k)
376  {
377  indices.resize(cell->get_fe().dofs_per_cell);
378  cell->get_dof_indices(indices);
379  add(k, indices);
380  }
381 }
382 
383 
384 template <typename ITERATOR>
385 inline
386 void
387 BlockList::initialize_mg(size_type n_blocks, const ITERATOR begin, const typename identity<ITERATOR>::type end)
388 {
389  index_sets.resize(n_blocks);
390  std::vector<size_type> indices;
391  size_type k = 0;
392  for (ITERATOR cell = begin; cell != end; ++cell, ++k)
393  {
394  indices.resize(cell->get_fe().dofs_per_cell);
395  cell->get_mg_dof_indices(indices);
396  add(k, indices);
397  }
398 }
399 
400 
401 template <typename ITERATOR>
402 inline
403 void
405  size_type n_blocks,
406  const ITERATOR begin,
407  const typename identity<ITERATOR>::type end,
408  const std::vector<bool> &selected_dofs,
409  size_type offset)
410 {
411  index_sets.resize(n_blocks);
412  std::vector<size_type> indices;
413  size_type k = 0;
414  for (ITERATOR cell = begin; cell != end; ++cell, ++k)
415  {
416  indices.resize(cell->get_fe().dofs_per_cell);
417  AssertDimension(selected_dofs.size(), cell->get_fe().dofs_per_cell);
418 
419  cell->get_dof_indices(indices);
420  add(k, indices, selected_dofs);
421  }
422 }
423 
424 
425 template <typename ITERATOR>
426 inline
427 void
429  size_type n_blocks,
430  const ITERATOR begin,
431  const typename identity<ITERATOR>::type end,
432  const std::vector<bool> &selected_dofs,
433  size_type offset)
434 {
435  index_sets.resize(n_blocks);
436  std::vector<size_type> indices;
437  size_type k = 0;
438  for (ITERATOR cell = begin; cell != end; ++cell, ++k)
439  {
440  indices.resize(cell->get_fe().dofs_per_cell);
441  AssertDimension(selected_dofs.size(), cell->get_fe().dofs_per_cell);
442 
443  cell->get_mg_dof_indices(indices);
444  add(k, indices, selected_dofs, offset);
445  }
446 }
447 
448 
449 
450 template <int dim, typename ITERATOR>
451 inline
452 bool
454  const ITERATOR cell,
455  bool same_level_only) const
456 {
457  switch (dim)
458  {
459  case 3:
460  if (cell->at_boundary(4)) break;
461  if (same_level_only && cell->neighbor(4)->level() != cell->level()) break;
462  if (cell->neighbor(4)->at_boundary(0)) break;
463  if (same_level_only && cell->neighbor(4)->neighbor(0)->level() != cell->level()) break;
464  if (cell->neighbor(4)->at_boundary(2)) break;
465  if (same_level_only && cell->neighbor(4)->neighbor(2)->level() != cell->level()) break;
466  if (cell->neighbor(4)->neighbor(0)->at_boundary(2)) break;
467  if (same_level_only && cell->neighbor(4)->neighbor(0)->neighbor(2)->level() != cell->level()) break;
468  // No break here
469  case 2:
470  if (cell->at_boundary(2)) break;
471  if (same_level_only && cell->neighbor(2)->level() != cell->level()) break;
472  if (cell->neighbor(2)->at_boundary(0)) break;
473  if (same_level_only && cell->neighbor(2)->neighbor(0)->level() != cell->level()) break;
474  case 1:
475  if (cell->at_boundary(0)) break;
476  if (same_level_only && cell->neighbor(0)->level() != cell->level()) break;
477  return true;
478  break;
479  default:
480  Assert(false, ExcNotImplemented());
481  break;
482  }
483  return false;
484 }
485 
486 
487 template <int dim, typename ITERATOR>
488 inline
489 unsigned int
491  const ITERATOR begin,
492  const typename identity<ITERATOR>::type end,
493  bool same_level_only) const
494 {
495  unsigned int count = 0;
496  for (ITERATOR cell = begin; cell != end; ++cell)
497  if (cell_generates_vertex_patch<dim>(cell, same_level_only))
498  ++count;
499  return count;
500 }
501 
502 
503 template <int dim, typename ITERATOR>
504 inline
505 void
507  size_type n_blocks,
508  const ITERATOR begin,
509  const typename identity<ITERATOR>::type end,
510  const std::vector<bool> &selected_dofs,
511  size_type offset)
512 {
513  Assert(selected_dofs.size() == 0, ExcNotImplemented());
514  Assert(offset==0, ExcNotImplemented());
515  const FiniteElement<dim> &fe = begin->get_fe();
517 
518  index_sets.resize(n_blocks);
519  std::vector<size_type> indices;
520  size_type k = 0;
521  for (ITERATOR cell = begin; cell != end; ++cell)
522  {
523  if (cell_generates_vertex_patch<dim>(cell, true))
524  {
525  indices.resize(cell->get_fe().dofs_per_cell);
526 
527  switch (dim)
528  {
529  case 3:
530  cell->neighbor(4)->get_mg_dof_indices(indices);
531  for (unsigned int i=0; i<fe.dofs_per_face; ++i)
532  {
536  }
537  add(k, indices);
538  cell->neighbor(4)->neighbor(0)->get_mg_dof_indices(indices);
539  for (unsigned int i=0; i<fe.dofs_per_face; ++i)
540  {
544  }
545  add(k, indices);
546  cell->neighbor(4)->neighbor(2)->get_mg_dof_indices(indices);
547  for (unsigned int i=0; i<fe.dofs_per_face; ++i)
548  {
552  }
553  add(k, indices);
554  cell->neighbor(4)->neighbor(2)->neighbor(0)->get_mg_dof_indices(indices);
555  for (unsigned int i=0; i<fe.dofs_per_face; ++i)
556  {
560  }
561  add(k, indices);
562  case 2:
563  cell->neighbor(2)->get_mg_dof_indices(indices);
564  for (unsigned int i=0; i<fe.dofs_per_face; ++i)
565  {
568  if (dim>2)
570  }
571  add(k, indices);
572  cell->neighbor(2)->neighbor(0)->get_mg_dof_indices(indices);
573  for (unsigned int i=0; i<fe.dofs_per_face; ++i)
574  {
577  if (dim>2)
579  }
580  add(k, indices);
581  // no break here
582  case 1:
583  cell->get_mg_dof_indices(indices);
584  for (unsigned int i=0; i<fe.dofs_per_face; ++i)
585  {
587  if (dim>1)
589  if (dim>2)
591  }
592  add(k, indices);
593  cell->neighbor(0)->get_mg_dof_indices(indices);
594  for (unsigned int i=0; i<fe.dofs_per_face; ++i)
595  {
597  if (dim>1)
599  if (dim>2)
601  }
602  add(k, indices);
603  break;
604  default:
605  Assert(false, ExcNotImplemented());
606  break;
607  }
608  ++k;
609  }
610  }
611 }
612 
613 
614 inline
617 {
618  return index_sets.size();
619 }
620 
621 
622 inline
624 BlockList::block_size(size_type block) const
625 {
626  return index_sets[block].size();
627 }
628 
629 
630 inline
632 BlockList::begin(size_type block) const
633 {
634  AssertIndexRange(block, index_sets.size());
635  return index_sets[block].begin();
636 }
637 
638 
639 inline
641 BlockList::end(size_type block) const
642 {
643  AssertIndexRange(block, index_sets.size());
644  return index_sets[block].end();
645 }
646 
647 
648 inline
650 BlockList::local_index(size_type block, size_type index) const
651 {
652  AssertIndexRange(block, index_sets.size());
653  const block_container &b = index_sets[block];
654  for (size_type i=0; i<b.size(); ++i)
655  if (b[i] == index)
656  return i;
658 }
659 
660 
661 DEAL_II_NAMESPACE_CLOSE
662 
663 #endif
const types::global_dof_index invalid_size_type
Definition: types.h:211
const_iterator begin(size_type block) const
Definition: block_list.h:632
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:858
size_type local_index(size_type block, size_type index) const
Definition: block_list.h:650
const_iterator end(size_type block) const
Definition: block_list.h:641
bool cell_generates_vertex_patch(const ITERATOR cell, bool same_level_only) const DEAL_II_DEPRECATED
Definition: block_list.h:453
#define AssertIndexRange(index, range)
Definition: exceptions.h:888
unsigned int count_vertex_patches(const ITERATOR begin, const typename identity< ITERATOR >::type end, bool same_level_only) const DEAL_II_DEPRECATED
Definition: block_list.h:490
STL namespace.
void add(size_type block, const std::vector< size_type > &indices)
Definition: block_list.h:320
void initialize(size_type n_blocks)
Definition: block_list.h:361
virtual unsigned int face_to_cell_index(const unsigned int face_dof_index, const unsigned int face, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false) const
const unsigned int dofs_per_face
Definition: fe_base.h:264
unsigned int global_dof_index
Definition: types.h:100
#define Assert(cond, exc)
Definition: exceptions.h:299
size_type size() const
Definition: block_list.h:616
std::vector< size_type > block_container
The container for each index set.
Definition: block_list.h:62
BlockCompressedSparsityPattern CompressedBlockSparsityPattern DEAL_II_DEPRECATED
void create_sparsity_pattern(SparsityPattern &sparsity, size_type n) const
Definition: block_list.h:302
const unsigned int dofs_per_vertex
Definition: fe_base.h:214
std::vector< block_container > index_sets
Definition: block_list.h:296
size_type block_size(size_type block) const
Definition: block_list.h:624
block_container::const_iterator const_iterator
The iterator for individual indices.
Definition: block_list.h:64
::ExceptionBase & ExcNotImplemented()
types::global_dof_index size_type
Definition: block_list.h:59
void initialize_mg(size_type n_blocks, const ITERATOR begin, const typename identity< ITERATOR >::type end) DEAL_II_DEPRECATED
Definition: block_list.h:387
const types::global_dof_index invalid_dof_index
Definition: types.h:217
void initialize_vertex_patches_mg(size_type n_blocks, const ITERATOR begin, const typename identity< ITERATOR >::type end, const std::vector< bool > &selected_dofs=std::vector< bool >(), size_type offset=0) DEAL_II_DEPRECATED
Definition: block_list.h:506