15#ifndef dealii_tensor_product_matrix_creator_h
16#define dealii_tensor_product_matrix_creator_h
64 template <
int dim,
typename Number>
65 std::pair<std::array<FullMatrix<Number>, dim>,
66 std::array<FullMatrix<Number>, dim>>
78 template <
int dim,
typename Number>
79 std::pair<std::array<FullMatrix<Number>, dim>,
80 std::array<FullMatrix<Number>, dim>>
101 template <
typename Number>
104 const unsigned int n,
114 template <
typename Number>
124 dof_handler.distribute_dofs(fe);
128 const unsigned int n_dofs_1D = fe.n_dofs_per_cell();
139 fe_values.reinit(tria.begin());
141 const auto lexicographic_to_hierarchic_numbering =
143 FETools::hierarchic_to_lexicographic_numbering<1>(
144 fe.tensor_degree()));
146 for (
const unsigned int q_index : fe_values.quadrature_point_indices())
147 for (
const unsigned int i : fe_values.dof_indices())
148 for (
const unsigned int j : fe_values.dof_indices())
151 (fe_values.shape_value(lexicographic_to_hierarchic_numbering[i],
153 fe_values.shape_value(lexicographic_to_hierarchic_numbering[
j],
155 fe_values.JxW(q_index));
158 (fe_values.shape_grad(lexicographic_to_hierarchic_numbering[i],
160 fe_values.shape_grad(lexicographic_to_hierarchic_numbering[
j],
162 fe_values.JxW(q_index));
172 template <
int dim,
typename Number>
173 std::pair<std::array<FullMatrix<Number>, dim>,
174 std::array<FullMatrix<Number>, dim>>
183 const auto create_reference_mass_and_stiffness_matrices =
184 internal::create_reference_mass_and_stiffness_matrices<Number>(
188 std::get<0>(create_reference_mass_and_stiffness_matrices);
190 std::get<1>(create_reference_mass_and_stiffness_matrices);
192 std::get<2>(create_reference_mass_and_stiffness_matrices);
204 std::array<FullMatrix<Number>, dim>
Ms;
205 std::array<FullMatrix<Number>, dim>
Ks;
207 for (
unsigned int d = 0; d < dim; ++d)
213 for (
unsigned int i = 0; i <
n_dofs_1D; ++i)
228 for (
unsigned int i = 0; i <
n_overlap; ++i)
265 for (
unsigned int i = 0; i <
n_overlap; ++i)
302 template <
int dim,
typename Number>
303 std::pair<std::array<FullMatrix<Number>, dim>,
304 std::array<FullMatrix<Number>, dim>>
316 for (
unsigned int d = 0; d < dim; ++d)
319 if ((cell->at_boundary(2 * d) ==
false) ||
320 cell->has_periodic_neighbor(2 * d))
329 const auto bid = cell->face(2 * d)->boundary_id();
349 if ((cell->at_boundary(2 * d + 1) ==
false) ||
350 cell->has_periodic_neighbor(2 * d + 1))
358 const auto bid = cell->face(2 * d + 1)->boundary_id();
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
#define AssertThrow(cond, exc)
@ update_values
Shape function values.
@ update_JxW_values
Transformed quadrature weights.
@ update_gradients
Shape function gradients.
void hyper_cube(Triangulation< dim, spacedim > &tria, const double left=0., const double right=1., const bool colorize=false)
void clear_row_and_column(const unsigned int n_dofs_1D_with_overlap, const unsigned int n, FullMatrix< Number > &matrix)
std::tuple< FullMatrix< Number >, FullMatrix< Number >, bool > create_reference_mass_and_stiffness_matrices(const FiniteElement< 1 > &fe, const Quadrature< 1 > &quadrature)
std::pair< std::array< FullMatrix< Number >, dim >, std::array< FullMatrix< Number >, dim > > create_laplace_tensor_product_matrix(const FiniteElement< 1 > &fe, const Quadrature< 1 > &quadrature, const ::ndarray< LaplaceBoundaryType, dim, 2 > &boundary_ids, const ::ndarray< double, dim, 3 > &cell_extent, const unsigned int n_overlap=1)
std::vector< Integer > invert_permutation(const std::vector< Integer > &permutation)