18 #ifndef __deal2__sparse_matrix_templates_h
19 #define __deal2__sparse_matrix_templates_h
22 #include <deal.II/base/config.h>
23 #include <deal.II/base/template_constraints.h>
24 #include <deal.II/base/parallel.h>
25 #include <deal.II/base/thread_management.h>
26 #include <deal.II/base/multithread_info.h>
27 #include <deal.II/base/utilities.h>
28 #include <deal.II/lac/sparse_matrix.h>
29 #include <deal.II/lac/trilinos_sparse_matrix.h>
30 #include <deal.II/lac/vector.h>
31 #include <deal.II/lac/full_matrix.h>
32 #include <deal.II/lac/compressed_simple_sparsity_pattern.h>
33 #include <deal.II/lac/vector_memory.h>
42 #include <deal.II/base/std_cxx1x/bind.h>
46 DEAL_II_NAMESPACE_OPEN
49 template <
typename number>
52 cols(0,
"SparseMatrix"),
59 template <
typename number>
63 cols(0,
"SparseMatrix"),
74 template <
typename number>
87 template <
typename number>
90 cols(0,
"SparseMatrix"),
99 template <
typename number>
103 cols(0,
"SparseMatrix"),
117 template <
typename number>
135 void zero_subrange (
const size_type begin,
139 std::memset (dst+begin,0,(end-begin)*
sizeof(T));
146 template <
typename number>
153 Assert (cols->compressed || cols->empty(), SparsityPattern::ExcNotCompressed());
162 const size_type matrix_size = cols->n_nonzero_elements();
164 internal::SparseMatrix::minimum_parallel_grain_size *
165 (cols->n_nonzero_elements()+m()) / m();
166 if (matrix_size>grain_size)
168 std_cxx1x::bind(&internal::SparseMatrix::template
169 zero_subrange<number>,
170 std_cxx1x::_1, std_cxx1x::_2,
173 else if (matrix_size > 0)
174 std::memset (&val[0], 0, matrix_size*
sizeof(number));
181 template <
typename number>
185 Assert (cols->n_rows() ==
id.m(),
187 Assert (cols->n_cols() ==
id.n(),
199 template <
typename number>
215 if (N > max_len || max_len == 0)
228 template <
typename number>
233 if (val)
delete[] val;
240 template <
typename number>
247 return cols->empty();
252 template <
typename number>
257 return cols->row_length(row);
262 template <
typename number>
267 return cols->n_nonzero_elements ();
272 template <
typename number>
279 const size_type nnz_alloc = n_nonzero_elements();
281 if (std::fabs(val[i]) > threshold)
288 template <
typename number>
293 Assert (cols->rows == cols->cols, ExcNotQuadratic());
299 number *val_ptr = &val[cols->rowstart[row]];
302 const size_type *colnum_ptr = &cols->colnums[cols->rowstart[row]+1];
303 const number *
const val_end_of_row = &val[cols->rowstart[row+1]];
306 while ((val_ptr != val_end_of_row) && (*colnum_ptr<row))
310 const number mean_value = (*val_ptr +
311 val[(*cols)(*colnum_ptr,row)]) / 2.0;
315 *val_ptr = mean_value;
316 set (*colnum_ptr, row, mean_value);
327 template <
typename number>
328 template <
typename somenumber>
334 Assert (cols == matrix.
cols, ExcDifferentSparsityPatterns());
336 std::copy (&matrix.
val[0], &matrix.
val[cols->n_nonzero_elements()],
344 template <
typename number>
345 template <
typename somenumber>
355 if (matrix(row,col) != 0)
356 set (row, col, matrix(row,col));
361 #ifdef DEAL_II_WITH_TRILINOS
363 template <
typename number>
373 std::vector < TrilinosScalar > value_cache;
374 std::vector<size_type> colnum_cache;
376 for (
size_type row = 0; row < matrix.
m(); ++row)
378 value_cache.resize(matrix.
n());
379 colnum_cache.resize(matrix.
n());
388 reinterpret_cast<TrilinosWrappers::types::int_type *>(&(colnum_cache[0])));
389 Assert (ierr==0, ExcTrilinosError(ierr));
392 value_cache.resize(ncols);
393 colnum_cache.resize(ncols);
407 template <
typename number>
408 template <
typename somenumber>
415 Assert (cols == matrix.
cols, ExcDifferentSparsityPatterns());
417 number *val_ptr = &val[0];
418 const somenumber *matrix_ptr = &matrix.
val[0];
419 const number *
const end_ptr = &val[cols->n_nonzero_elements()];
421 while (val_ptr != end_ptr)
422 *val_ptr++ += factor **matrix_ptr++;
439 template <
typename number,
442 void vmult_on_subrange (
const size_type begin_row,
444 const number *values,
445 const std::size_t *rowstart,
451 const number *val_ptr = &values[rowstart[begin_row]];
452 const size_type *colnum_ptr = &colnums[rowstart[begin_row]];
453 typename OutVector::iterator dst_ptr = dst.begin() + begin_row;
456 for (
size_type row=begin_row; row<end_row; ++row)
458 typename OutVector::value_type s = 0.;
459 const number *
const val_end_of_row = &values[rowstart[row+1]];
460 while (val_ptr != val_end_of_row)
461 s += *val_ptr++ * src(*colnum_ptr++);
465 for (
size_type row=begin_row; row<end_row; ++row)
467 typename OutVector::value_type s = *dst_ptr;
468 const number *
const val_end_of_row = &values[rowstart[row+1]];
469 while (val_ptr != val_end_of_row)
470 s += *val_ptr++ * src(*colnum_ptr++);
478 template <
typename number>
479 template <
typename number2>
484 const number2 *values,
485 const bool elide_zero_values,
486 const bool col_indices_are_sorted)
495 if (elide_zero_values ==
false && col_indices_are_sorted ==
true &&
502 Assert (col_indices[i] > col_indices[i-1],
503 ExcMessage(
"List of indices is unsorted or contains duplicates."));
507 &cols->colnums[cols->rowstart[row]];
508 const size_type row_length_1 = cols->row_length(row)-1;
509 number *val_ptr = &val[cols->rowstart[row]];
520 const size_type diag = diag_pos - col_indices;
522 if (diag != n_cols && *diag_pos == row)
524 val_ptr[0] += *(values+(diag_pos-col_indices));
532 while (this_cols[counter]<col_indices[i] && counter<row_length_1)
535 Assert (this_cols[counter] == col_indices[i] || values[i] == 0,
536 ExcInvalidIndex(row,col_indices[i]));
538 val_ptr[counter] += values[i];
542 for (
size_type i=post_diag; i<n_cols; ++i)
544 while (this_cols[counter]<col_indices[i] && counter<row_length_1)
547 Assert (this_cols[counter] == col_indices[i] || values[i] == 0,
548 ExcInvalidIndex(row,col_indices[i]));
550 val_ptr[counter] += values[i];
553 Assert (counter < cols->row_length(row),
554 ExcMessage(
"Specified invalid column indices in add function."));
561 while (this_cols[counter]<col_indices[i] && counter<row_length_1)
564 Assert (this_cols[counter] == col_indices[i] || values[i] == 0,
565 ExcInvalidIndex(row,col_indices[i]));
567 val_ptr[counter] += values[i];
569 Assert (counter < cols->row_length(row),
570 ExcMessage(
"Specified invalid column indices in add function."));
578 const size_type *
const my_cols = cols->colnums;
580 const size_type next_row_index = cols->rowstart[row+1];
584 const number value = values[j];
588 if (elide_zero_values==
true && value == 0)
599 if (index < next_row_index && my_cols[index] == col_indices[j])
602 index = cols->operator()(row, col_indices[j]);
610 Assert (value == 0., ExcInvalidIndex(row,col_indices[j]));
622 template <
typename number>
623 template <
typename number2>
628 const number2 *values,
629 const bool elide_zero_values)
632 Assert (row < m(), ExcInvalidIndex1(row));
637 const size_type *my_cols = cols->colnums;
638 std::size_t index = cols->rowstart[row], next_index = index;
639 const std::size_t next_row_index = cols->rowstart[row+1];
641 if (elide_zero_values ==
true)
645 const number value = values[j];
655 if (index != next_row_index && my_cols[index] == col_indices[j])
658 next_index = cols->operator()(row, col_indices[j]);
666 Assert (
false, ExcInvalidIndex(row,col_indices[j]));
681 const number value = values[j];
684 if (index != next_row_index && my_cols[index] == col_indices[j])
685 goto set_value_checked;
687 next_index = cols->operator()(row, col_indices[j]);
691 Assert (value == 0., ExcInvalidIndex(row,col_indices[j]));
705 template <
typename number>
706 template <
class OutVector,
class InVector>
709 const InVector &src)
const
719 std_cxx1x::bind (&internal::SparseMatrix::vmult_on_subrange
720 <number,InVector,OutVector>,
721 std_cxx1x::_1, std_cxx1x::_2,
725 std_cxx1x::cref(src),
728 internal::SparseMatrix::minimum_parallel_grain_size);
733 template <
typename number>
734 template <
class OutVector,
class InVector>
737 const InVector &src)
const
750 for (
size_type j=cols->rowstart[i]; j<cols->rowstart[i+1] ; j++)
753 dst(p) += val[j] * src(i);
760 template <
typename number>
761 template <
class OutVector,
class InVector>
764 const InVector &src)
const
774 std_cxx1x::bind (&internal::SparseMatrix::vmult_on_subrange
775 <number,InVector,OutVector>,
776 std_cxx1x::_1, std_cxx1x::_2,
780 std_cxx1x::cref(src),
783 internal::SparseMatrix::minimum_parallel_grain_size);
788 template <
typename number>
789 template <
class OutVector,
class InVector>
792 const InVector &src)
const
802 for (
size_type j=cols->rowstart[i]; j<cols->rowstart[i+1] ; j++)
805 dst(p) += val[j] * src(i);
825 template <
typename number,
827 number matrix_norm_sqr_on_subrange (
const size_type begin_row,
829 const number *values,
830 const std::size_t *rowstart,
836 for (
size_type i=begin_row; i<end_row; ++i)
839 for (
size_type j=rowstart[i]; j<rowstart[i+1] ; j++)
840 s += values[j] * v(colnums[j]);
850 template <
typename number>
851 template <
typename somenumber>
861 parallel::accumulate_from_subranges<number>
862 (std_cxx1x::bind (&internal::SparseMatrix::matrix_norm_sqr_on_subrange
864 std_cxx1x::_1, std_cxx1x::_2,
865 val, cols->rowstart, cols->colnums,
868 internal::SparseMatrix::minimum_parallel_grain_size);
888 template <
typename number,
890 number matrix_scalar_product_on_subrange (
const size_type begin_row,
892 const number *values,
893 const std::size_t *rowstart,
900 for (
size_type i=begin_row; i<end_row; ++i)
903 for (
size_type j=rowstart[i]; j<rowstart[i+1] ; j++)
904 s += values[j] * v(colnums[j]);
914 template <
typename number>
915 template <
typename somenumber>
926 parallel::accumulate_from_subranges<number>
927 (std_cxx1x::bind (&internal::SparseMatrix::matrix_scalar_product_on_subrange
929 std_cxx1x::_1, std_cxx1x::_2,
930 val, cols->rowstart, cols->colnums,
934 internal::SparseMatrix::minimum_parallel_grain_size);
939 template <
typename number>
940 template <
typename numberB,
typename numberC>
945 const bool rebuild_sparsity_C)
const
947 const bool use_vector = V.
size() == n() ?
true :
false;
957 if (rebuild_sparsity_C ==
true)
962 ExcMessage (
"Can't use the same sparsity pattern for "
963 "different matrices if it is to be rebuilt."));
965 ExcMessage (
"Can't use the same sparsity pattern for "
966 "different matrices if it is to be rebuilt."));
981 for (
size_type i = 0; i < csp.n_rows(); ++i)
986 for (; rows != end_rows; ++rows)
1002 csp.add_entries (i, new_cols, end_new_cols,
true);
1005 sp_C.copy_from (csp);
1018 unsigned int max_n_cols_B = 0;
1020 max_n_cols_B = std::max (max_n_cols_B, sp_B.
row_length(i));
1021 std::vector<numberC> new_entries(max_n_cols_B);
1031 for (; rows != end_rows; ++rows)
1041 C.
add (i, *new_cols, A_val *
1043 (use_vector ? V(col) : 1));
1050 numberC *new_ptr = &new_entries[0];
1051 const numberB *B_val_ptr =
1053 const numberB *
const end_cols = &B.
val[sp_B.
rowstart[col+1]];
1054 for (; B_val_ptr != end_cols; ++B_val_ptr)
1055 *new_ptr++ = A_val **B_val_ptr * (use_vector ? V(col) : 1);
1057 C.
add (i, new_ptr-&new_entries[0], new_cols, &new_entries[0],
1066 template <
typename number>
1067 template <
typename numberB,
typename numberC>
1072 const bool rebuild_sparsity_C)
const
1074 const bool use_vector = V.
size() == m() ?
true :
false;
1084 if (rebuild_sparsity_C ==
true)
1089 ExcMessage (
"Can't use the same sparsity pattern for "
1090 "different matrices if it is to be rebuilt."));
1092 ExcMessage (
"Can't use the same sparsity pattern for "
1093 "different matrices if it is to be rebuilt."));
1099 sp_C.reinit (0,0,0);
1123 for (; rows != end_rows; ++rows)
1132 csp.add_entries (row, new_cols, end_new_cols,
true);
1135 sp_C.copy_from (csp);
1148 unsigned int max_n_cols_B = 0;
1150 max_n_cols_B = std::max (max_n_cols_B, sp_B.
row_length(i));
1151 std::vector<numberC> new_entries(max_n_cols_B);
1165 const numberB *
const end_cols = &B.
val[sp_B.
rowstart[i+1]];
1167 for (; rows != end_rows; ++rows)
1174 C.
add (row, i, A_val *
1176 (use_vector ? V(i) : 1));
1181 numberC *new_ptr = &new_entries[0];
1182 const numberB *B_val_ptr =
1184 for (; B_val_ptr != end_cols; ++B_val_ptr)
1185 *new_ptr++ = A_val **B_val_ptr * (use_vector ? V(i) : 1);
1187 C.
add (row, new_ptr-&new_entries[0], new_cols, &new_entries[0],
1195 template <
typename number>
1204 for (
size_type row=0; row<n_rows; ++row)
1205 for (
size_type j=cols->rowstart[row]; j<cols->rowstart[row+1] ; ++j)
1213 template <
typename number>
1220 const number *val_ptr = &val[cols->rowstart[0]];
1224 for (
size_type row=0; row<n_rows; ++row)
1227 const number *
const val_end_of_row = &val[cols->rowstart[row+1]];
1228 while (val_ptr != val_end_of_row)
1238 template <
typename number>
1247 for (
const number *ptr = &val[0];
1248 ptr != &val[cols->rowstart[n_rows]]; ++ptr)
1251 return std::sqrt (norm_sqr);
1271 template <
typename number,
1274 number residual_sqr_on_subrange (
const size_type begin_row,
1276 const number *values,
1277 const std::size_t *rowstart,
1285 for (
size_type i=begin_row; i<end_row; ++i)
1288 for (
size_type j=rowstart[i]; j<rowstart[i+1] ; j++)
1289 s -= values[j] * u(colnums[j]);
1299 template <
typename number>
1300 template <
typename somenumber>
1312 Assert (&u != &dst, ExcSourceEqualsDestination());
1315 std::sqrt (parallel::accumulate_from_subranges<number>
1316 (std_cxx1x::bind (&internal::SparseMatrix::residual_sqr_on_subrange
1318 std_cxx1x::_1, std_cxx1x::_2,
1319 val, cols->rowstart, cols->colnums,
1322 std_cxx1x::ref(dst)),
1324 internal::SparseMatrix::minimum_parallel_grain_size));
1329 template <
typename number>
1330 template <
typename somenumber>
1334 const number om)
const
1343 somenumber *dst_ptr = dst.
begin();
1344 const somenumber *src_ptr = src.
begin();
1345 const std::size_t *rowstart_ptr = &cols->rowstart[0];
1359 for (
size_type i=0; i<n; ++i, ++dst_ptr, ++src_ptr, ++rowstart_ptr)
1360 *dst_ptr = om **src_ptr / val[*rowstart_ptr];
1362 for (
size_type i=0; i<n; ++i, ++dst_ptr, ++src_ptr, ++rowstart_ptr)
1363 *dst_ptr = *src_ptr / val[*rowstart_ptr];
1368 template <
typename number>
1369 template <
typename somenumber>
1374 const std::vector<std::size_t> &pos_right_of_diagonal)
const
1387 const std::size_t *rowstart_ptr = &cols->rowstart[0];
1388 somenumber *dst_ptr = &dst(0);
1393 if (pos_right_of_diagonal.size() != 0)
1395 Assert (pos_right_of_diagonal.size() == dst.
size(),
1399 for (
size_type row=0; row<n; ++row, ++dst_ptr, ++rowstart_ptr)
1401 *dst_ptr = src(row);
1402 const std::size_t first_right_of_diagonal_index =
1403 pos_right_of_diagonal[row];
1404 Assert (first_right_of_diagonal_index <= *(rowstart_ptr+1),
1407 for (
size_type j=(*rowstart_ptr)+1; j<first_right_of_diagonal_index; ++j)
1408 s += val[j] * dst(cols->colnums[j]);
1413 *dst_ptr /= val[*rowstart_ptr];
1416 rowstart_ptr = &cols->rowstart[0];
1418 for ( ; rowstart_ptr!=&cols->rowstart[n]; ++rowstart_ptr, ++dst_ptr)
1419 *dst_ptr *= om*(2.-om)*val[*rowstart_ptr];
1422 rowstart_ptr = &cols->rowstart[n-1];
1423 dst_ptr = &dst(n-1);
1424 for (
int row=n-1; row>=0; --row, --rowstart_ptr, --dst_ptr)
1426 const size_type end_row = *(rowstart_ptr+1);
1427 const size_type first_right_of_diagonal_index
1428 = pos_right_of_diagonal[row];
1430 for (
size_type j=first_right_of_diagonal_index; j<end_row; ++j)
1431 s += val[j] * dst(cols->colnums[j]);
1435 *dst_ptr /= val[*rowstart_ptr];
1444 for (
size_type row=0; row<n; ++row, ++dst_ptr, ++rowstart_ptr)
1446 *dst_ptr = src(row);
1454 const size_type first_right_of_diagonal_index
1456 &cols->colnums[*(rowstart_ptr+1)],
1462 for (
size_type j=(*rowstart_ptr)+1; j<first_right_of_diagonal_index; ++j)
1463 s += val[j] * dst(cols->colnums[j]);
1468 *dst_ptr /= val[*rowstart_ptr];
1471 rowstart_ptr = &cols->rowstart[0];
1473 for (
size_type row=0; row<n; ++row, ++rowstart_ptr, ++dst_ptr)
1474 *dst_ptr *= (2.-om)*val[*rowstart_ptr];
1477 rowstart_ptr = &cols->rowstart[n-1];
1478 dst_ptr = &dst(n-1);
1479 for (
int row=n-1; row>=0; --row, --rowstart_ptr, --dst_ptr)
1481 const size_type end_row = *(rowstart_ptr+1);
1482 const size_type first_right_of_diagonal_index
1484 &cols->colnums[end_row],
1485 static_cast<size_type>(row)) -
1488 for (
size_type j=first_right_of_diagonal_index; j<end_row; ++j)
1489 s += val[j] * dst(cols->colnums[j]);
1492 *dst_ptr /= val[*rowstart_ptr];
1497 template <
typename number>
1498 template <
typename somenumber>
1502 const number om)
const
1512 template <
typename number>
1513 template <
typename somenumber>
1517 const number om)
const
1527 template <
typename number>
1528 template <
typename somenumber>
1531 const number om)
const
1540 somenumber s = dst(row);
1541 for (
size_type j=cols->rowstart[row]; j<cols->rowstart[row+1]; ++j)
1545 s -= val[j] * dst(col);
1549 dst(row) = s * om / val[cols->rowstart[row]];
1554 template <
typename number>
1555 template <
typename somenumber>
1558 const number om)
const
1568 somenumber s = dst(row);
1569 for (
size_type j=cols->rowstart[row]; j<cols->rowstart[row+1]; ++j)
1570 if (cols->colnums[j] > row)
1571 s -= val[j] * dst(cols->colnums[j]);
1574 dst(row) = s * om / val[cols->rowstart[row]];
1584 template <
typename number>
1585 template <
typename somenumber>
1588 const std::vector<size_type> &permutation,
1589 const std::vector<size_type> &inverse_permutation,
1590 const number om)
const
1597 Assert (m() == permutation.size(),
1599 Assert (m() == inverse_permutation.size(),
1602 for (
size_type urow=0; urow<m(); ++urow)
1604 const size_type row = permutation[urow];
1605 somenumber s = dst(row);
1607 for (
size_type j=cols->rowstart[row]; j<cols->rowstart[row+1]; ++j)
1610 if (inverse_permutation[col] < urow)
1612 s -= val[j] * dst(col);
1617 dst(row) = s * om / val[cols->rowstart[row]];
1622 template <
typename number>
1623 template <
typename somenumber>
1626 const std::vector<size_type> &permutation,
1627 const std::vector<size_type> &inverse_permutation,
1628 const number om)
const
1635 Assert (m() == permutation.size(),
1637 Assert (m() == inverse_permutation.size(),
1643 const size_type row = permutation[urow];
1644 somenumber s = dst(row);
1645 for (
size_type j=cols->rowstart[row]; j<cols->rowstart[row+1]; ++j)
1648 if (inverse_permutation[col] > urow)
1649 s -= val[j] * dst(col);
1653 dst(row) = s * om / val[cols->rowstart[row]];
1659 template <
typename number>
1660 template <
typename somenumber>
1664 const number om)
const
1684 precondition_Jacobi (*w, *w, om);
1690 template <
typename number>
1691 template <
typename somenumber>
1695 const number om)
const
1705 somenumber s = b(row);
1706 for (
size_type j=cols->rowstart[row]; j<cols->rowstart[row+1]; ++j)
1708 s -= val[j] * v(cols->colnums[j]);
1711 v(row) += s * om / val[cols->rowstart[row]];
1717 template <
typename number>
1718 template <
typename somenumber>
1722 const number om)
const
1730 for (
int row=m()-1; row>=0; --row)
1732 somenumber s = b(row);
1733 for (
size_type j=cols->rowstart[row]; j<cols->rowstart[row+1]; ++j)
1735 s -= val[j] * v(cols->colnums[j]);
1738 v(row) += s * om / val[cols->rowstart[row]];
1744 template <
typename number>
1745 template <
typename somenumber>
1749 const number om)
const
1757 template <
typename number>
1758 template <
typename somenumber>
1761 const number om)
const
1778 for (j=cols->rowstart[i]; j<cols->rowstart[i+1] ; j++)
1783 if (i>j) s += val[j] * dst(p);
1788 dst(i) /= val[cols->rowstart[i]];
1791 for (
int i=n-1; i>=0; i--)
1794 for (j=cols->rowstart[i]; j<cols->rowstart[i+1] ; j++)
1799 if (static_cast<size_type>(i)<j) s += val[j] * dst(p);
1803 dst(i) -= s * om / val[cols->rowstart[i]];
1809 template <
typename number>
1819 template <
typename number>
1821 const unsigned int precision,
1822 const bool scientific,
1823 const unsigned int width_,
1824 const char *zero_string,
1825 const double denominator)
const
1830 unsigned int width = width_;
1832 std::ios::fmtflags old_flags = out.flags();
1833 unsigned int old_precision = out.precision (precision);
1837 out.setf (std::ios::scientific, std::ios::floatfield);
1839 width = precision+7;
1843 out.setf (std::ios::fixed, std::ios::floatfield);
1845 width = precision+2;
1852 out << std::setw(width)
1853 << val[cols->operator()(i,j)] * denominator <<
' ';
1855 out << std::setw(width) << zero_string <<
' ';
1861 out.precision(old_precision);
1862 out.flags (old_flags);
1867 template <
typename number>
1869 const double threshold)
const
1879 else if (std::fabs(val[cols->operator()(i,j)]) > threshold)
1890 template <
typename number>
1898 out <<
'[' << max_len <<
"][";
1900 out.write (reinterpret_cast<const char *>(&val[0]),
1901 reinterpret_cast<const char *>(&val[max_len])
1902 - reinterpret_cast<const char *>(&val[0]));
1910 template <
typename number>
1930 val =
new number[max_len];
1933 in.read (reinterpret_cast<char *>(&val[0]),
1934 reinterpret_cast<char *>(&val[max_len])
1935 - reinterpret_cast<char *>(&val[0]));
1941 template <
typename number>
1947 template <
typename number>
1951 return max_len*
static_cast<std::size_t
>(
sizeof(number)) +
sizeof(*this);
1955 DEAL_II_NAMESPACE_CLOSE
somenumber matrix_norm_square(const Vector< somenumber > &v) const
Iterator lower_bound(Iterator first, Iterator last, const T &val)
void precondition_Jacobi(Vector< somenumber > &dst, const Vector< somenumber > &src, const number omega=1.) const
static const number & conjugate(const number &x)
void Tvmult(OutVector &dst, const InVector &src) const
void vmult(OutVector &dst, const InVector &src) const
void PSOR(Vector< somenumber > &v, const std::vector< size_type > &permutation, const std::vector< size_type > &inverse_permutation, const number om=1.) const
#define AssertDimension(dim1, dim2)
void Tvmult_add(OutVector &dst, const InVector &src) const
numbers::NumberTraits< number >::real_type real_type
void mmult(SparseMatrix< numberC > &C, const SparseMatrix< numberB > &B, const Vector< number > &V=Vector< number >(), const bool rebuild_sparsity_pattern=true) const
const Epetra_CrsMatrix & trilinos_matrix() const
void SSOR(Vector< somenumber > &v, const number omega=1.) const
size_type get_row_length(const size_type row) const
::ExceptionBase & ExcMessage(std::string arg1)
void set(const size_type i, const size_type j, const number value)
void block_read(std::istream &in)
size_type n_actually_nonzero_elements(const double threshold=0.) const
size_type n_nonzero_elements() const
#define AssertThrow(cond, exc)
static real_type abs(const number &x)
void Tmmult(SparseMatrix< numberC > &C, const SparseMatrix< numberB > &B, const Vector< number > &V=Vector< number >(), const bool rebuild_sparsity_pattern=true) const
void SOR(Vector< somenumber > &v, const number om=1.) const
void apply_to_subranges(const RangeType &begin, const typename identity< RangeType >::type &end, const Function &f, const unsigned int grainsize)
real_type linfty_norm() const
bool is_finite(const double x)
void print_pattern(std::ostream &out, const double threshold=0.) const
::ExceptionBase & ExcIO()
real_type l1_norm() const
void precondition_TSOR(Vector< somenumber > &dst, const Vector< somenumber > &src, const number om=1.) const
SparseMatrix< number > & operator=(const SparseMatrix< number > &)
unsigned int global_dof_index
#define Assert(cond, exc)
std::size_t memory_consumption() const
const SparsityPattern & get_sparsity_pattern() const
static bool equal(const T *p1, const T *p2)
void precondition_SSOR(Vector< somenumber > &dst, const Vector< somenumber > &src, const number omega=1., const std::vector< std::size_t > &pos_right_of_diagonal=std::vector< std::size_t >()) const
unsigned int row_length(const size_type row) const
void precondition_SOR(Vector< somenumber > &dst, const Vector< somenumber > &src, const number om=1.) const
::ExceptionBase & ExcInvalidConstructorCall()
size_type n_nonzero_elements() const
somenumber residual(Vector< somenumber > &dst, const Vector< somenumber > &x, const Vector< somenumber > &b) const
somenumber matrix_scalar_product(const Vector< somenumber > &u, const Vector< somenumber > &v) const
void add(const size_type i, const size_type j, const number value)
void TPSOR(Vector< somenumber > &v, const std::vector< size_type > &permutation, const std::vector< size_type > &inverse_permutation, const number om=1.) const
SparseMatrix< number > & copy_from(const SparseMatrix< somenumber > &source)
virtual void reinit(const SparsityPattern &sparsity)
void Jacobi_step(Vector< somenumber > &v, const Vector< somenumber > &b, const number om=1.) const
void vmult_add(OutVector &dst, const InVector &src) const
types::global_dof_index size_type
::ExceptionBase & ExcNumberNotFinite()
void TSOR(Vector< somenumber > &v, const number om=1.) const
void SSOR_step(Vector< somenumber > &v, const Vector< somenumber > &b, const number om=1.) const
SmartPointer< const SparsityPattern, SparseMatrix< number > > cols
void block_write(std::ostream &out) const
::ExceptionBase & ExcNotImplemented()
::ExceptionBase & ExcNotInitialized()
void print_formatted(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const unsigned int width=0, const char *zero_string=" ", const double denominator=1.) const
::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
::ExceptionBase & ExcScalarAssignmentOnlyForZeroValue()
void SOR_step(Vector< somenumber > &v, const Vector< somenumber > &b, const number om=1.) const
::ExceptionBase & ExcInternalError()
::ExceptionBase & ExcDivideByZero()
real_type frobenius_norm() const
VectorizedArray< Number > max(const ::VectorizedArray< Number > &x, const ::VectorizedArray< Number > &y)
unsigned int row_length(const size_type row) const
void compress(::VectorOperation::values)
real_type linfty_norm() const
static const size_type invalid_entry
void TSOR_step(Vector< somenumber > &v, const Vector< somenumber > &b, const number om=1.) const