17 #ifndef dealii_matrix_free_operators_h 18 #define dealii_matrix_free_operators_h 21 #include <deal.II/base/exceptions.h> 22 #include <deal.II/base/subscriptor.h> 23 #include <deal.II/base/vectorization.h> 25 #include <deal.II/lac/diagonal_matrix.h> 26 #include <deal.II/lac/la_parallel_vector.h> 28 #include <deal.II/matrix_free/fe_evaluation.h> 29 #include <deal.II/matrix_free/matrix_free.h> 31 #include <deal.II/multigrid/mg_constrained_dofs.h> 34 DEAL_II_NAMESPACE_OPEN
44 template <
typename VectorType>
45 typename std::enable_if<IsBlockVector<VectorType>::value,
47 n_blocks(
const VectorType &vector)
49 return vector.n_blocks();
52 template <
typename VectorType>
53 typename std::enable_if<!IsBlockVector<VectorType>::value,
55 n_blocks(
const VectorType &)
60 template <
typename VectorType>
61 typename std::enable_if<IsBlockVector<VectorType>::value,
62 typename VectorType::BlockType &>::type
63 subblock(VectorType &vector,
unsigned int block_no)
65 return vector.block(block_no);
68 template <
typename VectorType>
69 typename std::enable_if<IsBlockVector<VectorType>::value,
70 const typename VectorType::BlockType &>::type
71 subblock(
const VectorType &vector,
unsigned int block_no)
74 return vector.block(block_no);
77 template <
typename VectorType>
78 typename std::enable_if<!IsBlockVector<VectorType>::value,
80 subblock(VectorType &vector,
unsigned int)
85 template <
typename VectorType>
86 typename std::enable_if<!IsBlockVector<VectorType>::value,
87 const VectorType &>::type
88 subblock(
const VectorType &vector,
unsigned int)
93 template <
typename VectorType>
94 typename std::enable_if<IsBlockVector<VectorType>::value,
void>::type
95 collect_sizes(VectorType &vector)
97 vector.collect_sizes();
100 template <
typename VectorType>
101 typename std::enable_if<!IsBlockVector<VectorType>::value,
void>::type
102 collect_sizes(
const VectorType &)
208 virtual ~
Base()
override =
default;
235 const std::vector<unsigned int> &selected_row_blocks =
236 std::vector<unsigned int>(),
237 const std::vector<unsigned int> &selected_column_blocks =
238 std::vector<unsigned int>());
256 const unsigned int level,
257 const std::vector<unsigned int> &selected_row_blocks =
258 std::vector<unsigned int>());
276 const std::vector<MGConstrainedDoFs> &mg_constrained_dofs,
277 const unsigned int level,
278 const std::vector<unsigned int> & selected_row_blocks =
279 std::vector<unsigned int>());
297 vmult_interface_down(VectorType &dst,
const VectorType &src)
const;
303 vmult_interface_up(VectorType &dst,
const VectorType &src)
const;
309 vmult(VectorType &dst,
const VectorType &src)
const;
315 Tvmult(VectorType &dst,
const VectorType &src)
const;
321 vmult_add(VectorType &dst,
const VectorType &src)
const;
327 Tvmult_add(VectorType &dst,
const VectorType &src)
const;
334 el(
const unsigned int row,
const unsigned int col)
const;
341 memory_consumption()
const;
347 initialize_dof_vector(VectorType &vec)
const;
357 compute_diagonal() = 0;
362 std::shared_ptr<const MatrixFree<dim, value_type>>
363 get_matrix_free()
const;
368 const std::shared_ptr<DiagonalMatrix<VectorType>> &
369 get_matrix_diagonal_inverse()
const;
374 const std::shared_ptr<DiagonalMatrix<VectorType>> &
375 get_matrix_diagonal()
const;
384 precondition_Jacobi(VectorType & dst,
385 const VectorType &src,
394 preprocess_constraints(VectorType &dst,
const VectorType &src)
const;
401 postprocess_constraints(VectorType &dst,
const VectorType &src)
const;
408 set_constrained_entries_to_one(VectorType &dst)
const;
414 apply_add(VectorType &dst,
const VectorType &src)
const = 0;
422 Tapply_add(VectorType &dst,
const VectorType &src)
const;
427 std::shared_ptr<const MatrixFree<dim, value_type>>
data;
462 mutable std::vector<std::vector<std::pair<value_type, value_type>>>
476 mult_add(VectorType & dst,
477 const VectorType &src,
488 adjust_ghost_range_if_necessary(
const VectorType &vec,
489 const bool is_row)
const;
530 template <
typename OperatorType>
559 initialize(
const OperatorType &operator_in);
564 template <
typename VectorType>
566 vmult(VectorType &dst,
const VectorType &src)
const;
571 template <
typename VectorType>
573 Tvmult(VectorType &dst,
const VectorType &src)
const;
578 template <
typename VectorType>
580 initialize_dof_vector(VectorType &vec)
const;
613 int n_components = 1,
614 typename Number =
double>
635 const unsigned int n_actual_components,
645 fill_inverse_JxW_values(
673 int n_q_points_1d = fe_degree + 1,
674 int n_components = 1,
699 compute_diagonal()
override;
708 apply_add(VectorType &dst,
const VectorType &src)
const override;
717 const VectorType & src,
718 const std::pair<unsigned int, unsigned int> &cell_range)
const;
737 int n_q_points_1d = fe_degree + 1,
738 int n_components = 1,
815 &scalar_coefficient);
826 std::shared_ptr<Table<2, VectorizedArray<value_type>>>
836 apply_add(VectorType &dst,
const VectorType &src)
const;
845 const VectorType & src,
846 const std::pair<unsigned int, unsigned int> &cell_range)
const;
855 const unsigned int &,
856 const std::pair<unsigned int, unsigned int> &cell_range)
const;
862 do_operation_on_cell(
865 const unsigned int cell)
const;
877 template <
int dim,
int fe_degree,
int n_components,
typename Number>
884 for (
unsigned int i = 0, c = 0; i < shapes_1d.
m(); ++i)
885 for (
unsigned int j = 0; j < shapes_1d.
n(); ++j, ++c)
888 const unsigned int stride = (fe_degree + 2) / 2;
890 for (
unsigned int i = 0; i < stride; ++i)
891 for (
unsigned int q = 0; q < (fe_degree + 2) / 2; ++q)
894 0.5 * (shapes_1d(i, q) + shapes_1d(i, fe_degree - q));
896 0.5 * (shapes_1d(i, q) - shapes_1d(i, fe_degree - q));
898 if (fe_degree % 2 == 0)
899 for (
unsigned int q = 0; q < (fe_degree + 2) / 2; ++q)
900 inverse_shape[fe_degree / 2 * stride + q] = shapes_1d(fe_degree / 2, q);
905 template <
int dim,
int fe_degree,
int n_components,
typename Number>
911 constexpr
unsigned int dofs_per_cell =
Utilities::pow(fe_degree + 1, dim);
912 Assert(inverse_jxw.size() > 0 && inverse_jxw.size() % dofs_per_cell == 0,
914 "Expected diagonal to be a multiple of scalar dof per cells"));
918 const unsigned int previous_size = inverse_jxw.size();
919 inverse_jxw.resize(dofs_per_cell);
923 inverse_jxw.resize_fast(previous_size);
924 for (
unsigned int q = 0; q < dofs_per_cell; ++q)
925 inverse_jxw[q] = 1. / inverse_jxw[q];
927 for (
unsigned int q = dofs_per_cell; q < previous_size;)
928 for (
unsigned int i = 0; i < dofs_per_cell; ++i, ++q)
929 inverse_jxw[q] = inverse_jxw[i];
934 template <
int dim,
int fe_degree,
int n_components,
typename Number>
938 const unsigned int n_actual_components,
942 constexpr
unsigned int dofs_per_component =
944 Assert(inverse_coefficients.size() > 0 &&
945 inverse_coefficients.size() % dofs_per_component == 0,
947 "Expected diagonal to be a multiple of scalar dof per cells"));
948 if (inverse_coefficients.size() != dofs_per_component)
950 inverse_coefficients.size());
961 const unsigned int shift_coefficient =
962 inverse_coefficients.size() > dofs_per_component ? dofs_per_component : 0;
965 for (
unsigned int d = 0; d < n_actual_components; ++d)
971 evaluator.template hessians<0, false, false>(in, temp_data_field);
974 evaluator.template hessians<1, false, false>(temp_data_field, out);
978 evaluator.template hessians<2, false, false>(out,
980 for (
unsigned int q = 0; q < dofs_per_component; ++q)
981 temp_data_field[q] *= inv_coefficient[q];
982 evaluator.template hessians<2, true, false>(temp_data_field,
986 for (
unsigned int q = 0; q < dofs_per_component; ++q)
987 out[q] *= inv_coefficient[q];
989 evaluator.template hessians<1, true, false>(out, temp_data_field);
993 for (
unsigned int q = 0; q < dofs_per_component; ++q)
994 temp_data_field[q] *= inv_coefficient[q];
996 evaluator.template hessians<0, true, false>(temp_data_field, out);
998 inv_coefficient += shift_coefficient;
1003 template <
int dim,
typename VectorType>
1006 , have_interface_matrices(false)
1011 template <
int dim,
typename VectorType>
1024 template <
int dim,
typename VectorType>
1037 template <
int dim,
typename VectorType>
1047 template <
int dim,
typename VectorType>
1050 const unsigned int col)
const 1057 return 1.0 / (*inverse_diagonal_entries)(row, row);
1062 template <
int dim,
typename VectorType>
1068 for (
unsigned int i = 0; i < BlockHelper::n_blocks(vec); ++i)
1071 if (!BlockHelper::subblock(vec, index)
1072 .partitioners_are_compatible(
1073 *
data->get_dof_info(index).vector_partitioner))
1074 data->initialize_dof_vector(BlockHelper::subblock(vec, index), index);
1076 Assert(BlockHelper::subblock(vec, index)
1077 .partitioners_are_globally_compatible(
1078 *
data->get_dof_info(index).vector_partitioner),
1081 BlockHelper::collect_sizes(vec);
1086 template <
int dim,
typename VectorType>
1091 const std::vector<unsigned int> &given_row_selection,
1092 const std::vector<unsigned int> &given_column_selection)
1098 if (given_row_selection.empty())
1099 for (
unsigned int i = 0; i < data_->n_components(); ++i)
1103 for (
unsigned int i = 0; i < given_row_selection.size(); ++i)
1106 for (
unsigned int j = 0; j < given_row_selection.size(); ++j)
1108 Assert(given_row_selection[j] != given_row_selection[i],
1109 ExcMessage(
"Given row indices must be unique"));
1114 if (given_column_selection.size() == 0)
1118 for (
unsigned int i = 0; i < given_column_selection.size(); ++i)
1121 for (
unsigned int j = 0; j < given_column_selection.size(); ++j)
1123 Assert(given_column_selection[j] != given_column_selection[i],
1124 ExcMessage(
"Given column indices must be unique"));
1139 template <
int dim,
typename VectorType>
1145 const unsigned int level,
1146 const std::vector<unsigned int> &given_row_selection)
1148 std::vector<MGConstrainedDoFs> mg_constrained_dofs_vector(
1149 1, mg_constrained_dofs);
1150 initialize(data_, mg_constrained_dofs_vector, level, given_row_selection);
1155 template <
int dim,
typename VectorType>
1159 const std::vector<MGConstrainedDoFs> & mg_constrained_dofs,
1160 const unsigned int level,
1161 const std::vector<unsigned int> & given_row_selection)
1168 if (given_row_selection.empty())
1169 for (
unsigned int i = 0; i < data_->n_components(); ++i)
1173 for (
unsigned int i = 0; i < given_row_selection.size(); ++i)
1176 for (
unsigned int j = 0; j < given_row_selection.size(); ++j)
1178 Assert(given_row_selection[j] != given_row_selection[i],
1179 ExcMessage(
"Given row indices must be unique"));
1196 if (data_->n_macro_cells() > 0)
1199 data_->get_cell_iterator(0, 0, j)->level());
1203 std::vector<types::global_dof_index> interface_indices;
1204 mg_constrained_dofs[j]
1205 .get_refinement_edge_indices(level)
1206 .fill_index_vector(interface_indices);
1212 for (
unsigned int i = 0; i < interface_indices.size(); ++i)
1213 if (locally_owned.
is_element(interface_indices[i]))
1219 data->get_vector_partitioner()->get_mpi_communicator()) > 0;
1225 template <
int dim,
typename VectorType>
1229 for (
unsigned int j = 0; j < BlockHelper::n_blocks(dst); ++j)
1231 const std::vector<unsigned int> &constrained_dofs =
1233 for (
unsigned int i = 0; i < constrained_dofs.size(); ++i)
1234 BlockHelper::subblock(dst, j).local_element(constrained_dofs[i]) = 1.;
1236 BlockHelper::subblock(dst, j).local_element(
1243 template <
int dim,
typename VectorType>
1254 template <
int dim,
typename VectorType>
1263 template <
int dim,
typename VectorType>
1266 const VectorType &src)
const 1273 template <
int dim,
typename VectorType>
1276 const VectorType &src,
1277 const bool is_row)
const 1280 for (
unsigned int i = 0; i < BlockHelper::n_blocks(src); ++i)
1282 const unsigned int mf_component =
1285 if (BlockHelper::subblock(src, i).get_partitioner().get() ==
1286 data->get_dof_info(mf_component).vector_partitioner.get())
1292 BlockHelper::subblock(src, i).get_partitioner()->local_size() ==
1293 data->get_dof_info(mf_component).vector_partitioner->local_size(),
1294 ExcMessage(
"The vector passed to the vmult() function does not have " 1295 "the correct size for compatibility with MatrixFree."));
1300 BlockHelper::subblock(src, i));
1301 BlockHelper::subblock(const_cast<VectorType &>(src), i)
1302 .reinit(
data->get_dof_info(mf_component).vector_partitioner);
1303 BlockHelper::subblock(const_cast<VectorType &>(src), i)
1304 .copy_locally_owned_data_from(copy_vec);
1310 template <
int dim,
typename VectorType>
1313 const VectorType &src)
const 1321 for (
unsigned int j = 0; j < BlockHelper::n_blocks(dst); ++j)
1326 BlockHelper::subblock(src, j).local_element(
1328 BlockHelper::subblock(dst, j).local_element(
1330 BlockHelper::subblock(const_cast<VectorType &>(src), j)
1338 template <
int dim,
typename VectorType>
1341 const VectorType &src,
1345 AssertDimension(BlockHelper::n_blocks(dst), BlockHelper::n_blocks(src));
1357 template <
int dim,
typename VectorType>
1360 const VectorType &src)
const 1362 for (
unsigned int j = 0; j < BlockHelper::n_blocks(dst); ++j)
1364 const std::vector<unsigned int> &constrained_dofs =
1366 for (
unsigned int i = 0; i < constrained_dofs.size(); ++i)
1367 BlockHelper::subblock(dst, j).local_element(constrained_dofs[i]) +=
1368 BlockHelper::subblock(src, j).local_element(constrained_dofs[i]);
1373 for (
unsigned int j = 0; j < BlockHelper::n_blocks(dst); ++j)
1377 BlockHelper::subblock(const_cast<VectorType &>(src), j)
1380 BlockHelper::subblock(dst, j).local_element(
1390 template <
int dim,
typename VectorType>
1393 const VectorType &src)
const 1407 for (
unsigned int j = 0; j < BlockHelper::n_blocks(dst); ++j)
1411 BlockHelper::subblock(src, j).local_element(
1413 BlockHelper::subblock(dst, j).local_element(
1415 BlockHelper::subblock(const_cast<VectorType &>(src), j)
1421 for (
unsigned int j = 0; j < BlockHelper::n_blocks(dst); ++j)
1427 BlockHelper::subblock(dst, j).local_element(c) = 0.;
1431 BlockHelper::subblock(const_cast<VectorType &>(src), j)
1432 .local_element(edge_constrained_indices[j][i]) =
1435 for (; c < BlockHelper::subblock(dst, j).local_size(); ++c)
1436 BlockHelper::subblock(dst, j).local_element(c) = 0.;
1442 template <
int dim,
typename VectorType>
1445 const VectorType &src)
const 1457 VectorType src_cpy(src);
1458 for (
unsigned int j = 0; j < BlockHelper::n_blocks(dst); ++j)
1464 BlockHelper::subblock(src_cpy, j).local_element(c) = 0.;
1467 for (; c < BlockHelper::subblock(src_cpy, j).local_size(); ++c)
1468 BlockHelper::subblock(src_cpy, j).local_element(c) = 0.;
1473 for (
unsigned int j = 0; j < BlockHelper::n_blocks(dst); ++j)
1475 BlockHelper::subblock(dst, j).local_element(
1481 template <
int dim,
typename VectorType>
1492 template <
int dim,
typename VectorType>
1503 template <
int dim,
typename VectorType>
1513 template <
int dim,
typename VectorType>
1514 const std::shared_ptr<DiagonalMatrix<VectorType>> &
1525 template <
int dim,
typename VectorType>
1526 const std::shared_ptr<DiagonalMatrix<VectorType>> &
1536 template <
int dim,
typename VectorType>
1539 const VectorType &src)
const 1546 template <
int dim,
typename VectorType>
1550 const VectorType & src,
1563 template <
typename OperatorType>
1566 , mf_base_operator(nullptr)
1571 template <
typename OperatorType>
1580 template <
typename OperatorType>
1589 template <
typename OperatorType>
1590 template <
typename VectorType>
1593 const VectorType &src)
const 1595 #ifndef DEAL_II_MSVC 1597 std::is_same<typename VectorType::value_type, value_type>::value,
1598 "The vector type must be based on the same value type as this" 1609 template <
typename OperatorType>
1610 template <
typename VectorType>
1613 const VectorType &src)
const 1615 #ifndef DEAL_II_MSVC 1617 std::is_same<typename VectorType::value_type, value_type>::value,
1618 "The vector type must be based on the same value type as this" 1629 template <
typename OperatorType>
1630 template <
typename VectorType>
1633 VectorType &vec)
const 1648 typename VectorType>
1651 :
Base<dim, VectorType>()
1660 typename VectorType>
1670 VectorType &inverse_diagonal_vector =
1675 inverse_diagonal_vector = Number(1.);
1676 apply_add(diagonal_vector, inverse_diagonal_vector);
1679 inverse_diagonal_vector = diagonal_vector;
1681 const unsigned int local_size = inverse_diagonal_vector.local_size();
1682 for (
unsigned int i = 0; i < local_size; ++i)
1683 inverse_diagonal_vector.local_element(i) =
1684 Number(1.) / inverse_diagonal_vector.local_element(i);
1686 inverse_diagonal_vector.update_ghost_values();
1687 diagonal_vector.update_ghost_values();
1696 typename VectorType>
1713 typename VectorType>
1719 const VectorType & src,
1720 const std::pair<unsigned int, unsigned int> &cell_range)
const 1725 for (
unsigned int cell = cell_range.first; cell < cell_range.second; ++cell)
1730 for (
unsigned int q = 0; q < phi.
n_q_points; ++q)
1744 typename VectorType>
1747 :
Base<dim, VectorType>()
1756 typename VectorType>
1771 typename VectorType>
1775 const std::shared_ptr<
1777 &scalar_coefficient_)
1788 typename VectorType>
1809 typename VectorType>
1817 unsigned int dummy = 0;
1820 VectorType &inverse_diagonal_vector =
1832 inverse_diagonal_vector = diagonal_vector;
1834 for (
unsigned int i = 0; i < inverse_diagonal_vector.local_size(); ++i)
1835 if (std::abs(inverse_diagonal_vector.local_element(i)) >
1836 std::sqrt(std::numeric_limits<Number>::epsilon()))
1837 inverse_diagonal_vector.local_element(i) =
1838 1. / inverse_diagonal_vector.local_element(i);
1840 inverse_diagonal_vector.local_element(i) = 1.;
1842 inverse_diagonal_vector.update_ghost_values();
1843 diagonal_vector.update_ghost_values();
1852 typename VectorType>
1863 namespace Implementation
1865 template <
typename Number>
1869 for (
unsigned int v = 0; v < VectorizedArray<Number>::n_array_elements;
1884 typename VectorType>
1893 const unsigned int cell)
const 1895 phi.evaluate(
false,
true,
false);
1898 for (
unsigned int q = 0; q < phi.n_q_points; ++q)
1901 ExcMessage(
"Coefficient must be non-negative"));
1903 phi.get_gradient(q),
1909 for (
unsigned int q = 0; q < phi.n_q_points; ++q)
1911 phi.submit_gradient(phi.get_gradient(q), q);
1914 phi.integrate(
false,
true);
1923 typename VectorType>
1929 const VectorType & src,
1930 const std::pair<unsigned int, unsigned int> &cell_range)
const 1935 for (
unsigned int cell = cell_range.first; cell < cell_range.second; ++cell)
1949 typename VectorType>
1955 const unsigned int &,
1956 const std::pair<unsigned int, unsigned int> &cell_range)
const 1961 for (
unsigned int cell = cell_range.first; cell < cell_range.second; ++cell)
1976 local_diagonal_vector[i];
1985 DEAL_II_NAMESPACE_CLOSE
typename OperatorType::value_type value_type
const std::shared_ptr< DiagonalMatrix< VectorType > > & get_matrix_diagonal() const
void fill_inverse_JxW_values(AlignedVector< VectorizedArray< Number >> &inverse_jxw) const
static const unsigned int invalid_unsigned_int
std::vector< unsigned int > selected_rows
void do_operation_on_cell(FEEvaluation< dim, fe_degree, n_q_points_1d, n_components, value_type > &phi, const unsigned int cell) const
void vmult_interface_up(VectorType &dst, const VectorType &src) const
#define AssertDimension(dim1, dim2)
void initialize_dof_vector(VectorType &vec) const
void Tvmult(VectorType &dst, const VectorType &src) const
typename VectorType::size_type size_type
constexpr unsigned int pow(const unsigned int base, const unsigned int iexp)
void evaluate(const bool evaluate_values, const bool evaluate_gradients, const bool evaluate_hessians=false)
void integrate(const bool integrate_values, const bool integrate_gradients)
typename VectorType::value_type value_type
const internal::MatrixFreeFunctions::ShapeInfo< VectorizedArray< Number > > & get_shape_info() const
static constexpr unsigned int static_dofs_per_cell
#define AssertIndexRange(index, range)
std::vector< std::vector< std::pair< value_type, value_type > > > edge_constrained_values
void local_apply_cell(const MatrixFree< dim, value_type > &data, VectorType &dst, const VectorType &src, const std::pair< unsigned int, unsigned int > &cell_range) const
void adjust_ghost_range_if_necessary(const VectorType &vec, const bool is_row) const
static::ExceptionBase & ExcNotInitialized()
std::shared_ptr< DiagonalMatrix< VectorType > > inverse_diagonal_entries
#define AssertThrow(cond, exc)
std::shared_ptr< DiagonalMatrix< VectorType > > diagonal_entries
void read_dof_values(const VectorType &src, const unsigned int first_index=0)
static constexpr unsigned int n_components
value_type get_value(const unsigned int q_point) const
void postprocess_constraints(VectorType &dst, const VectorType &src) const
typename OperatorType::size_type size_type
void resize(const size_type size_in)
AlignedVector< Number > shape_values
const unsigned int dofs_per_component
void submit_value(const value_type val_in, const unsigned int q_point)
const FEEvaluationBase< dim, n_components, Number > & fe_eval
virtual void compute_diagonal()
std::vector< unsigned int > selected_columns
void Tvmult(VectorType &dst, const VectorType &src) const
static::ExceptionBase & ExcMessage(std::string arg1)
const unsigned int n_q_points
virtual std::size_t memory_consumption() const
void preprocess_constraints(VectorType &dst, const VectorType &src) const
void apply(const AlignedVector< VectorizedArray< Number >> &inverse_coefficient, const unsigned int n_actual_components, const VectorizedArray< Number > *in_array, VectorizedArray< Number > *out_array) const
#define Assert(cond, exc)
std::shared_ptr< const MatrixFree< dim, value_type > > get_matrix_free() const
void set_coefficient(const std::shared_ptr< Table< 2, VectorizedArray< value_type >>> &scalar_coefficient)
void fill_JxW_values(AlignedVector< VectorizedArray< Number >> &JxW_values) const
std::vector< std::vector< unsigned int > > edge_constrained_indices
AlignedVector< VectorizedArray< Number > > inverse_shape
SymmetricTensor< rank_, dim, Number > transpose(const SymmetricTensor< rank_, dim, Number > &t)
void set_constrained_entries_to_one(VectorType &dst) const
void vmult(VectorType &dst, const VectorType &src) const
void reinit(const unsigned int cell_batch_index)
void distribute_local_to_global(VectorType &dst, const unsigned int first_index=0) const
virtual void apply_add(VectorType &dst, const VectorType &src) const =0
std::shared_ptr< Table< 2, VectorizedArray< value_type > > > scalar_coefficient
void local_diagonal_cell(const MatrixFree< dim, value_type > &data, VectorType &dst, const unsigned int &, const std::pair< unsigned int, unsigned int > &cell_range) const
void vmult_add(VectorType &dst, const VectorType &src) const
virtual void compute_diagonal() override
CellwiseInverseMassMatrix(const FEEvaluationBase< dim, n_components, Number > &fe_eval)
const VectorizedArray< Number > * begin_dof_values() const
value_type el(const unsigned int row, const unsigned int col) const
std::shared_ptr< Table< 2, VectorizedArray< value_type > > > get_coefficient()
void local_apply_cell(const MatrixFree< dim, value_type > &data, VectorType &dst, const VectorType &src, const std::pair< unsigned int, unsigned int > &cell_range) const
void mult_add(VectorType &dst, const VectorType &src, const bool transpose) const
size_type index_within_set(const size_type global_index) const
void Tvmult_add(VectorType &dst, const VectorType &src) const
std::shared_ptr< const MatrixFree< dim, value_type > > data
void precondition_Jacobi(VectorType &dst, const VectorType &src, const value_type omega) const
virtual void apply_add(VectorType &dst, const VectorType &src) const override
bool is_element(const size_type index) const
static::ExceptionBase & ExcNotImplemented()
void vmult_interface_down(VectorType &dst, const VectorType &src) const
void initialize(const OperatorType &operator_in)
SmartPointer< const OperatorType > mf_base_operator
virtual void Tapply_add(VectorType &dst, const VectorType &src) const
void initialize_dof_vector(VectorType &vec) const
void vmult(VectorType &dst, const VectorType &src) const
const std::shared_ptr< DiagonalMatrix< VectorType > > & get_matrix_diagonal_inverse() const
T max(const T &t, const MPI_Comm &mpi_communicator)
bool have_interface_matrices
void initialize(std::shared_ptr< const MatrixFree< dim, value_type >> data, const std::vector< unsigned int > &selected_row_blocks=std::vector< unsigned int >(), const std::vector< unsigned int > &selected_column_blocks=std::vector< unsigned int >())
static::ExceptionBase & ExcInternalError()
virtual void apply_add(VectorType &dst, const VectorType &src) const