Reference documentation for deal.II version 8.4.1
block_sparsity_pattern.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2000 - 2015 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 at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii__block_sparsity_pattern_h
17 #define dealii__block_sparsity_pattern_h
18 
19 
20 #include <deal.II/base/config.h>
21 #include <deal.II/base/exceptions.h>
22 #include <deal.II/base/table.h>
23 #include <deal.II/base/subscriptor.h>
24 #include <deal.II/base/smartpointer.h>
25 #include <deal.II/lac/sparsity_pattern.h>
26 #include <deal.II/lac/trilinos_sparsity_pattern.h>
27 #include <deal.II/lac/compressed_sparsity_pattern.h>
28 #include <deal.II/lac/compressed_set_sparsity_pattern.h>
29 #include <deal.II/lac/compressed_simple_sparsity_pattern.h>
30 #include <deal.II/lac/block_indices.h>
31 
32 DEAL_II_NAMESPACE_OPEN
33 
34 
35 template <typename number> class BlockSparseMatrix;
38 #ifdef DEAL_II_WITH_TRILINOS
39 namespace TrilinosWrappers
40 {
42 }
43 #endif
44 
82 template <typename SparsityPatternType>
84 {
85 public:
90 
100 
107 
113  BlockSparsityPatternBase (const size_type n_block_rows,
114  const size_type n_block_columns);
115 
123 
128 
142  void reinit (const size_type n_block_rows,
143  const size_type n_block_columns);
144 
151 
158  void collect_sizes ();
159 
163  SparsityPatternType &
164  block (const size_type row,
165  const size_type column);
166 
167 
172  const SparsityPatternType &
173  block (const size_type row,
174  const size_type column) const;
175 
180  const BlockIndices &
181  get_row_indices () const;
182 
187  const BlockIndices &
188  get_column_indices () const;
189 
194  void compress ();
195 
199  size_type n_block_rows () const;
200 
204  size_type n_block_cols () const;
205 
212  bool empty () const;
213 
219  size_type max_entries_per_row () const;
220 
230  void add (const size_type i, const size_type j);
231 
242  template <typename ForwardIterator>
243  void add_entries (const size_type row,
244  ForwardIterator begin,
245  ForwardIterator end,
246  const bool indices_are_sorted = false);
247 
252  size_type n_rows () const;
253 
259  size_type n_cols () const;
260 
264  bool exists (const size_type i, const size_type j) const;
265 
270  unsigned int row_length (const size_type row) const;
271 
283  size_type n_nonzero_elements () const;
284 
290  void print (std::ostream &out) const;
291 
299  void print_gnuplot (std::ostream &out) const;
300 
309  DeclException4 (ExcIncompatibleRowNumbers,
310  int, int, int, int,
311  << "The blocks [" << arg1 << ',' << arg2 << "] and ["
312  << arg3 << ',' << arg4 << "] have differing row numbers.");
316  DeclException4 (ExcIncompatibleColNumbers,
317  int, int, int, int,
318  << "The blocks [" << arg1 << ',' << arg2 << "] and ["
319  << arg3 << ',' << arg4 << "] have differing column numbers.");
323  DeclException0 (ExcInvalidConstructorCall);
325 
326 protected:
327 
331  size_type rows;
332 
336  size_type columns;
337 
342 
348 
354 
355 private:
360  std::vector<size_type > counter_within_block;
361 
366  std::vector<std::vector<size_type > > block_column_indices;
367 
372  template <typename number>
373  friend class BlockSparseMatrix;
374 };
375 
376 
377 
389 class BlockSparsityPattern : public BlockSparsityPatternBase<SparsityPattern>
390 {
391 public:
392 
399 
406  const size_type n_columns);
407 
411  void reinit (const size_type n_block_rows,
412  const size_type n_block_columns);
413 
426  void reinit (const BlockIndices &row_indices,
427  const BlockIndices &col_indices,
428  const std::vector<std::vector<unsigned int> > &row_lengths);
429 
430 
435  bool is_compressed () const;
436 
441  std::size_t memory_consumption () const;
442 
448  void copy_from (const BlockDynamicSparsityPattern &dsp);
449 };
450 
451 
452 
495 class BlockDynamicSparsityPattern : public BlockSparsityPatternBase<DynamicSparsityPattern>
496 {
497 public:
498 
505 
512  const size_type n_columns);
513 
520  BlockDynamicSparsityPattern (const std::vector<size_type> &row_block_sizes,
521  const std::vector<size_type> &col_block_sizes);
522 
531  BlockDynamicSparsityPattern (const std::vector<IndexSet> &partitioning);
532 
538  const BlockIndices &col_indices);
539 
540 
549  void reinit (const std::vector<size_type> &row_block_sizes,
550  const std::vector<size_type> &col_block_sizes);
551 
556  void reinit(const std::vector<IndexSet> &partitioning);
557 
563  void reinit (const BlockIndices &row_indices, const BlockIndices &col_indices);
564 
569  size_type column_number (const size_type row,
570  const unsigned int index) const;
571 
576 };
577 
581 typedef BlockDynamicSparsityPattern BlockCompressedSparsityPattern DEAL_II_DEPRECATED;
582 
586 typedef BlockDynamicSparsityPattern BlockCompressedSetSparsityPattern DEAL_II_DEPRECATED;
587 
591 typedef BlockDynamicSparsityPattern BlockCompressedSimpleSparsityPattern DEAL_II_DEPRECATED;
592 
593 
594 
595 
596 #ifdef DEAL_II_WITH_TRILINOS
597 
598 
618 namespace TrilinosWrappers
619 {
620  class BlockSparsityPattern :
621  public ::BlockSparsityPatternBase<SparsityPattern>
622  {
623  public:
624 
631 
637  BlockSparsityPattern (const size_type n_rows,
638  const size_type n_columns);
639 
646  BlockSparsityPattern (const std::vector<size_type> &row_block_sizes,
647  const std::vector<size_type> &col_block_sizes);
648 
659  BlockSparsityPattern (const std::vector<Epetra_Map> &parallel_partitioning) DEAL_II_DEPRECATED;
660 
668  BlockSparsityPattern (const std::vector<IndexSet> &parallel_partitioning,
669  const MPI_Comm &communicator = MPI_COMM_WORLD);
670 
682  BlockSparsityPattern (const std::vector<IndexSet> &row_parallel_partitioning,
683  const std::vector<IndexSet> &column_parallel_partitioning,
684  const std::vector<IndexSet> &writeable_rows,
685  const MPI_Comm &communicator = MPI_COMM_WORLD);
686 
696  void reinit (const std::vector<size_type> &row_block_sizes,
697  const std::vector<size_type> &col_block_sizes);
698 
706  void reinit (const std::vector<Epetra_Map> &parallel_partitioning) DEAL_II_DEPRECATED;
707 
712  void reinit (const std::vector<IndexSet> &parallel_partitioning,
713  const MPI_Comm &communicator = MPI_COMM_WORLD);
714 
720  void reinit (const std::vector<IndexSet> &row_parallel_partitioning,
721  const std::vector<IndexSet> &column_parallel_partitioning,
722  const MPI_Comm &communicator = MPI_COMM_WORLD);
723 
732  void reinit (const std::vector<IndexSet> &row_parallel_partitioning,
733  const std::vector<IndexSet> &column_parallel_partitioning,
734  const std::vector<IndexSet> &writeable_rows,
735  const MPI_Comm &communicator = MPI_COMM_WORLD);
736 
741  };
742 }
743 
744 #endif
745 
746 
748 /*---------------------- Template functions -----------------------------------*/
749 
750 
751 
752 template <typename SparsityPatternType>
753 inline
754 SparsityPatternType &
756  const size_type column)
757 {
758  Assert (row<rows, ExcIndexRange(row,0,rows));
759  Assert (column<columns, ExcIndexRange(column,0,columns));
760  return *sub_objects[row][column];
761 }
762 
763 
764 
765 template <typename SparsityPatternType>
766 inline
767 const SparsityPatternType &
769  const size_type column) const
770 {
771  Assert (row<rows, ExcIndexRange(row,0,rows));
772  Assert (column<columns, ExcIndexRange(column,0,columns));
773  return *sub_objects[row][column];
774 }
775 
776 
777 
778 template <typename SparsityPatternType>
779 inline
780 const BlockIndices &
782 {
783  return row_indices;
784 }
785 
786 
787 
788 template <typename SparsityPatternType>
789 inline
790 const BlockIndices &
792 {
793  return column_indices;
794 }
795 
796 
797 
798 template <typename SparsityPatternType>
799 inline
800 void
802  const size_type j)
803 {
804  // if you get an error here, are
805  // you sure you called
806  // <tt>collect_sizes()</tt> before?
807  const std::pair<size_type,size_type>
808  row_index = row_indices.global_to_local (i),
809  col_index = column_indices.global_to_local (j);
810  sub_objects[row_index.first][col_index.first]->add (row_index.second,
811  col_index.second);
812 }
813 
814 
815 
816 template <typename SparsityPatternType>
817 template <typename ForwardIterator>
818 void
820  ForwardIterator begin,
821  ForwardIterator end,
822  const bool indices_are_sorted)
823 {
824  // Resize scratch arrays
825  if (block_column_indices.size() < this->n_block_cols())
826  {
827  block_column_indices.resize (this->n_block_cols());
828  counter_within_block.resize (this->n_block_cols());
829  }
830 
831  const size_type n_cols = static_cast<size_type>(end - begin);
832 
833  // Resize sub-arrays to n_cols. This
834  // is a bit wasteful, but we resize
835  // only a few times (then the maximum
836  // row length won't increase that
837  // much any more). At least we know
838  // that all arrays are going to be of
839  // the same size, so we can check
840  // whether the size of one is large
841  // enough before actually going
842  // through all of them.
843  if (block_column_indices[0].size() < n_cols)
844  for (size_type i=0; i<this->n_block_cols(); ++i)
845  block_column_indices[i].resize(n_cols);
846 
847  // Reset the number of added elements
848  // in each block to zero.
849  for (size_type i=0; i<this->n_block_cols(); ++i)
850  counter_within_block[i] = 0;
851 
852  // Go through the column indices to
853  // find out which portions of the
854  // values should be set in which
855  // block of the matrix. We need to
856  // touch all the data, since we can't
857  // be sure that the data of one block
858  // is stored contiguously (in fact,
859  // indices will be intermixed when it
860  // comes from an element matrix).
861  for (ForwardIterator it = begin; it != end; ++it)
862  {
863  const size_type col = *it;
864 
865  const std::pair<size_type , size_type>
866  col_index = this->column_indices.global_to_local(col);
867 
868  const size_type local_index = counter_within_block[col_index.first]++;
869 
870  block_column_indices[col_index.first][local_index] = col_index.second;
871  }
872 
873 #ifdef DEBUG
874  // If in debug mode, do a check whether
875  // the right length has been obtained.
876  size_type length = 0;
877  for (size_type i=0; i<this->n_block_cols(); ++i)
878  length += counter_within_block[i];
879  Assert (length == n_cols, ExcInternalError());
880 #endif
881 
882  // Now we found out about where the
883  // individual columns should start and
884  // where we should start reading out
885  // data. Now let's write the data into
886  // the individual blocks!
887  const std::pair<size_type , size_type>
888  row_index = this->row_indices.global_to_local (row);
889  for (size_type block_col=0; block_col<n_block_cols(); ++block_col)
890  {
891  if (counter_within_block[block_col] == 0)
892  continue;
893  sub_objects[row_index.first][block_col]->
894  add_entries (row_index.second,
895  block_column_indices[block_col].begin(),
896  block_column_indices[block_col].begin()+counter_within_block[block_col],
897  indices_are_sorted);
898  }
899 }
900 
901 
902 
903 template <typename SparsityPatternType>
904 inline
905 bool
907  const size_type j) const
908 {
909  // if you get an error here, are
910  // you sure you called
911  // <tt>collect_sizes()</tt> before?
912  const std::pair<size_type , size_type>
913  row_index = row_indices.global_to_local (i),
914  col_index = column_indices.global_to_local (j);
915  return sub_objects[row_index.first][col_index.first]->exists (row_index.second,
916  col_index.second);
917 }
918 
919 
920 
921 template <typename SparsityPatternType>
922 inline
923 unsigned int
925 row_length (const size_type row) const
926 {
927  const std::pair<size_type , size_type>
928  row_index = row_indices.global_to_local (row);
929 
930  unsigned int c = 0;
931 
932  for (size_type b=0; b<rows; ++b)
933  c += sub_objects[row_index.first][b]->row_length (row_index.second);
934 
935  return c;
936 }
937 
938 
939 
940 template <typename SparsityPatternType>
941 inline
944 {
945  return columns;
946 }
947 
948 
949 
950 template <typename SparsityPatternType>
951 inline
954 {
955  return rows;
956 }
957 
958 
959 inline
962  const unsigned int index) const
963 {
964  // .first= ith block, .second = jth row in that block
965  const std::pair<size_type ,size_type >
966  row_index = row_indices.global_to_local (row);
967 
968  Assert(index<row_length(row), ExcIndexRange(index, 0, row_length(row)));
969 
970  size_type c = 0;
971  size_type block_columns = 0; //sum of n_cols for all blocks to the left
972  for (unsigned int b=0; b<columns; ++b)
973  {
974  unsigned int rowlen = sub_objects[row_index.first][b]->row_length (row_index.second);
975  if (index<c+rowlen)
976  return block_columns+sub_objects[row_index.first][b]->column_number(row_index.second, index-c);
977  c += rowlen;
978  block_columns += sub_objects[row_index.first][b]->n_cols();
979  }
980 
981  Assert(false, ExcInternalError());
982  return 0;
983 }
984 
985 
986 inline
987 void
989  const size_type n_block_rows,
990  const size_type n_block_columns)
991 {
993  n_block_rows, n_block_columns);
994 }
995 
996 
997 DEAL_II_NAMESPACE_CLOSE
998 
999 #endif
const BlockIndices & get_row_indices() const
void reinit(const std::vector< size_type > &row_block_sizes, const std::vector< size_type > &col_block_sizes)
void print(std::ostream &out) const
static const size_type invalid_entry
Auxiliary class aiding in the handling of block structures like in BlockVector or FESystem...
Definition: block_indices.h:54
std::size_t memory_consumption() const
size_type n_nonzero_elements() const
SparsityPatternType & block(const size_type row, const size_type column)
std::vector< size_type > counter_within_block
void reinit(const size_type n_block_rows, const size_type n_block_columns)
types::global_dof_index size_type
BlockDynamicSparsityPattern BlockCompressedSparsityPattern DEAL_II_DEPRECATED
bool exists(const size_type i, const size_type j) const
void reinit(const size_type n_block_rows, const size_type n_block_columns)
unsigned int row_length(const size_type row) const
static const size_type invalid_entry
void add(const size_type i, const size_type j)
unsigned int global_dof_index
Definition: types.h:88
#define Assert(cond, exc)
Definition: exceptions.h:294
void add_entries(const size_type row, ForwardIterator begin, ForwardIterator end, const bool indices_are_sorted=false)
DeclException4(ExcIncompatibleRowNumbers, int, int, int, int,<< "The blocks ["<< arg1<< ','<< arg2<< "] and ["<< arg3<< ','<< arg4<< "] have differing row numbers.")
BlockSparsityPatternBase & operator=(const BlockSparsityPatternBase &)
size_type max_entries_per_row() const
std::pair< unsigned int, size_type > global_to_local(const size_type i) const
std::vector< std::vector< size_type > > block_column_indices
const BlockIndices & get_column_indices() const
Table< 2, SmartPointer< SparsityPatternType, BlockSparsityPatternBase< SparsityPatternType > > > sub_objects
void print_gnuplot(std::ostream &out) const
size_type column_number(const size_type row, const unsigned int index) const
Definition: table.h:33
void copy_from(const BlockDynamicSparsityPattern &dsp)
DeclException0(ExcInvalidConstructorCall)