From: MichaƂ Wichrowski Date: Thu, 17 Apr 2025 11:40:20 +0000 (+0200) Subject: fixed compile X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=5cd5978ce5d7f17fe420c2af4ff45ec0f4336a73;p=dealii.git fixed compile --- diff --git a/include/deal.II/numerics/matrix_creator.templates.h b/include/deal.II/numerics/matrix_creator.templates.h index 203e9dde98..1d56931c8e 100644 --- a/include/deal.II/numerics/matrix_creator.templates.h +++ b/include/deal.II/numerics/matrix_creator.templates.h @@ -2270,148 +2270,6 @@ namespace MatrixCreator coefficient, constraints); } - - - - FullMatrix - create_1d_cell_mass_matrix(const FiniteElement<1> &fe, - const double &h, - const std::pair include_endpoints, - std::vector numbering) - { - if (dynamic_cast *>(&fe) == nullptr && - numbering.size() == 0) - { - Assert( - include_endpoints.first == true && include_endpoints.second == true, - ExcMessage( - "You tried to genereate a 1D mass matrix with excluding boundary " - "dofs for a non-DGQ element without providing a numbering.")); - } - - - if (numbering.size() == 0) - { - numbering.resize(fe.dofs_per_cell); - std::iota(numbering.begin(), numbering.end(), 0); - } - const unsigned int degree = fe.degree; - const unsigned int n_dofs_per_cell = fe.dofs_per_cell; - const double &JxW = h; - QGauss<1> quadrature(degree + 1); - - FullMatrix cell_mass_matrix(n_dofs_per_cell, n_dofs_per_cell); - cell_mass_matrix = 0; - - unsigned int start_dof = include_endpoints.first ? 0 : 1; - unsigned int end_dof = - include_endpoints.second ? n_dofs_per_cell : n_dofs_per_cell - 1; - const unsigned int shift = include_endpoints.first ? 0 : 1; - - for (unsigned int i = start_dof; i < end_dof; ++i) - for (unsigned int j = start_dof; j < end_dof; ++j) - for (unsigned int q = 0; q < quadrature.size(); ++q) - cell_mass_matrix(i - shift, j - shift) += - (fe.shape_value(numbering[i], quadrature.point(q)) * - fe.shape_value(numbering[j], quadrature.point(q))) * - JxW * quadrature.weight(q); - - return cell_mass_matrix; - } - - - FullMatrix - create_1d_cell_derivative_matrix( - const FiniteElement<1> &fe, - const double &h, - const std::pair include_endpoints, - std::vector numbering) - { - if (numbering.size() == 0) - { - numbering.resize(fe.dofs_per_cell); - std::iota(numbering.begin(), numbering.end(), 0); - } - - const unsigned int degree = fe.degree; - const unsigned int n_dofs_per_cell = fe.dofs_per_cell; - const double &JxW = h; - QGauss<1> quadrature(degree + 1); - - FullMatrix cell_matrix(n_dofs_per_cell, n_dofs_per_cell); - cell_matrix = 0; - - unsigned int start_dof = include_endpoints.first ? 0 : 1; - unsigned int end_dof = - include_endpoints.second ? n_dofs_per_cell : n_dofs_per_cell - 1; - const unsigned int shift = include_endpoints.first ? 0 : 1; - - for (unsigned int i = start_dof; i < end_dof; ++i) - for (unsigned int j = start_dof; j < end_dof; ++j) - for (unsigned int q = 0; q < quadrature.size(); ++q) - cell_matrix(i - shift, j - shift) += - (fe.shape_grad(numbering[i], quadrature.point(q)) / h * - fe.shape_grad(numbering[j], quadrature.point(q))) / - h * JxW * quadrature.weight(q); - - return cell_matrix; - } - - FullMatrix - create_1D_discretization_matrix(FullMatrix &cell_matrix, - const unsigned int &n_cells, - const unsigned int &overlap, - const std::pair include_endpoints) - { - const unsigned int n_dofs_per_cell = cell_matrix.n(); - - Assert(cell_matrix.m() == n_dofs_per_cell, - ExcMessage( - "The provided cell mass matrix must be a square matrix.")); - AssertThrow( - n_cells <= 10, - ExcMessage( - "create_1D_discretization_matrix() returns a full matrix and is not meant to be used with a larger number of cells. ")); - Assert(n_cells > 0, - ExcMessage("You are trying to get a mass matrix of zero cells.")); - Assert(overlap < n_dofs_per_cell, - ExcMessage("The overlap must be smaller than the number of dofs.")); - - unsigned int n_total_dofs = - n_cells * n_dofs_per_cell - overlap * (n_cells - 1); - - if (!include_endpoints.first) - n_total_dofs -= 1; - if (!include_endpoints.second) - n_total_dofs -= 1; - - FullMatrix result_matrix(n_total_dofs, n_total_dofs); - result_matrix = 0; - - for (unsigned int cell = 0; cell < n_cells; ++cell) - { - const unsigned int dof_shift = - cell * overlap + !include_endpoints.first; - - const unsigned int start_dof = - (cell == 0 && !include_endpoints.first) ? 1 : 0; - - const unsigned int end_dof = - (cell == n_cells - 1 && !include_endpoints.second) ? - n_dofs_per_cell - 1 : - n_dofs_per_cell; - for (unsigned int i = start_dof; i < end_dof; ++i) - for (unsigned int j = start_dof; j < end_dof; ++j) - { - result_matrix(i + cell * n_dofs_per_cell - dof_shift, - j + cell * n_dofs_per_cell - dof_shift) += - cell_matrix(i, j); - } - } - return result_matrix; - } - - } // namespace MatrixCreator DEAL_II_NAMESPACE_CLOSE diff --git a/source/numerics/matrix_creator.cc b/source/numerics/matrix_creator.cc index 6e5187f529..425abd8eac 100644 --- a/source/numerics/matrix_creator.cc +++ b/source/numerics/matrix_creator.cc @@ -24,4 +24,149 @@ DEAL_II_NAMESPACE_OPEN #endif #include "numerics/matrix_creator.inst" +#if defined(SPLIT_INSTANTIATIONS_INDEX) && SPLIT_INSTANTIATIONS_INDEX == 0 +namespace MatrixCreator +{ + + FullMatrix + create_1d_cell_mass_matrix(const FiniteElement<1> &fe, + const double &h, + const std::pair include_endpoints, + std::vector numbering) + { + if (dynamic_cast *>(&fe) == nullptr && + numbering.size() == 0) + { + Assert( + include_endpoints.first == true && include_endpoints.second == true, + ExcMessage( + "You tried to genereate a 1D mass matrix with excluding boundary " + "dofs for a non-DGQ element without providing a numbering.")); + } + + + if (numbering.size() == 0) + { + numbering.resize(fe.dofs_per_cell); + std::iota(numbering.begin(), numbering.end(), 0); + } + const unsigned int degree = fe.degree; + const unsigned int n_dofs_per_cell = fe.dofs_per_cell; + const double &JxW = h; + QGauss<1> quadrature(degree + 1); + + FullMatrix cell_mass_matrix(n_dofs_per_cell, n_dofs_per_cell); + cell_mass_matrix = 0; + + unsigned int start_dof = include_endpoints.first ? 0 : 1; + unsigned int end_dof = + include_endpoints.second ? n_dofs_per_cell : n_dofs_per_cell - 1; + const unsigned int shift = include_endpoints.first ? 0 : 1; + + for (unsigned int i = start_dof; i < end_dof; ++i) + for (unsigned int j = start_dof; j < end_dof; ++j) + for (unsigned int q = 0; q < quadrature.size(); ++q) + cell_mass_matrix(i - shift, j - shift) += + (fe.shape_value(numbering[i], quadrature.point(q)) * + fe.shape_value(numbering[j], quadrature.point(q))) * + JxW * quadrature.weight(q); + + return cell_mass_matrix; + } + + + FullMatrix + create_1d_cell_derivative_matrix( + const FiniteElement<1> &fe, + const double &h, + const std::pair include_endpoints, + std::vector numbering) + { + if (numbering.size() == 0) + { + numbering.resize(fe.dofs_per_cell); + std::iota(numbering.begin(), numbering.end(), 0); + } + + const unsigned int degree = fe.degree; + const unsigned int n_dofs_per_cell = fe.dofs_per_cell; + const double &JxW = h; + QGauss<1> quadrature(degree + 1); + + FullMatrix cell_matrix(n_dofs_per_cell, n_dofs_per_cell); + cell_matrix = 0; + + unsigned int start_dof = include_endpoints.first ? 0 : 1; + unsigned int end_dof = + include_endpoints.second ? n_dofs_per_cell : n_dofs_per_cell - 1; + const unsigned int shift = include_endpoints.first ? 0 : 1; + + for (unsigned int i = start_dof; i < end_dof; ++i) + for (unsigned int j = start_dof; j < end_dof; ++j) + for (unsigned int q = 0; q < quadrature.size(); ++q) + cell_matrix(i - shift, j - shift) += + (fe.shape_grad(numbering[i], quadrature.point(q)) / h * + fe.shape_grad(numbering[j], quadrature.point(q))) / + h * JxW * quadrature.weight(q); + + return cell_matrix; + } + + FullMatrix + create_1D_discretization_matrix(FullMatrix &cell_matrix, + const unsigned int &n_cells, + const unsigned int &overlap, + const std::pair include_endpoints) + { + const unsigned int n_dofs_per_cell = cell_matrix.n(); + + Assert(cell_matrix.m() == n_dofs_per_cell, + ExcMessage( + "The provided cell mass matrix must be a square matrix.")); + AssertThrow( + n_cells <= 10, + ExcMessage( + "create_1D_discretization_matrix() returns a full matrix and is not meant to be used with a larger number of cells. ")); + Assert(n_cells > 0, + ExcMessage("You are trying to get a mass matrix of zero cells.")); + Assert(overlap < n_dofs_per_cell, + ExcMessage("The overlap must be smaller than the number of dofs.")); + + unsigned int n_total_dofs = + n_cells * n_dofs_per_cell - overlap * (n_cells - 1); + + if (!include_endpoints.first) + n_total_dofs -= 1; + if (!include_endpoints.second) + n_total_dofs -= 1; + + FullMatrix result_matrix(n_total_dofs, n_total_dofs); + result_matrix = 0; + + for (unsigned int cell = 0; cell < n_cells; ++cell) + { + const unsigned int dof_shift = + cell * overlap + !include_endpoints.first; + + const unsigned int start_dof = + (cell == 0 && !include_endpoints.first) ? 1 : 0; + + const unsigned int end_dof = + (cell == n_cells - 1 && !include_endpoints.second) ? + n_dofs_per_cell - 1 : + n_dofs_per_cell; + for (unsigned int i = start_dof; i < end_dof; ++i) + for (unsigned int j = start_dof; j < end_dof; ++j) + { + result_matrix(i + cell * n_dofs_per_cell - dof_shift, + j + cell * n_dofs_per_cell - dof_shift) += + cell_matrix(i, j); + } + } + return result_matrix; + } +} // namespace MatrixCreator + +#endif // SPLIT_INSTANTIATIONS_INDEX == 0 + DEAL_II_NAMESPACE_CLOSE