Reference documentation for deal.II version 8.1.0
dof_renumbering.h
1 // ---------------------------------------------------------------------
2 // @f$Id: dof_renumbering.h 31766 2013-11-22 20:22:01Z heister @f$
3 //
4 // Copyright (C) 2003 - 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__dof_renumbering_h
18 #define __deal2__dof_renumbering_h
19 
20 
21 #include <deal.II/base/config.h>
23 #include <deal.II/base/point.h>
24 #include <deal.II/dofs/dof_handler.h>
25 #include <deal.II/hp/dof_handler.h>
26 #include <deal.II/multigrid/mg_dof_handler.h>
27 
28 #include <vector>
29 
31 
392 namespace DoFRenumbering
393 {
403  template <class Iterator, int dim>
405  {
410  :
411  dir(dir)
412  {}
416  bool operator () (const Iterator &c1, const Iterator &c2) const
417  {
418  const Point<dim> diff = c2->center() - c1->center();
419  return (diff*dir > 0);
420  }
421 
422  private:
427  };
428 
429 
437  template <int dim>
439  {
444  :
445  dir(dir)
446  {}
451  const std::pair<Point<dim>,types::global_dof_index> &c2) const
452  {
453  const Point<dim> diff = c2.first-c1.first;
454  return (diff*dir > 0 || (diff*dir==0 && c1.second<c2.second));
455  }
456 
457  private:
462  };
463 
477  namespace boost
478  {
490  template <class DH>
491  void
492  Cuthill_McKee (DH &dof_handler,
493  const bool reversed_numbering = false,
494  const bool use_constraints = false);
495 
501  template <class DH>
502  void
503  compute_Cuthill_McKee (std::vector<types::global_dof_index> &new_dof_indices,
504  const DH &,
505  const bool reversed_numbering = false,
506  const bool use_constraints = false);
507 
520  template <class DH>
521  void
522  king_ordering (DH &dof_handler,
523  const bool reversed_numbering = false,
524  const bool use_constraints = false);
525 
530  template <class DH>
531  void
532  compute_king_ordering (std::vector<types::global_dof_index> &new_dof_indices,
533  const DH &,
534  const bool reversed_numbering = false,
535  const bool use_constraints = false);
536 
548  template <class DH>
549  void
550  minimum_degree (DH &dof_handler,
551  const bool reversed_numbering = false,
552  const bool use_constraints = false);
553 
558  template <class DH>
559  void
560  compute_minimum_degree (std::vector<types::global_dof_index> &new_dof_indices,
561  const DH &,
562  const bool reversed_numbering = false,
563  const bool use_constraints = false);
564  }
565 
577  template <class DH>
578  void
579  Cuthill_McKee (DH &dof_handler,
580  const bool reversed_numbering = false,
581  const bool use_constraints = false,
582  const std::vector<types::global_dof_index> &starting_indices = std::vector<types::global_dof_index>());
583 
589  template <class DH>
590  void
591  compute_Cuthill_McKee (std::vector<types::global_dof_index> &new_dof_indices,
592  const DH &,
593  const bool reversed_numbering = false,
594  const bool use_constraints = false,
595  const std::vector<types::global_dof_index> &starting_indices = std::vector<types::global_dof_index>());
596 
610  template <class DH>
611  void
612  Cuthill_McKee (DH &dof_handler,
613  const unsigned int level,
614  const bool reversed_numbering = false,
615  const std::vector<types::global_dof_index> &starting_indices = std::vector<types::global_dof_index> ());
616 
652  template <int dim, int spacedim>
653  void
655  const std::vector<unsigned int> &target_component
656  = std::vector<unsigned int>());
657 
658 
663  template <int dim>
664  void
665  component_wise (hp::DoFHandler<dim> &dof_handler,
666  const std::vector<unsigned int> &target_component = std::vector<unsigned int> ());
667 
674  template <class DH>
675  void
676  component_wise (DH &dof_handler,
677  const unsigned int level,
678  const std::vector<unsigned int> &target_component = std::vector<unsigned int>());
679 
680 
686  template <int dim>
687  void
688  component_wise (MGDoFHandler<dim> &dof_handler,
689  const std::vector<unsigned int> &target_component = std::vector<unsigned int>());
690 
696  template <int dim, int spacedim, class ITERATOR, class ENDITERATOR>
698  compute_component_wise (std::vector<types::global_dof_index> &new_dof_indices,
699  const ITERATOR &start,
700  const ENDITERATOR &end,
701  const std::vector<unsigned int> &target_component,
702  bool is_level_operation);
703 
720  template <int dim, int spacedim>
721  void
722  block_wise (DoFHandler<dim,spacedim> &dof_handler);
723 
724 
740  template <int dim>
741  void
742  block_wise (hp::DoFHandler<dim> &dof_handler);
743 
750  template <int dim>
751  void
752  block_wise (MGDoFHandler<dim> &dof_handler,
753  const unsigned int level);
754 
755 
761  template <int dim>
762  void
763  block_wise (MGDoFHandler<dim> &dof_handler);
764 
770  template <int dim, int spacedim, class ITERATOR, class ENDITERATOR>
772  compute_block_wise (std::vector<types::global_dof_index> &new_dof_indices,
773  const ITERATOR &start,
774  const ENDITERATOR &end);
775 
791  template <int dim>
792  void
793  hierarchical (DoFHandler<dim> &dof_handler);
794 
802  template <class DH>
803  void
804  cell_wise (DH &dof_handler,
805  const std::vector<typename DH::active_cell_iterator> &cell_order);
806 
812  template <class DH>
813  void
814  compute_cell_wise (std::vector<types::global_dof_index> &renumbering,
815  std::vector<types::global_dof_index> &inverse_renumbering,
816  const DH &dof_handler,
817  const std::vector<typename DH::active_cell_iterator> &cell_order);
818 
823  template <class DH>
824  void
825  cell_wise (DH &dof_handler,
826  const unsigned int level,
827  const std::vector<typename DH::level_cell_iterator> &cell_order);
828 
834  template <class DH>
835  void
836  compute_cell_wise (std::vector<types::global_dof_index> &renumbering,
837  std::vector<types::global_dof_index> &inverse_renumbering,
838  const DH &dof_handler,
839  const unsigned int level,
840  const std::vector<typename DH::level_cell_iterator> &cell_order);
841 
876  template <class DH>
877  void
878  downstream (DH &dof_handler,
879  const Point<DH::space_dimension> &direction,
880  const bool dof_wise_renumbering = false);
881 
882 
887  template <class DH>
888  void
889  downstream (DH &dof_handler,
890  const unsigned int level,
891  const Point<DH::space_dimension> &direction,
892  const bool dof_wise_renumbering = false);
893 
897  template <class DH>
898  void
899  downstream_dg (DH &dof,
901 
902  template <class DH>
903  void
904  downstream_dg (DH &dof,
905  const Point<DH::space_dimension> &direction)
906  {
907  downstream(dof, direction);
908  }
909 
910 
914  template <class DH>
915  void
916  downstream_dg (DH &dof,
917  unsigned int level,
919 
920  template <class DH>
921  void
922  downstream_dg (DH &dof,
923  unsigned int level,
924  const Point<DH::space_dimension> &direction)
925  {
926  downstream(dof, level, direction);
927  }
928 
934  template <class DH>
935  void
936  compute_downstream (std::vector<types::global_dof_index> &new_dof_indices,
937  std::vector<types::global_dof_index> &reverse,
938  const DH &dof_handler,
939  const Point<DH::space_dimension> &direction,
940  const bool dof_wise_renumbering);
941 
947  template <class DH>
948  void
949  compute_downstream (std::vector<types::global_dof_index> &new_dof_indices,
950  std::vector<types::global_dof_index> &reverse,
951  const DH &dof_handler,
952  const unsigned int level,
953  const Point<DH::space_dimension> &direction,
954  const bool dof_wise_renumbering);
955 
964  template <class DH>
965  void
966  clockwise_dg (DH &dof_handler,
967  const Point<DH::space_dimension> &center,
968  const bool counter = false);
969 
974  template <class DH>
975  void
976  clockwise_dg (DH &dof_handler,
977  const unsigned int level,
978  const Point<DH::space_dimension> &center,
979  const bool counter = false);
980 
986  template <class DH>
987  void
988  compute_clockwise_dg (std::vector<types::global_dof_index> &new_dof_indices,
989  const DH &dof_handler,
990  const Point<DH::space_dimension> &center,
991  const bool counter);
992 
1011  template <class DH>
1012  void
1013  sort_selected_dofs_back (DH &dof_handler,
1014  const std::vector<bool> &selected_dofs);
1015 
1026  template <class DH>
1027  void
1028  sort_selected_dofs_back (DH &dof_handler,
1029  const std::vector<bool> &selected_dofs,
1030  const unsigned int level);
1031 
1040  template <class DH>
1041  void
1042  compute_sort_selected_dofs_back (std::vector<types::global_dof_index> &new_dof_indices,
1043  const DH &dof_handler,
1044  const std::vector<bool> &selected_dofs);
1045 
1054  template <class DH>
1055  void
1056  compute_sort_selected_dofs_back (std::vector<types::global_dof_index> &new_dof_indices,
1057  const DH &dof_handler,
1058  const std::vector<bool> &selected_dofs,
1059  const unsigned int level);
1060 
1064  template <class DH>
1065  void
1066  random (DH &dof_handler);
1067 
1073  template <class DH>
1074  void
1075  compute_random (std::vector<types::global_dof_index> &new_dof_indices,
1076  const DH &dof_handler);
1077 
1105  template <class DH>
1106  void
1107  subdomain_wise (DH &dof_handler);
1108 
1114  template <class DH>
1115  void
1116  compute_subdomain_wise (std::vector<types::global_dof_index> &new_dof_indices,
1117  const DH &dof_handler);
1118 
1128  DeclException0 (ExcRenumberingIncomplete);
1134  DeclException0 (ExcInvalidComponentOrder);
1142  DeclException0 (ExcNotDGFEM);
1143 }
1144 
1145 /* ------------------------- inline functions -------------- */
1146 
1147 #ifndef DOXYGEN
1148 namespace DoFRenumbering
1149 {
1150  template <class DH>
1151  void
1152  inline
1153  downstream (DH &dof,
1154  const Point<DH::space_dimension> &direction,
1155  const bool dof_wise_renumbering)
1156  {
1157  std::vector<types::global_dof_index> renumbering(dof.n_dofs());
1158  std::vector<types::global_dof_index> reverse(dof.n_dofs());
1159  compute_downstream(renumbering, reverse, dof, direction,
1160  dof_wise_renumbering);
1161 
1162  dof.renumber_dofs(renumbering);
1163  }
1164 }
1165 #endif
1166 
1167 
1168 DEAL_II_NAMESPACE_CLOSE
1169 
1170 #endif
types::global_dof_index compute_block_wise(std::vector< types::global_dof_index > &new_dof_indices, const ITERATOR &start, const ENDITERATOR &end)
void sort_selected_dofs_back(DH &dof_handler, const std::vector< bool > &selected_dofs)
void random(DH &dof_handler)
void compute_subdomain_wise(std::vector< types::global_dof_index > &new_dof_indices, const DH &dof_handler)
void Cuthill_McKee(DH &dof_handler, const bool reversed_numbering=false, const bool use_constraints=false)
void compute_random(std::vector< types::global_dof_index > &new_dof_indices, const DH &dof_handler)
void Cuthill_McKee(DH &dof_handler, const bool reversed_numbering=false, const bool use_constraints=false, const std::vector< types::global_dof_index > &starting_indices=std::vector< types::global_dof_index >())
void hierarchical(DoFHandler< dim > &dof_handler)
void component_wise(DoFHandler< dim, spacedim > &dof_handler, const std::vector< unsigned int > &target_component=std::vector< unsigned int >())
void compute_clockwise_dg(std::vector< types::global_dof_index > &new_dof_indices, const DH &dof_handler, const Point< DH::space_dimension > &center, const bool counter)
void compute_sort_selected_dofs_back(std::vector< types::global_dof_index > &new_dof_indices, const DH &dof_handler, const std::vector< bool > &selected_dofs)
void compute_cell_wise(std::vector< types::global_dof_index > &renumbering, std::vector< types::global_dof_index > &inverse_renumbering, const DH &dof_handler, const std::vector< typename DH::active_cell_iterator > &cell_order)
void downstream(DH &dof_handler, const Point< DH::space_dimension > &direction, const bool dof_wise_renumbering=false)
void compute_minimum_degree(std::vector< types::global_dof_index > &new_dof_indices, const DH &, const bool reversed_numbering=false, const bool use_constraints=false)
void compute_downstream(std::vector< types::global_dof_index > &new_dof_indices, std::vector< types::global_dof_index > &reverse, const DH &dof_handler, const Point< DH::space_dimension > &direction, const bool dof_wise_renumbering)
void minimum_degree(DH &dof_handler, const bool reversed_numbering=false, const bool use_constraints=false)
DeclException0(ExcRenumberingIncomplete)
unsigned int global_dof_index
Definition: types.h:100
void clockwise_dg(DH &dof_handler, const Point< DH::space_dimension > &center, const bool counter=false)
void compute_king_ordering(std::vector< types::global_dof_index > &new_dof_indices, const DH &, const bool reversed_numbering=false, const bool use_constraints=false)
void compute_Cuthill_McKee(std::vector< types::global_dof_index > &new_dof_indices, const DH &, const bool reversed_numbering=false, const bool use_constraints=false)
void downstream_dg(DH &dof, const Point< DH::space_dimension > &direction) DEAL_II_DEPRECATED
ComparePointwiseDownstream(const Point< dim > &dir)
BlockCompressedSparsityPattern CompressedBlockSparsityPattern DEAL_II_DEPRECATED
void compute_Cuthill_McKee(std::vector< types::global_dof_index > &new_dof_indices, const DH &, const bool reversed_numbering=false, const bool use_constraints=false, const std::vector< types::global_dof_index > &starting_indices=std::vector< types::global_dof_index >())
void king_ordering(DH &dof_handler, const bool reversed_numbering=false, const bool use_constraints=false)
bool operator()(const Iterator &c1, const Iterator &c2) const
void block_wise(DoFHandler< dim, spacedim > &dof_handler)
types::global_dof_index compute_component_wise(std::vector< types::global_dof_index > &new_dof_indices, const ITERATOR &start, const ENDITERATOR &end, const std::vector< unsigned int > &target_component, bool is_level_operation)
CompareDownstream(const Point< dim > &dir)
bool operator()(const std::pair< Point< dim >, types::global_dof_index > &c1, const std::pair< Point< dim >, types::global_dof_index > &c2) const
void cell_wise(DH &dof_handler, const std::vector< typename DH::active_cell_iterator > &cell_order)
void subdomain_wise(DH &dof_handler)