00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012 #ifndef __deal2__block_list_h
00013 #define __deal2__block_list_h
00014
00015
00016 #include <deal.II/base/subscriptor.h>
00017 #include <deal.II/base/template_constraints.h>
00018 #include <deal.II/lac/sparsity_pattern.h>
00019 #include <deal.II/fe/fe.h>
00020
00021 #include <vector>
00022
00023 DEAL_II_NAMESPACE_OPEN
00024
00047 class BlockList :
00048 public Subscriptor
00049 {
00050 public:
00052 typedef std::vector<unsigned int> block_container;
00054 typedef block_container::const_iterator const_iterator;
00055
00068 void create_sparsity_pattern(SparsityPattern& sparsity, unsigned int n) const;
00069
00076 void add(unsigned int block, const std::vector<unsigned int>& indices);
00077
00086 void add(unsigned int block,
00087 const std::vector<unsigned int>& indices,
00088 const std::vector<bool>& selected_indices,
00089 unsigned int offset = 0);
00090
00095 void initialize(unsigned int n_blocks);
00096
00112 template <typename ITERATOR>
00113 void initialize(unsigned int n_blocks,
00114 const ITERATOR begin,
00115 const typename identity<ITERATOR>::type end);
00116
00135 template <typename ITERATOR>
00136 void initialize_mg(unsigned int n_blocks,
00137 const ITERATOR begin,
00138 const typename identity<ITERATOR>::type end);
00139
00171 template <typename ITERATOR>
00172 void initialize(unsigned int n_blocks,
00173 const ITERATOR begin,
00174 const typename identity<ITERATOR>::type end,
00175 const std::vector<bool>& selected_dofs,
00176 unsigned int offset = 0);
00208 template <typename ITERATOR>
00209 void initialize_mg(unsigned int n_blocks,
00210 const ITERATOR begin,
00211 const typename identity<ITERATOR>::type end,
00212 const std::vector<bool>& selected_dofs,
00213 unsigned int offset = 0);
00214
00226 template <int dim, typename ITERATOR>
00227 void initialize_vertex_patches_mg(unsigned int n_blocks,
00228 const ITERATOR begin,
00229 const typename identity<ITERATOR>::type end,
00230 const std::vector<bool>& selected_dofs = std::vector<bool>(),
00231 unsigned int offset = 0);
00232
00240 template <int dim, typename ITERATOR>
00241 unsigned int count_vertex_patches(
00242 const ITERATOR begin,
00243 const typename identity<ITERATOR>::type end,
00244 bool same_level_only) const;
00245
00251 template <int dim, typename ITERATOR>
00252 bool cell_generates_vertex_patch(const ITERATOR cell,
00253 bool same_level_only) const;
00254
00258 unsigned int size() const;
00259
00263 unsigned int block_size(unsigned int block) const;
00264
00268 const_iterator begin(unsigned int block) const;
00272 const_iterator end(unsigned int block) const;
00280 unsigned int local_index(unsigned int block, unsigned int index) const;
00281
00282 private:
00286 std::vector<block_container> index_sets;
00287 };
00288
00289
00290 inline
00291 void
00292 BlockList::create_sparsity_pattern(SparsityPattern& sparsity, unsigned int n) const
00293 {
00294 std::vector<unsigned int> sizes(size());
00295 for (unsigned int b=0;b<size();++b)
00296 sizes[b] = block_size(b);
00297
00298 sparsity.reinit(size(), n, sizes);
00299 for (unsigned int b=0;b<size();++b)
00300 {
00301 for (const_iterator i = begin(b); i != end(b); ++i)
00302 sparsity.add(b,*i);
00303 }
00304 sparsity.compress();
00305 }
00306
00307
00308 inline
00309 void
00310 BlockList::add(const unsigned int block, const std::vector<unsigned int>& indices)
00311 {
00312 AssertIndexRange(block, index_sets.size());
00313
00314 for (unsigned int i=0;i<indices.size();++i)
00315 {
00316 const unsigned int k = indices[i];
00317 if (k==numbers::invalid_unsigned_int)
00318 continue;
00319 if (std::find(index_sets[block].begin(), index_sets[block].end(), k)
00320 == index_sets[block].end())
00321 index_sets[block].push_back(k);
00322 }
00323 }
00324
00325
00326 inline
00327 void
00328 BlockList::add(
00329 const unsigned int block,
00330 const std::vector<unsigned int>& indices,
00331 const std::vector<bool>& selected,
00332 unsigned int offset)
00333 {
00334 AssertIndexRange(block, index_sets.size());
00335 AssertDimension(indices.size(), selected.size());
00336
00337 for (unsigned int i=0;i<indices.size();++i)
00338 {
00339 const unsigned int k = indices[i];
00340 if (k==numbers::invalid_unsigned_int)
00341 continue;
00342 if (selected[i] && std::find(index_sets[block].begin(), index_sets[block].end(), k-offset)
00343 == index_sets[block].end())
00344 index_sets[block].push_back(k-offset);
00345 }
00346 }
00347
00348
00349 inline
00350 void
00351 BlockList::initialize(unsigned int n_blocks)
00352 {
00353 index_sets.resize(n_blocks);
00354 }
00355
00356
00357 template <typename ITERATOR>
00358 inline
00359 void
00360 BlockList::initialize(unsigned int n_blocks, const ITERATOR begin, const typename identity<ITERATOR>::type end)
00361 {
00362 index_sets.resize(n_blocks);
00363 std::vector<unsigned int> indices;
00364 unsigned int k = 0;
00365 for (ITERATOR cell = begin; cell != end; ++cell, ++k)
00366 {
00367 indices.resize(cell->get_fe().dofs_per_cell);
00368 cell->get_dof_indices(indices);
00369 add(k, indices);
00370 }
00371 }
00372
00373
00374 template <typename ITERATOR>
00375 inline
00376 void
00377 BlockList::initialize_mg(unsigned int n_blocks, const ITERATOR begin, const typename identity<ITERATOR>::type end)
00378 {
00379 index_sets.resize(n_blocks);
00380 std::vector<unsigned int> indices;
00381 unsigned int k = 0;
00382 for (ITERATOR cell = begin; cell != end; ++cell, ++k)
00383 {
00384 indices.resize(cell->get_fe().dofs_per_cell);
00385 cell->get_mg_dof_indices(indices);
00386 add(k, indices);
00387 }
00388 }
00389
00390
00391 template <typename ITERATOR>
00392 inline
00393 void
00394 BlockList::initialize(
00395 unsigned int n_blocks,
00396 const ITERATOR begin,
00397 const typename identity<ITERATOR>::type end,
00398 const std::vector<bool>& selected_dofs,
00399 unsigned int offset)
00400 {
00401 index_sets.resize(n_blocks);
00402 std::vector<unsigned int> indices;
00403 unsigned int k = 0;
00404 for (ITERATOR cell = begin; cell != end; ++cell, ++k)
00405 {
00406 indices.resize(cell->get_fe().dofs_per_cell);
00407 AssertDimension(selected_dofs.size(), cell->get_fe().dofs_per_cell);
00408
00409 cell->get_dof_indices(indices);
00410 add(k, indices, selected_dofs);
00411 }
00412 }
00413
00414
00415 template <typename ITERATOR>
00416 inline
00417 void
00418 BlockList::initialize_mg(
00419 unsigned int n_blocks,
00420 const ITERATOR begin,
00421 const typename identity<ITERATOR>::type end,
00422 const std::vector<bool>& selected_dofs,
00423 unsigned int offset)
00424 {
00425 index_sets.resize(n_blocks);
00426 std::vector<unsigned int> indices;
00427 unsigned int k = 0;
00428 for (ITERATOR cell = begin; cell != end; ++cell, ++k)
00429 {
00430 indices.resize(cell->get_fe().dofs_per_cell);
00431 AssertDimension(selected_dofs.size(), cell->get_fe().dofs_per_cell);
00432
00433 cell->get_mg_dof_indices(indices);
00434 add(k, indices, selected_dofs, offset);
00435 }
00436 }
00437
00438
00439
00440 template <int dim, typename ITERATOR>
00441 inline
00442 bool
00443 BlockList::cell_generates_vertex_patch(
00444 const ITERATOR cell,
00445 bool same_level_only) const
00446 {
00447 switch(dim)
00448 {
00449 case 3:
00450 if (cell->at_boundary(4)) break;
00451 if (same_level_only && cell->neighbor(4)->level() != cell->level()) break;
00452 if (cell->neighbor(4)->at_boundary(0)) break;
00453 if (same_level_only && cell->neighbor(4)->neighbor(0)->level() != cell->level()) break;
00454 if (cell->neighbor(4)->at_boundary(2)) break;
00455 if (same_level_only && cell->neighbor(4)->neighbor(2)->level() != cell->level()) break;
00456 if (cell->neighbor(4)->neighbor(0)->at_boundary(2)) break;
00457 if (same_level_only && cell->neighbor(4)->neighbor(0)->neighbor(2)->level() != cell->level()) break;
00458
00459 case 2:
00460 if (cell->at_boundary(2)) break;
00461 if (same_level_only && cell->neighbor(2)->level() != cell->level()) break;
00462 if (cell->neighbor(2)->at_boundary(0)) break;
00463 if (same_level_only && cell->neighbor(2)->neighbor(0)->level() != cell->level()) break;
00464 case 1:
00465 if (cell->at_boundary(0)) break;
00466 if (same_level_only && cell->neighbor(0)->level() != cell->level()) break;
00467 return true;
00468 break;
00469 default:
00470 Assert(false, ExcNotImplemented());
00471 break;
00472 }
00473 return false;
00474 }
00475
00476
00477 template <int dim, typename ITERATOR>
00478 inline
00479 unsigned int
00480 BlockList::count_vertex_patches(
00481 const ITERATOR begin,
00482 const typename identity<ITERATOR>::type end,
00483 bool same_level_only) const
00484 {
00485 unsigned int count = 0;
00486 for (ITERATOR cell = begin; cell != end; ++cell)
00487 if (cell_generates_vertex_patch<dim>(cell, same_level_only))
00488 ++count;
00489 return count;
00490 }
00491
00492
00493 template <int dim, typename ITERATOR>
00494 inline
00495 void
00496 BlockList::initialize_vertex_patches_mg(
00497 unsigned int n_blocks,
00498 const ITERATOR begin,
00499 const typename identity<ITERATOR>::type end,
00500 const std::vector<bool>& selected_dofs,
00501 unsigned int offset)
00502 {
00503 Assert(selected_dofs.size() == 0, ExcNotImplemented());
00504 Assert(offset==0, ExcNotImplemented());
00505 const FiniteElement<dim>& fe = begin->get_fe();
00506 Assert(fe.dofs_per_vertex == 0, ExcNotImplemented());
00507
00508 index_sets.resize(n_blocks);
00509 std::vector<unsigned int> indices;
00510 unsigned int k = 0;
00511 for (ITERATOR cell = begin; cell != end; ++cell)
00512 {
00513 if (cell_generates_vertex_patch<dim>(cell, true))
00514 {
00515 indices.resize(cell->get_fe().dofs_per_cell);
00516
00517 switch(dim)
00518 {
00519 case 3:
00520 cell->neighbor(4)->get_mg_dof_indices(indices);
00521 for (unsigned int i=0;i<fe.dofs_per_face;++i)
00522 {
00523 indices[fe.face_to_cell_index(i,1)] = numbers::invalid_unsigned_int;
00524 indices[fe.face_to_cell_index(i,3)] = numbers::invalid_unsigned_int;
00525 indices[fe.face_to_cell_index(i,4)] = numbers::invalid_unsigned_int;
00526 }
00527 add(k, indices);
00528 cell->neighbor(4)->neighbor(0)->get_mg_dof_indices(indices);
00529 for (unsigned int i=0;i<fe.dofs_per_face;++i)
00530 {
00531 indices[fe.face_to_cell_index(i,0)] = numbers::invalid_unsigned_int;
00532 indices[fe.face_to_cell_index(i,3)] = numbers::invalid_unsigned_int;
00533 indices[fe.face_to_cell_index(i,4)] = numbers::invalid_unsigned_int;
00534 }
00535 add(k, indices);
00536 cell->neighbor(4)->neighbor(2)->get_mg_dof_indices(indices);
00537 for (unsigned int i=0;i<fe.dofs_per_face;++i)
00538 {
00539 indices[fe.face_to_cell_index(i,1)] = numbers::invalid_unsigned_int;
00540 indices[fe.face_to_cell_index(i,2)] = numbers::invalid_unsigned_int;
00541 indices[fe.face_to_cell_index(i,4)] = numbers::invalid_unsigned_int;
00542 }
00543 add(k, indices);
00544 cell->neighbor(4)->neighbor(2)->neighbor(0)->get_mg_dof_indices(indices);
00545 for (unsigned int i=0;i<fe.dofs_per_face;++i)
00546 {
00547 indices[fe.face_to_cell_index(i,0)] = numbers::invalid_unsigned_int;
00548 indices[fe.face_to_cell_index(i,2)] = numbers::invalid_unsigned_int;
00549 indices[fe.face_to_cell_index(i,4)] = numbers::invalid_unsigned_int;
00550 }
00551 add(k, indices);
00552 case 2:
00553 cell->neighbor(2)->get_mg_dof_indices(indices);
00554 for (unsigned int i=0;i<fe.dofs_per_face;++i)
00555 {
00556 indices[fe.face_to_cell_index(i,1)] = numbers::invalid_unsigned_int;
00557 indices[fe.face_to_cell_index(i,2)] = numbers::invalid_unsigned_int;
00558 if (dim>2)
00559 indices[fe.face_to_cell_index(i,5)] = numbers::invalid_unsigned_int;
00560 }
00561 add(k, indices);
00562 cell->neighbor(2)->neighbor(0)->get_mg_dof_indices(indices);
00563 for (unsigned int i=0;i<fe.dofs_per_face;++i)
00564 {
00565 indices[fe.face_to_cell_index(i,0)] = numbers::invalid_unsigned_int;
00566 indices[fe.face_to_cell_index(i,2)] = numbers::invalid_unsigned_int;
00567 if (dim>2)
00568 indices[fe.face_to_cell_index(i,5)] = numbers::invalid_unsigned_int;
00569 }
00570 add(k, indices);
00571
00572 case 1:
00573 cell->get_mg_dof_indices(indices);
00574 for (unsigned int i=0;i<fe.dofs_per_face;++i)
00575 {
00576 indices[fe.face_to_cell_index(i,1)] = numbers::invalid_unsigned_int;
00577 if (dim>1)
00578 indices[fe.face_to_cell_index(i,3)] = numbers::invalid_unsigned_int;
00579 if (dim>2)
00580 indices[fe.face_to_cell_index(i,5)] = numbers::invalid_unsigned_int;
00581 }
00582 add(k, indices);
00583 cell->neighbor(0)->get_mg_dof_indices(indices);
00584 for (unsigned int i=0;i<fe.dofs_per_face;++i)
00585 {
00586 indices[fe.face_to_cell_index(i,0)] = numbers::invalid_unsigned_int;
00587 if (dim>1)
00588 indices[fe.face_to_cell_index(i,3)] = numbers::invalid_unsigned_int;
00589 if (dim>2)
00590 indices[fe.face_to_cell_index(i,5)] = numbers::invalid_unsigned_int;
00591 }
00592 add(k, indices);
00593 break;
00594 default:
00595 Assert(false, ExcNotImplemented());
00596 break;
00597 }
00598 ++k;
00599 }
00600 }
00601 }
00602
00603
00604 inline
00605 unsigned int
00606 BlockList::size() const
00607 {
00608 return index_sets.size();
00609 }
00610
00611
00612 inline
00613 unsigned int
00614 BlockList::block_size(unsigned int block) const
00615 {
00616 return index_sets[block].size();
00617 }
00618
00619
00620 inline
00621 BlockList::const_iterator
00622 BlockList::begin(unsigned int block) const
00623 {
00624 AssertIndexRange(block, index_sets.size());
00625 return index_sets[block].begin();
00626 }
00627
00628
00629 inline
00630 BlockList::const_iterator
00631 BlockList::end(unsigned int block) const
00632 {
00633 AssertIndexRange(block, index_sets.size());
00634 return index_sets[block].end();
00635 }
00636
00637
00638 inline
00639 unsigned int
00640 BlockList::local_index(unsigned int block, unsigned int index) const
00641 {
00642 AssertIndexRange(block, index_sets.size());
00643 const block_container& b = index_sets[block];
00644 for (unsigned i=0;i<b.size();++i)
00645 if (b[i] == index)
00646 return i;
00647 return numbers::invalid_unsigned_int;
00648 }
00649
00650
00651 DEAL_II_NAMESPACE_CLOSE
00652
00653 #endif
00654