namespace MPI
{
class SparseMatrix;
+ class BlockSparseMatrix;
class Vector;
+ class BlockVector;
}
}
#endif
PETScWrappers::Vector &solution,
PETScWrappers::MPI::Vector &right_hand_side,
const bool eliminate_columns = true);
+
+ /**
+ * Same as above but for BlockSparseMatrix.
+ */
+ static void
+ apply_boundary_values (const std::map<unsigned int,double> &boundary_values,
+ PETScWrappers::MPI::BlockSparseMatrix &matrix,
+ PETScWrappers::MPI::BlockVector &solution,
+ PETScWrappers::MPI::BlockVector &right_hand_side,
+ const bool eliminate_columns = true);
+
#endif
#ifdef DEAL_II_USE_TRILINOS
# include <lac/petsc_sparse_matrix.h>
# include <lac/petsc_parallel_vector.h>
# include <lac/petsc_vector.h>
+# include <lac/petsc_parallel_block_sparse_matrix.h>
#endif
#ifdef DEAL_II_USE_TRILINOS
}
+void
+MatrixTools::
+apply_boundary_values (const std::map<unsigned int,double> &boundary_values,
+ PETScWrappers::MPI::BlockSparseMatrix &matrix,
+ PETScWrappers::MPI::BlockVector &solution,
+ PETScWrappers::MPI::BlockVector &right_hand_side,
+ const bool eliminate_columns)
+{
+ Assert (matrix.n() == right_hand_side.size(),
+ ExcDimensionMismatch(matrix.n(), right_hand_side.size()));
+ Assert (matrix.n() == solution.size(),
+ ExcDimensionMismatch(matrix.n(), solution.size()));
+ Assert (matrix.n_block_rows() == matrix.n_block_cols(),
+ ExcNotQuadratic());
+
+ const unsigned int n_blocks = matrix.n_block_rows();
+
+ matrix.compress();
+
+ // We need to find the subdivision
+ // into blocks for the boundary values.
+ // To this end, generate a vector of
+ // maps with the respective indices.
+ std::vector<std::map<unsigned int,double> > block_boundary_values(n_blocks);
+ {
+ int offset = 0, block=0;
+ for (std::map<unsigned int,double>::const_iterator
+ dof = boundary_values.begin();
+ dof != boundary_values.end();
+ ++dof)
+ {
+ if (dof->first >= matrix.block(block,0).m() + offset)
+ {
+ offset += matrix.block(block,0).m();
+ block++;
+ }
+ const unsigned int index = dof->first - offset;
+ block_boundary_values[block].insert(std::pair<unsigned int, double> (index,dof->second));
+ }
+ }
+
+ // Now call the non-block variants on
+ // the diagonal subblocks and the
+ // solution/rhs.
+ for (unsigned int block=0; block<n_blocks; ++block)
+ PETScWrappers::apply_boundary_values(block_boundary_values[block],
+ matrix.block(block,block),
+ solution.block(block),
+ right_hand_side.block(block),
+ eliminate_columns);
+
+ // Finally, we need to do something
+ // about the off-diagonal matrices. This
+ // is luckily not difficult. Just clear
+ // the whole row.
+ for (unsigned int block_m=0; block_m<n_blocks; ++block_m)
+ {
+ const std::pair<unsigned int, unsigned int> local_range
+ = matrix.block(block_m,0).local_range();
+
+ std::vector<unsigned int> constrained_rows;
+ for (std::map<unsigned int,double>::const_iterator
+ dof = block_boundary_values[block_m].begin();
+ dof != block_boundary_values[block_m].end();
+ ++dof)
+ if ((dof->first >= local_range.first) &&
+ (dof->first < local_range.second))
+ constrained_rows.push_back (dof->first);
+
+ for (unsigned int block_n=0; block_n<n_blocks; ++block_n)
+ if (block_m != block_n)
+ matrix.block(block_m,block_n).clear_rows(constrained_rows);
+ }
+}
+
#endif
(WB 2008/09/23)
</p>
+ <li>
+ <p>
+ New: Added the CompressedSimpleSparsityPattern as an alternative to
+ CompressedSparsityPattern and CompressedSetSparsityPattern, which
+ should be faster in most cases but its memory usage is somewhere
+ inbetween.
+ <br>
+ (Timo Heister 2008/09/03)
+ </p>
+
<li>
<p>
Improved: The CompressedSparsityPattern can now elide some operations
<h3>deal.II</h3>
<ol>
+ <li>
+ <p>
+ New: MatrixTools::apply_boundary_values() also for PETScWrappers::MPI::BlockSparseMatrix.
+ <br>
+ (Timo Heister 2008/10/27)
+ </p>
+
<li>
<p>
New: When calling function DoFTools::make_sparsity_pattern with a ConstraintMatrix, it is now possible to set a bool argument keep_constrained_dofs. When this flag is set to false, constrained rows and columns will not be part of the sparsity pattern, which increases the performance of matrix operations and decrease memory consumption in case there are many constraints.