Reference documentation for deal.II version Git 92ca227 2018-05-25 19:47:08 -0500
constraint_matrix.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2017 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 
17 #ifndef dealii_constraint_matrix_h
18 #define dealii_constraint_matrix_h
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/exceptions.h>
23 #include <deal.II/base/index_set.h>
24 #include <deal.II/base/subscriptor.h>
25 #include <deal.II/base/template_constraints.h>
26 
27 #include <deal.II/lac/vector.h>
28 #include <deal.II/lac/vector_element_access.h>
29 
30 #include <boost/range/iterator_range.hpp>
31 
32 #include <set>
33 #include <utility>
34 #include <vector>
35 
36 
37 DEAL_II_NAMESPACE_OPEN
38 
39 template <int dim, class T>
40 class Table;
41 template <typename>
42 class FullMatrix;
43 class SparsityPattern;
47 template <typename number>
48 class SparseMatrix;
49 template <typename number>
50 class BlockSparseMatrix;
51 
52 namespace internals
53 {
54  class GlobalRowsFromLocal;
55 }
56 
57 
58 // TODO[WB]: We should have a function of the kind
59 // ConstraintMatrix::add_constraint (const size_type constrained_dof,
60 // const std::vector<std::pair<size_type, double> > &entries,
61 // const double inhomogeneity = 0);
62 // rather than building up constraints piecemeal through add_line/add_entry
63 // etc. This would also eliminate the possibility of accidentally changing
64 // existing constraints into something pointless, see the discussion on the
65 // mailing list on "Tiny bug in interpolate_boundary_values" in Sept. 2010.
66 
143 {
144 public:
149 
156  {
162 
169 
176  };
177 
195  explicit ConstraintMatrix(const IndexSet &local_constraints = IndexSet());
196 
200  explicit ConstraintMatrix(const ConstraintMatrix &constraint_matrix);
201 
205  ConstraintMatrix(ConstraintMatrix &&constraint_matrix) = default;
206 
217  operator=(const ConstraintMatrix &) = delete;
218 
223  operator=(ConstraintMatrix &&constraint_matrix) = default;
224 
231  void
232  copy_from(const ConstraintMatrix &other);
233 
240  void
241  reinit(const IndexSet &local_constraints = IndexSet());
242 
249  bool
250  can_store_line(const size_type line_index) const;
251 
258  const IndexSet &
259  get_local_lines() const;
260 
280  void
281  add_selected_constraints(const ConstraintMatrix &constraints_in,
282  const IndexSet & filter);
283 
293  void
294  add_line(const size_type line);
295 
308  void
309  add_lines(const std::vector<bool> &lines);
310 
323  void
324  add_lines(const std::set<size_type> &lines);
325 
338  void
339  add_lines(const IndexSet &lines);
340 
352  void
353  add_entry(const size_type line, const size_type column, const double value);
354 
360  void
361  add_entries(const size_type line,
362  const std::vector<std::pair<size_type, double>> &col_val_pairs);
363 
370  void
371  set_inhomogeneity(const size_type line, const double value);
372 
393  void
394  close();
395 
420  void
421  merge(
422  const ConstraintMatrix & other_constraints,
423  const MergeConflictBehavior merge_conflict_behavior = no_conflicts_allowed,
424  const bool allow_different_local_lines = false);
425 
438  void
439  shift(const size_type offset);
440 
448  void
449  clear();
450 
464  size_type
465  n_constraints() const;
466 
476  bool
477  is_constrained(const size_type index) const;
478 
490  bool
491  is_identity_constrained(const size_type index) const;
492 
499  bool
500  are_identity_constrained(const size_type index1,
501  const size_type index2) const;
502 
513  size_type
515 
520  bool
521  is_inhomogeneously_constrained(const size_type index) const;
522 
528  bool
529  has_inhomogeneities() const;
530 
535  const std::vector<std::pair<size_type, double>> *
536  get_constraint_entries(const size_type line) const;
537 
542  double
543  get_inhomogeneity(const size_type line) const;
544 
565  void
566  print(std::ostream &out) const;
567 
580  void
581  write_dot(std::ostream &) const;
582 
587  std::size_t
588  memory_consumption() const;
589 
596  void
597  resolve_indices(std::vector<types::global_dof_index> &indices) const;
598 
626  void
627  condense(SparsityPattern &sparsity) const;
628 
632  void
633  condense(BlockSparsityPattern &sparsity) const;
634 
639  void
640  condense(DynamicSparsityPattern &sparsity) const;
641 
646  void
647  condense(BlockDynamicSparsityPattern &sparsity) const;
648 
656  template <typename number>
657  void
658  condense(SparseMatrix<number> &matrix) const;
659 
663  template <typename number>
664  void
665  condense(BlockSparseMatrix<number> &matrix) const;
666 
678  template <class VectorType>
679  void
680  condense(VectorType &vec) const;
681 
688  template <class VectorType>
689  void
690  condense(const VectorType &vec_ghosted, VectorType &output) const;
691 
704  template <typename number, class VectorType>
705  void
706  condense(SparseMatrix<number> &matrix, VectorType &vector) const;
707 
712  template <typename number, class BlockVectorType>
713  void
714  condense(BlockSparseMatrix<number> &matrix, BlockVectorType &vector) const;
715 
722  template <class VectorType>
723  void
724  set_zero(VectorType &vec) const;
725 
777  template <class InVector, class OutVector>
778  void
779  distribute_local_to_global(const InVector & local_vector,
780  const std::vector<size_type> &local_dof_indices,
781  OutVector & global_vector) const;
782 
826  template <typename VectorType, typename LocalType>
827  void
828  distribute_local_to_global(const Vector<LocalType> & local_vector,
829  const std::vector<size_type> &local_dof_indices,
830  VectorType & global_vector,
831  const FullMatrix<LocalType> & local_matrix) const;
832 
852  template <typename VectorType, typename LocalType>
853  void
855  const Vector<LocalType> & local_vector,
856  const std::vector<size_type> &local_dof_indices_row,
857  const std::vector<size_type> &local_dof_indices_col,
858  VectorType & global_vector,
859  const FullMatrix<LocalType> & local_matrix,
860  bool diagonal = false) const;
861 
865  template <class VectorType>
866  void
867  distribute_local_to_global(const size_type index,
868  const double value,
869  VectorType & global_vector) const;
870 
899  template <typename ForwardIteratorVec,
900  typename ForwardIteratorInd,
901  class VectorType>
902  void
903  distribute_local_to_global(ForwardIteratorVec local_vector_begin,
904  ForwardIteratorVec local_vector_end,
905  ForwardIteratorInd local_indices_begin,
906  VectorType & global_vector) const;
907 
955  template <typename MatrixType>
956  void
959  const std::vector<size_type> & local_dof_indices,
960  MatrixType & global_matrix) const;
961 
989  template <typename MatrixType>
990  void
993  const std::vector<size_type> & row_indices,
994  const std::vector<size_type> & col_indices,
995  MatrixType & global_matrix) const;
996 
1013  template <typename MatrixType>
1014  void
1016  const FullMatrix<typename MatrixType::value_type> &local_matrix,
1017  const std::vector<size_type> & row_indices,
1018  const ConstraintMatrix & column_constraint_matrix,
1019  const std::vector<size_type> & column_indices,
1020  MatrixType & global_matrix) const;
1021 
1038  template <typename MatrixType, typename VectorType>
1039  void
1041  const FullMatrix<typename MatrixType::value_type> &local_matrix,
1042  const Vector<typename VectorType::value_type> & local_vector,
1043  const std::vector<size_type> & local_dof_indices,
1044  MatrixType & global_matrix,
1045  VectorType & global_vector,
1046  bool use_inhomogeneities_for_rhs = false) const;
1047 
1101  template <typename SparsityPatternType>
1102  void
1104  const std::vector<size_type> &local_dof_indices,
1105  SparsityPatternType & sparsity_pattern,
1106  const bool keep_constrained_entries = true,
1107  const Table<2, bool> & dof_mask = default_empty_table) const;
1108 
1112  template <typename SparsityPatternType>
1113  void
1115  const std::vector<size_type> &row_indices,
1116  const std::vector<size_type> &col_indices,
1117  SparsityPatternType & sparsity_pattern,
1118  const bool keep_constrained_entries = true,
1119  const Table<2, bool> & dof_mask = default_empty_table) const;
1120 
1140  template <typename ForwardIteratorVec,
1141  typename ForwardIteratorInd,
1142  class VectorType>
1143  void
1144  get_dof_values(const VectorType & global_vector,
1145  ForwardIteratorInd local_indices_begin,
1146  ForwardIteratorVec local_vector_begin,
1147  ForwardIteratorVec local_vector_end) const;
1148 
1170  template <class VectorType>
1171  void
1172  distribute(VectorType &vec) const;
1173 
1184  {
1189  typedef std::vector<std::pair<size_type, double>> Entries;
1190 
1196  size_type index;
1197 
1205  Entries entries;
1206 
1211 
1219  bool
1220  operator<(const ConstraintLine &) const;
1221 
1227  bool
1228  operator==(const ConstraintLine &) const;
1229 
1234  std::size_t
1235  memory_consumption() const;
1236 
1240  template <class Archive>
1241  void
1242  serialize(Archive &ar, const unsigned int)
1243  {
1244  ar &index &entries &inhomogeneity;
1245  }
1246  };
1247 
1248 
1252  typedef std::vector<ConstraintLine>::const_iterator const_iterator;
1253 
1254 
1258  typedef boost::iterator_range<const_iterator> LineRange;
1259 
1260 
1269  const LineRange
1270  get_lines() const;
1271 
1272 
1300  bool
1301  is_consistent_in_parallel(const std::vector<IndexSet> &locally_owned_dofs,
1302  const IndexSet & locally_active_dofs,
1303  const MPI_Comm mpi_communicator,
1304  const bool verbose = false) const;
1305 
1306 
1325  size_type,
1326  << "The specified line " << arg1 << " does not exist.");
1333  size_type,
1334  size_type,
1335  double,
1336  double,
1337  << "The entry for the indices " << arg1 << " and " << arg2
1338  << " already exists, but the values " << arg3 << " (old) and "
1339  << arg4 << " (new) differ "
1340  << "by " << (arg4 - arg3) << ".");
1347  int,
1348  int,
1349  << "You tried to constrain DoF " << arg1 << " to DoF " << arg2
1350  << ", but that one is also constrained. This is not allowed!");
1357  size_type,
1358  << "Degree of freedom " << arg1
1359  << " is constrained from both object in a merge operation.");
1366  size_type,
1367  << "In the given argument a degree of freedom is constrained "
1368  << "to another DoF with number " << arg1
1369  << ", which however is constrained by this object. This is not"
1370  << " allowed.");
1377  size_type,
1378  << "The index set given to this constraint matrix indicates "
1379  << "constraints for degree of freedom " << arg1
1380  << " should not be stored by this object, but a constraint "
1381  << "is being added.");
1382 
1389  size_type,
1390  size_type,
1391  << "The index set given to this constraint matrix indicates "
1392  << "constraints using degree of freedom " << arg2
1393  << " should not be stored by this object, but a constraint "
1394  << "for degree of freedom " << arg1 << " uses it.");
1395 
1402  int,
1403  int,
1404  << "While distributing the constraint for DoF " << arg1
1405  << ", it turns out that one of the processors "
1406  << "who own the " << arg2 << " degrees of freedom that x_"
1407  << arg1 << " is constrained against does not know about "
1408  << "the constraint on x_" << arg1
1409  << ". Did you not initialize the ConstraintMatrix "
1410  << "with the appropriate locally_relevant set so "
1411  << "that every processor who owns a DoF that constrains "
1412  << "another DoF also knows about this constraint?");
1413 
1414 private:
1426  std::vector<ConstraintLine> lines;
1427 
1460  std::vector<size_type> lines_cache;
1461 
1468 
1472  bool sorted;
1473 
1478  size_type
1479  calculate_line_index(const size_type line) const;
1480 
1485  static bool
1486  check_zero_weight(const std::pair<size_type, double> &p);
1487 
1493 
1498  template <typename MatrixType, typename VectorType>
1499  void
1501  const FullMatrix<typename MatrixType::value_type> &local_matrix,
1502  const Vector<typename VectorType::value_type> & local_vector,
1503  const std::vector<size_type> & local_dof_indices,
1504  MatrixType & global_matrix,
1505  VectorType & global_vector,
1506  bool use_inhomogeneities_for_rhs,
1507  std::integral_constant<bool, false>) const;
1508 
1513  template <typename MatrixType, typename VectorType>
1514  void
1516  const FullMatrix<typename MatrixType::value_type> &local_matrix,
1517  const Vector<typename VectorType::value_type> & local_vector,
1518  const std::vector<size_type> & local_dof_indices,
1519  MatrixType & global_matrix,
1520  VectorType & global_vector,
1521  bool use_inhomogeneities_for_rhs,
1522  std::integral_constant<bool, true>) const;
1523 
1528  template <typename SparsityPatternType>
1529  void
1530  add_entries_local_to_global(const std::vector<size_type> &local_dof_indices,
1531  SparsityPatternType & sparsity_pattern,
1532  const bool keep_constrained_entries,
1533  const Table<2, bool> &dof_mask,
1534  std::integral_constant<bool, false>) const;
1535 
1540  template <typename SparsityPatternType>
1541  void
1542  add_entries_local_to_global(const std::vector<size_type> &local_dof_indices,
1543  SparsityPatternType & sparsity_pattern,
1544  const bool keep_constrained_entries,
1545  const Table<2, bool> &dof_mask,
1546  std::integral_constant<bool, true>) const;
1547 
1555  void
1556  make_sorted_row_list(const std::vector<size_type> & local_dof_indices,
1557  internals::GlobalRowsFromLocal &global_rows) const;
1558 
1566  void
1567  make_sorted_row_list(const std::vector<size_type> &local_dof_indices,
1568  std::vector<size_type> & active_dofs) const;
1569 
1573  template <typename MatrixScalar, typename VectorScalar>
1574  typename ProductType<VectorScalar, MatrixScalar>::type
1575  resolve_vector_entry(const size_type i,
1576  const internals::GlobalRowsFromLocal &global_rows,
1577  const Vector<VectorScalar> & local_vector,
1578  const std::vector<size_type> & local_dof_indices,
1579  const FullMatrix<MatrixScalar> &local_matrix) const;
1580 };
1581 
1582 
1583 
1584 /* ---------------- template and inline functions ----------------- */
1585 
1586 inline ConstraintMatrix::ConstraintMatrix(const IndexSet &local_constraints) :
1587  lines(),
1588  local_lines(local_constraints),
1589  sorted(false)
1590 {
1591  // make sure the IndexSet is compressed. Otherwise this can lead to crashes
1592  // that are hard to find (only happen in release mode).
1593  // see tests/mpi/constraint_matrix_crash_01
1595 }
1596 
1597 
1598 
1600  const ConstraintMatrix &constraint_matrix) :
1601  Subscriptor(),
1602  lines(constraint_matrix.lines),
1603  lines_cache(constraint_matrix.lines_cache),
1604  local_lines(constraint_matrix.local_lines),
1605  sorted(constraint_matrix.sorted)
1606 {}
1607 
1608 
1609 inline void
1610 ConstraintMatrix::add_line(const size_type line)
1611 {
1612  Assert(sorted == false, ExcMatrixIsClosed());
1613 
1614  // the following can happen when we compute with distributed meshes and dof
1615  // handlers and we constrain a degree of freedom whose number we don't have
1616  // locally. if we don't abort here the program will try to allocate several
1617  // terabytes of memory to resize the various arrays below :-)
1619  const size_type line_index = calculate_line_index(line);
1620 
1621  // check whether line already exists; it may, in which case we can just quit
1622  if (is_constrained(line))
1623  return;
1624 
1625  // if necessary enlarge vector of existing entries for cache
1626  if (line_index >= lines_cache.size())
1627  lines_cache.resize(
1628  std::max(2 * static_cast<size_type>(lines_cache.size()), line_index + 1),
1630 
1631  // push a new line to the end of the list
1632  lines.emplace_back();
1633  lines.back().index = line;
1634  lines.back().inhomogeneity = 0.;
1635  lines_cache[line_index] = lines.size() - 1;
1636 }
1637 
1638 
1639 
1640 inline void
1641 ConstraintMatrix::add_entry(const size_type line,
1642  const size_type column,
1643  const double value)
1644 {
1645  Assert(sorted == false, ExcMatrixIsClosed());
1646  Assert(line != column,
1647  ExcMessage("Can't constrain a degree of freedom to itself"));
1648 
1649  // Ensure that the current line is present in the cache:
1650  const size_type line_index = calculate_line_index(line);
1651  Assert(line_index < lines_cache.size(),
1652  ExcMessage("The current ConstraintMatrix does not contain the line "
1653  "for the current entry. Call ConstraintMatrix::add_line "
1654  "before calling this function."));
1655 
1656  // if in debug mode, check whether an entry for this column already exists
1657  // and if it's the same as the one entered at present
1658  //
1659  // in any case: exit the function if an entry for this column already
1660  // exists, since we don't want to enter it twice
1662  ExcInternalError());
1664  ExcColumnNotStoredHere(line, column));
1665  ConstraintLine *line_ptr = &lines[lines_cache[line_index]];
1666  Assert(line_ptr->index == line, ExcInternalError());
1667  for (ConstraintLine::Entries::const_iterator p = line_ptr->entries.begin();
1668  p != line_ptr->entries.end();
1669  ++p)
1670  if (p->first == column)
1671  {
1672  Assert(std::fabs(p->second - value) < 1.e-14,
1673  ExcEntryAlreadyExists(line, column, p->second, value));
1674  return;
1675  }
1676 
1677  line_ptr->entries.emplace_back(column, value);
1678 }
1679 
1680 
1681 
1682 inline void
1683 ConstraintMatrix::set_inhomogeneity(const size_type line, const double value)
1684 {
1685  const size_type line_index = calculate_line_index(line);
1686  Assert(line_index < lines_cache.size() &&
1687  lines_cache[line_index] != numbers::invalid_size_type,
1688  ExcMessage("call add_line() before calling set_inhomogeneity()"));
1689  Assert(lines_cache[line_index] < lines.size(), ExcInternalError());
1690  ConstraintLine *line_ptr = &lines[lines_cache[line_index]];
1691  line_ptr->inhomogeneity = value;
1692 }
1693 
1694 
1695 
1698 {
1699  return lines.size();
1700 }
1701 
1702 
1703 
1704 inline bool
1705 ConstraintMatrix::is_constrained(const size_type index) const
1706 {
1707  const size_type line_index = calculate_line_index(index);
1708  return ((line_index < lines_cache.size()) &&
1709  (lines_cache[line_index] != numbers::invalid_size_type));
1710 }
1711 
1712 
1713 
1714 inline bool
1716 {
1717  // check whether the entry is constrained. could use is_constrained, but
1718  // that means computing the line index twice
1719  const size_type line_index = calculate_line_index(index);
1720  if (line_index >= lines_cache.size() ||
1721  lines_cache[line_index] == numbers::invalid_size_type)
1722  return false;
1723  else
1724  {
1725  Assert(lines_cache[line_index] < lines.size(), ExcInternalError());
1726  return !(lines[lines_cache[line_index]].inhomogeneity == 0);
1727  }
1728 }
1729 
1730 
1731 
1732 inline const std::vector<std::pair<types::global_dof_index, double>> *
1733 ConstraintMatrix::get_constraint_entries(const size_type line) const
1734 {
1735  // check whether the entry is constrained. could use is_constrained, but
1736  // that means computing the line index twice
1737  const size_type line_index = calculate_line_index(line);
1738  if (line_index >= lines_cache.size() ||
1739  lines_cache[line_index] == numbers::invalid_size_type)
1740  return nullptr;
1741  else
1742  return &lines[lines_cache[line_index]].entries;
1743 }
1744 
1745 
1746 
1747 inline double
1748 ConstraintMatrix::get_inhomogeneity(const size_type line) const
1749 {
1750  // check whether the entry is constrained. could use is_constrained, but
1751  // that means computing the line index twice
1752  const size_type line_index = calculate_line_index(line);
1753  if (line_index >= lines_cache.size() ||
1754  lines_cache[line_index] == numbers::invalid_size_type)
1755  return 0;
1756  else
1757  return lines[lines_cache[line_index]].inhomogeneity;
1758 }
1759 
1760 
1761 
1763 ConstraintMatrix::calculate_line_index(const size_type line) const
1764 {
1765  // IndexSet is unused (serial case)
1766  if (!local_lines.size())
1767  return line;
1768 
1770 
1771  return local_lines.index_within_set(line);
1772 }
1773 
1774 
1775 
1776 inline bool
1777 ConstraintMatrix::can_store_line(size_type line_index) const
1778 {
1779  return !local_lines.size() || local_lines.is_element(line_index);
1780 }
1781 
1782 
1783 
1784 inline const IndexSet &
1786 {
1787  return local_lines;
1788 }
1789 
1790 
1791 
1792 template <class VectorType>
1793 inline void
1795  const double value,
1796  VectorType &global_vector) const
1797 {
1798  Assert(lines.empty() || sorted == true, ExcMatrixNotClosed());
1799 
1800  if (is_constrained(index) == false)
1801  global_vector(index) += value;
1802  else
1803  {
1804  const ConstraintLine &position =
1806  for (size_type j = 0; j < position.entries.size(); ++j)
1807  global_vector(position.entries[j].first) +=
1808  value * position.entries[j].second;
1809  }
1810 }
1811 
1812 
1813 template <typename ForwardIteratorVec,
1814  typename ForwardIteratorInd,
1815  class VectorType>
1816 inline void
1818  ForwardIteratorVec local_vector_begin,
1819  ForwardIteratorVec local_vector_end,
1820  ForwardIteratorInd local_indices_begin,
1821  VectorType & global_vector) const
1822 {
1823  Assert(lines.empty() || sorted == true, ExcMatrixNotClosed());
1824  for (; local_vector_begin != local_vector_end;
1825  ++local_vector_begin, ++local_indices_begin)
1826  {
1827  if (is_constrained(*local_indices_begin) == false)
1828  internal::ElementAccess<VectorType>::add(
1829  *local_vector_begin, *local_indices_begin, global_vector);
1830  else
1831  {
1832  const ConstraintLine &position =
1833  lines[lines_cache[calculate_line_index(*local_indices_begin)]];
1834  for (size_type j = 0; j < position.entries.size(); ++j)
1835  internal::ElementAccess<VectorType>::add(
1836  (*local_vector_begin) * position.entries[j].second,
1837  position.entries[j].first,
1838  global_vector);
1839  }
1840  }
1841 }
1842 
1843 
1844 template <class InVector, class OutVector>
1845 inline void
1847  const InVector & local_vector,
1848  const std::vector<size_type> &local_dof_indices,
1849  OutVector & global_vector) const
1850 {
1851  Assert(local_vector.size() == local_dof_indices.size(),
1852  ExcDimensionMismatch(local_vector.size(), local_dof_indices.size()));
1853  distribute_local_to_global(local_vector.begin(),
1854  local_vector.end(),
1855  local_dof_indices.begin(),
1856  global_vector);
1857 }
1858 
1859 
1860 
1861 template <typename ForwardIteratorVec,
1862  typename ForwardIteratorInd,
1863  class VectorType>
1864 inline void
1865 ConstraintMatrix::get_dof_values(const VectorType & global_vector,
1866  ForwardIteratorInd local_indices_begin,
1867  ForwardIteratorVec local_vector_begin,
1868  ForwardIteratorVec local_vector_end) const
1869 {
1870  Assert(lines.empty() || sorted == true, ExcMatrixNotClosed());
1871  for (; local_vector_begin != local_vector_end;
1872  ++local_vector_begin, ++local_indices_begin)
1873  {
1874  if (is_constrained(*local_indices_begin) == false)
1875  *local_vector_begin = global_vector(*local_indices_begin);
1876  else
1877  {
1878  const ConstraintLine &position =
1879  lines[lines_cache[calculate_line_index(*local_indices_begin)]];
1880  typename VectorType::value_type value = position.inhomogeneity;
1881  for (size_type j = 0; j < position.entries.size(); ++j)
1882  value += (global_vector(position.entries[j].first) *
1883  position.entries[j].second);
1884  *local_vector_begin = value;
1885  }
1886  }
1887 }
1888 
1889 
1890 template <typename MatrixType>
1891 class BlockMatrixBase;
1892 template <typename SparsityPatternType>
1894 template <typename number>
1895 class BlockSparseMatrixEZ;
1896 
1915 template <typename MatrixType>
1917 {
1918 private:
1919  struct yes_type
1920  {
1921  char c[1];
1922  };
1923  struct no_type
1924  {
1925  char c[2];
1926  };
1927 
1933  template <typename T>
1934  static yes_type
1936 
1941  template <typename T>
1942  static yes_type
1944 
1949  template <typename T>
1950  static yes_type
1952 
1957  static no_type
1959 
1960 public:
1966  static const bool value =
1967  (sizeof(check_for_block_matrix((MatrixType *)nullptr)) == sizeof(yes_type));
1968 };
1969 
1970 
1971 // instantiation of the static member
1972 template <typename MatrixType>
1974 
1975 
1976 template <typename MatrixType>
1977 inline void
1979  const FullMatrix<typename MatrixType::value_type> &local_matrix,
1980  const std::vector<size_type> & local_dof_indices,
1981  MatrixType & global_matrix) const
1982 {
1983  // create a dummy and hand on to the function actually implementing this
1984  // feature in the cm.templates.h file.
1985  Vector<typename MatrixType::value_type> dummy(0);
1987  local_matrix,
1988  dummy,
1989  local_dof_indices,
1990  global_matrix,
1991  dummy,
1992  false,
1993  std::integral_constant<bool, IsBlockMatrix<MatrixType>::value>());
1994 }
1995 
1996 
1997 
1998 template <typename MatrixType, typename VectorType>
1999 inline void
2001  const FullMatrix<typename MatrixType::value_type> &local_matrix,
2002  const Vector<typename VectorType::value_type> & local_vector,
2003  const std::vector<size_type> & local_dof_indices,
2004  MatrixType & global_matrix,
2005  VectorType & global_vector,
2006  bool use_inhomogeneities_for_rhs) const
2007 {
2008  // enter the internal function with the respective block information set,
2009  // the actual implementation follows in the cm.templates.h file.
2011  local_matrix,
2012  local_vector,
2013  local_dof_indices,
2014  global_matrix,
2015  global_vector,
2016  use_inhomogeneities_for_rhs,
2017  std::integral_constant<bool, IsBlockMatrix<MatrixType>::value>());
2018 }
2019 
2020 
2021 
2022 template <typename SparsityPatternType>
2023 inline void
2025  const std::vector<size_type> &local_dof_indices,
2026  SparsityPatternType & sparsity_pattern,
2027  const bool keep_constrained_entries,
2028  const Table<2, bool> & dof_mask) const
2029 {
2030  // enter the internal function with the respective block information set,
2031  // the actual implementation follows in the cm.templates.h file.
2033  local_dof_indices,
2034  sparsity_pattern,
2035  keep_constrained_entries,
2036  dof_mask,
2037  std::integral_constant<bool, IsBlockMatrix<SparsityPatternType>::value>());
2038 }
2039 
2040 
2041 DEAL_II_NAMESPACE_CLOSE
2042 
2043 #endif
const types::global_dof_index invalid_size_type
Definition: types.h:182
double get_inhomogeneity(const size_type line) const
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:422
std::vector< size_type > lines_cache
static::ExceptionBase & ExcDoFIsConstrainedFromBothObjects(size_type arg1)
static::ExceptionBase & ExcEntryAlreadyExists(size_type arg1, size_type arg2, double arg3, double arg4)
std::vector< ConstraintLine >::const_iterator const_iterator
void merge(const ConstraintMatrix &other_constraints, const MergeConflictBehavior merge_conflict_behavior=no_conflicts_allowed, const bool allow_different_local_lines=false)
bool can_store_line(const size_type line_index) const
bool is_constrained(const size_type index) const
static bool check_zero_weight(const std::pair< size_type, double > &p)
void condense(SparsityPattern &sparsity) const
boost::iterator_range< const_iterator > LineRange
const std::vector< std::pair< size_type, double > > * get_constraint_entries(const size_type line) const
const IndexSet & get_local_lines() const
void distribute(VectorType &vec) const
bool has_inhomogeneities() const
static::ExceptionBase & ExcDoFIsConstrainedToConstrainedDoF(size_type arg1)
static const bool value
void print(std::ostream &out) const
void add_entries(const size_type line, const std::vector< std::pair< size_type, double >> &col_val_pairs)
types::global_dof_index size_type
size_type size() const
Definition: index_set.h:1608
static yes_type check_for_block_matrix(const BlockMatrixBase< T > *)
size_type max_constraint_indirections() const
bool are_identity_constrained(const size_type index1, const size_type index2) const
static::ExceptionBase & ExcLineInexistant(size_type arg1)
void add_line(const size_type line)
void write_dot(std::ostream &) const
void distribute_local_to_global(const InVector &local_vector, const std::vector< size_type > &local_dof_indices, OutVector &global_vector) const
static::ExceptionBase & ExcMessage(std::string arg1)
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:410
void add_entry(const size_type line, const size_type column, const double value)
unsigned int global_dof_index
Definition: types.h:89
ConstraintMatrix & operator=(const ConstraintMatrix &)=delete
#define Assert(cond, exc)
Definition: exceptions.h:1227
static::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static::ExceptionBase & ExcColumnNotStoredHere(size_type arg1, size_type arg2)
ConstraintMatrix(const IndexSet &local_constraints=IndexSet())
bool operator<(const ConstraintLine &) const
bool is_inhomogeneously_constrained(const size_type index) const
#define DeclException0(Exception0)
Definition: exceptions.h:387
bool operator==(const ConstraintLine &) const
static const Table< 2, bool > default_empty_table
void compress() const
Definition: index_set.h:1616
size_type calculate_line_index(const size_type line) const
void add_lines(const std::vector< bool > &lines)
static::ExceptionBase & ExcIncorrectConstraint(int arg1, int arg2)
std::vector< std::pair< size_type, double > > Entries
static::ExceptionBase & ExcMatrixNotClosed()
static::ExceptionBase & ExcDoFConstrainedToConstrainedDoF(int arg1, int arg2)
std::size_t memory_consumption() const
std::size_t memory_consumption() const
void copy_from(const ConstraintMatrix &other)
void shift(const size_type offset)
void set_zero(VectorType &vec) const
void resolve_indices(std::vector< types::global_dof_index > &indices) const
size_type n_constraints() const
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
Definition: exceptions.h:446
ProductType< VectorScalar, MatrixScalar >::type resolve_vector_entry(const size_type i, const internals::GlobalRowsFromLocal &global_rows, const Vector< VectorScalar > &local_vector, const std::vector< size_type > &local_dof_indices, const FullMatrix< MatrixScalar > &local_matrix) const
const LineRange get_lines() const
size_type index_within_set(const size_type global_index) const
Definition: index_set.h:1830
static::ExceptionBase & ExcRowNotStoredHere(size_type arg1)
bool is_element(const size_type index) const
Definition: index_set.h:1672
Definition: table.h:37
void reinit(const IndexSet &local_constraints=IndexSet())
void add_selected_constraints(const ConstraintMatrix &constraints_in, const IndexSet &filter)
bool is_consistent_in_parallel(const std::vector< IndexSet > &locally_owned_dofs, const IndexSet &locally_active_dofs, const MPI_Comm mpi_communicator, const bool verbose=false) const
void set_inhomogeneity(const size_type line, const double value)
void add_entries_local_to_global(const std::vector< size_type > &local_dof_indices, SparsityPatternType &sparsity_pattern, const bool keep_constrained_entries=true, const Table< 2, bool > &dof_mask=default_empty_table) const
std::vector< ConstraintLine > lines
void make_sorted_row_list(const std::vector< size_type > &local_dof_indices, internals::GlobalRowsFromLocal &global_rows) const
void get_dof_values(const VectorType &global_vector, ForwardIteratorInd local_indices_begin, ForwardIteratorVec local_vector_begin, ForwardIteratorVec local_vector_end) const
static::ExceptionBase & ExcMatrixIsClosed()
static::ExceptionBase & ExcInternalError()
void serialize(Archive &ar, const unsigned int)
bool is_identity_constrained(const size_type index) const