1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2018 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/lac/petsc_block_sparse_matrix.h>
17 
18 #ifdef DEAL_II_WITH_PETSC
19 
20 DEAL_II_NAMESPACE_OPEN
21 
22 namespace PETScWrappers
23 {
24   namespace MPI
25   {
26     BlockSparseMatrix &
operator =(const BlockSparseMatrix & m)27     BlockSparseMatrix::operator=(const BlockSparseMatrix &m)
28     {
29       BaseClass::operator=(m);
30 
31       return *this;
32     }
33 
34 
35     void
reinit(const size_type n_block_rows,const size_type n_block_columns)36     BlockSparseMatrix::reinit(const size_type n_block_rows,
37                               const size_type n_block_columns)
38     {
39       // first delete previous content of
40       // the subobjects array
41       clear();
42 
43       // then resize. set sizes of blocks to
44       // zero. user will later have to call
45       // collect_sizes for this
46       this->sub_objects.reinit(n_block_rows, n_block_columns);
47       this->row_block_indices.reinit(n_block_rows, 0);
48       this->column_block_indices.reinit(n_block_columns, 0);
49 
50       // and reinitialize the blocks
51       for (size_type r = 0; r < this->n_block_rows(); ++r)
52         for (size_type c = 0; c < this->n_block_cols(); ++c)
53           {
54             BlockType *p            = new BlockType();
55             this->sub_objects[r][c] = p;
56           }
57     }
58 
59     void
reinit(const std::vector<IndexSet> & rows,const std::vector<IndexSet> & cols,const BlockDynamicSparsityPattern & bdsp,const MPI_Comm & com)60     BlockSparseMatrix::reinit(const std::vector<IndexSet> &      rows,
61                               const std::vector<IndexSet> &      cols,
62                               const BlockDynamicSparsityPattern &bdsp,
63                               const MPI_Comm &                   com)
64     {
65       Assert(rows.size() == bdsp.n_block_rows(), ExcMessage("invalid size"));
66       Assert(cols.size() == bdsp.n_block_cols(), ExcMessage("invalid size"));
67 
68 
69       clear();
70       this->sub_objects.reinit(bdsp.n_block_rows(), bdsp.n_block_cols());
71 
72       std::vector<types::global_dof_index> row_sizes;
73       for (unsigned int r = 0; r < bdsp.n_block_rows(); ++r)
74         row_sizes.push_back(bdsp.block(r, 0).n_rows());
75       this->row_block_indices.reinit(row_sizes);
76 
77       std::vector<types::global_dof_index> col_sizes;
78       for (unsigned int c = 0; c < bdsp.n_block_cols(); ++c)
79         col_sizes.push_back(bdsp.block(0, c).n_cols());
80       this->column_block_indices.reinit(col_sizes);
81 
82       for (unsigned int r = 0; r < this->n_block_rows(); ++r)
83         for (unsigned int c = 0; c < this->n_block_cols(); ++c)
84           {
85             Assert(rows[r].size() == bdsp.block(r, c).n_rows(),
86                    ExcMessage("invalid size"));
87             Assert(cols[c].size() == bdsp.block(r, c).n_cols(),
88                    ExcMessage("invalid size"));
89 
90             BlockType *p = new BlockType();
91             p->reinit(rows[r], cols[c], bdsp.block(r, c), com);
92             this->sub_objects[r][c] = p;
93           }
94 
95       collect_sizes();
96     }
97 
98     void
reinit(const std::vector<IndexSet> & sizes,const BlockDynamicSparsityPattern & bdsp,const MPI_Comm & com)99     BlockSparseMatrix::reinit(const std::vector<IndexSet> &      sizes,
100                               const BlockDynamicSparsityPattern &bdsp,
101                               const MPI_Comm &                   com)
102     {
103       reinit(sizes, sizes, bdsp, com);
104     }
105 
106 
107 
108     void
collect_sizes()109     BlockSparseMatrix::collect_sizes()
110     {
111       BaseClass::collect_sizes();
112     }
113 
114     std::vector<IndexSet>
locally_owned_domain_indices() const115     BlockSparseMatrix::locally_owned_domain_indices() const
116     {
117       std::vector<IndexSet> index_sets;
118 
119       for (unsigned int i = 0; i < this->n_block_cols(); ++i)
120         index_sets.push_back(this->block(0, i).locally_owned_domain_indices());
121 
122       return index_sets;
123     }
124 
125     std::vector<IndexSet>
locally_owned_range_indices() const126     BlockSparseMatrix::locally_owned_range_indices() const
127     {
128       std::vector<IndexSet> index_sets;
129 
130       for (unsigned int i = 0; i < this->n_block_rows(); ++i)
131         index_sets.push_back(this->block(i, 0).locally_owned_range_indices());
132 
133       return index_sets;
134     }
135 
136     const MPI_Comm &
get_mpi_communicator() const137     BlockSparseMatrix::get_mpi_communicator() const
138     {
139       return block(0, 0).get_mpi_communicator();
140     }
141 
142   } // namespace MPI
143 } // namespace PETScWrappers
144 
145 
146 
147 DEAL_II_NAMESPACE_CLOSE
148 
149 #endif
150