diff --git a/examples/poisson.cc b/examples/poisson.cc index 8266ec63..f39f3da1 100644 --- a/examples/poisson.cc +++ b/examples/poisson.cc @@ -11,7 +11,6 @@ // ----------------------------------------------------------------------------- -#include #include #include @@ -27,6 +26,7 @@ #include #include +#include #include #include @@ -409,8 +409,8 @@ class Poisson Triangulation tria; #ifdef HEX - MappingQ1 mapping; - FE_DGQ dg_fe; + MappingQ1 mapping; + FE_AggloDGP dg_fe; #else MappingFE mapping; FE_SimplexDGP dg_fe; diff --git a/include/fe_agglodgp.h b/include/fe_agglodgp.h new file mode 100644 index 00000000..97f12e8e --- /dev/null +++ b/include/fe_agglodgp.h @@ -0,0 +1,477 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2002 - 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +#ifndef dealii_fe_agglodgp_h +#define dealii_fe_agglodgp_h + +#include + +#include + +#include + +DEAL_II_NAMESPACE_OPEN + +/** + * @addtogroup fe + * @{ + */ + +/** + * Discontinuous finite elements based on Legendre polynomials. + * + * This finite element implements complete polynomial spaces, that is, + * dim-dimensional polynomials of degree p. For example, in 2d the element + * FE_AggloDGP(1) would represent the span of the functions $\{1,\hat x,\hat + * y\}$, which is in contrast to the element FE_DGQ(1) that is formed by the + * span of + * $\{1,\hat x,\hat y,\hat x\hat y\}$. Since the DGP space has only three + * unknowns for each quadrilateral, it is immediately clear that this element + * can not be continuous. + * + * The basis functions used in this element for the space described above are + * chosen to form a Legendre basis on the unit square, i.e., in particular + * they are $L_2$-orthogonal and normalized on the reference cell (but not + * necessarily on the real cell). As a consequence, the first basis function + * of this element is always the function that is constant and equal to one, + * regardless of the polynomial degree of the element. In addition, as a + * result of the orthogonality of the basis functions, the @ref GlossMassMatrix "mass matrix" is + * diagonal if the grid cells are parallelograms. Note that this is in + * contrast to the FE_DGPMonomial class that actually uses the monomial basis + * listed above as basis functions, without transformation from reference to + * real cell. + * + * The shape functions are defined in the class PolynomialSpace. The + * polynomials used inside PolynomialSpace are Polynomials::Legendre up to + * degree p given in FE_AggloDGP. For the ordering of the basis + * functions, refer to PolynomialSpace, remembering that the Legendre + * polynomials are ordered by ascending degree. + * + * @note This element is not defined by finding shape functions within the + * given function space that interpolate a particular set of points. + * Consequently, there are no support points to which a given function could + * be interpolated; finding a finite element function that approximates a + * given function is therefore only possible through projection, rather than + * interpolation. Secondly, the shape functions of this element do not jointly + * add up to one. As a consequence of this, adding or subtracting a constant + * value -- such as one would do to make a function have mean value zero -- + * can not be done by simply subtracting the constant value from each degree + * of freedom. Rather, one needs to use the fact that the first basis function + * is constant equal to one and simply subtract the constant from the value of + * the degree of freedom corresponding to this first shape function on each + * cell. + * + * + * @note This class is only partially implemented for the codimension one case + * (spacedim != dim ), since no passage of information between meshes + * of different refinement level is possible because the embedding and + * projection matrices are not computed in the class constructor. + * + *

Transformation properties

+ * + * It is worth noting that under a (bi-, tri-)linear mapping, the space + * described by this element does not contain $P(k)$, even if we use a basis + * of polynomials of degree $k$. Consequently, for example, on meshes with + * non-affine cells, a linear function can not be exactly represented by + * elements of type FE_AggloDGP(1) or FE_DGPMonomial(1). + * + * This can be understood by the following 2-d example: consider the cell with + * vertices at $(0,0),(1,0),(0,1),(s,s)$: + * @image html dgp_doesnt_contain_p.png + * + * For this cell, a bilinear transformation $F$ produces the relations $x=\hat + * x+\hat x\hat y$ and $y=\hat y+\hat x\hat y$ that correlate reference + * coordinates $\hat x,\hat y$ and coordinates in real space $x,y$. Under this + * mapping, the constant function is clearly mapped onto itself, but the two + * other shape functions of the $P_1$ space, namely $\phi_1(\hat x,\hat + * y)=\hat x$ and $\phi_2(\hat x,\hat y)=\hat y$ are mapped onto + * $\phi_1(x,y)=\frac{x-t}{t(s-1)},\phi_2(x,y)=t$ where + * $t=\frac{y}{s-x+sx+y-sy}$. + * + * For the simple case that $s=1$, i.e. if the real cell is the unit square, + * the expressions can be simplified to $t=y$ and + * $\phi_1(x,y)=x,\phi_2(x,y)=y$. However, for all other cases, the functions + * $\phi_1(x,y),\phi_2(x,y)$ are not linear any more, and neither is any + * linear combination of them. Consequently, the linear functions are not + * within the range of the mapped $P_1$ polynomials. + * + *

Visualization of shape functions

In 2d, the shape functions of + * this element look as follows. + * + *

$P_0$ element

+ * + * + * + * + * + *
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P1/P1_DGP_shape0000.png + *
$P_0$ element, + * shape function 0
+ * + *

$P_1$ element

+ * + * + * + * + * + * + * + * + * + * + * + *
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P1/P1_DGP_shape0000.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P1/P1_DGP_shape0001.png + *
$P_1$ element, shape function 0 $P_1$ element, shape function 1
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P1/P1_DGP_shape0002.png + *
$P_1$ element, + * shape function 2
+ * + * + *

$P_2$ element

+ * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + *
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P2/P2_DGP_shape0000.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P2/P2_DGP_shape0001.png + *
$P_2$ element, shape function 0 $P_2$ element, shape function 1
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P2/P2_DGP_shape0002.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P2/P2_DGP_shape0003.png + *
$P_2$ element, shape function 2 $P_2$ element, shape function 3
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P2/P2_DGP_shape0004.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P2/P2_DGP_shape0005.png + *
$P_2$ element, shape function 4 $P_2$ element, shape function 5
+ * + * + *

$P_3$ element

+ * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + *
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P3/P3_DGP_shape0000.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P3/P3_DGP_shape0001.png + *
$P_3$ element, shape function 0 $P_3$ element, shape function 1
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P3/P3_DGP_shape0002.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P3/P3_DGP_shape0003.png + *
$P_3$ element, shape function 2 $P_3$ element, shape function 3
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P3/P3_DGP_shape0004.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P3/P3_DGP_shape0005.png + *
$P_3$ element, shape function 4 $P_3$ element, shape function 5
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P3/P3_DGP_shape0006.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P3/P3_DGP_shape0007.png + *
$P_3$ element, shape function 6 $P_3$ element, shape function 7
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P3/P3_DGP_shape0008.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P3/P3_DGP_shape0009.png + *
$P_3$ element, shape function 8 $P_3$ element, shape function 9
+ * + * + *

$P_4$ element

+ * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + * + *
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0000.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0001.png + *
$P_4$ element, shape function 0 $P_4$ element, shape function 1
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0002.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0003.png + *
$P_4$ element, shape function 2 $P_4$ element, shape function 3
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0004.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0005.png + *
$P_4$ element, shape function 4 $P_4$ element, shape function 5
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0006.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0007.png + *
$P_4$ element, shape function 6 $P_4$ element, shape function 7
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0008.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0009.png + *
$P_4$ element, shape function 8 $P_4$ element, shape function 9
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0010.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0011.png + *
$P_4$ element, shape function 10 $P_4$ element, shape function 11
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0012.png + * + * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0013.png + *
$P_4$ element, shape function 12 $P_4$ element, shape function 13
+ * @image html http://www.dealii.org/images/shape-functions/DGP/P4/P4_DGP_shape0014.png + *
$P_4$ element, + * shape function 14
+ */ +template +class FE_AggloDGP : public FE_Poly +{ +public: + /** + * Constructor for tensor product polynomials of degree @p p. + */ + FE_AggloDGP(const unsigned int p); + + /** + * Return a string that uniquely identifies a finite element. This class + * returns FE_AggloDGP(degree), with @p dim and @p degree replaced + * by appropriate values. + */ + virtual std::string + get_name() const override; + + /** + * @name Functions to support hp + * @{ + */ + + /** + * If, on a vertex, several finite elements are active, the hp-code first + * assigns the degrees of freedom of each of these FEs different global + * indices. It then calls this function to find out which of them should get + * identical values, and consequently can receive the same global DoF index. + * This function therefore returns a list of identities between DoFs of the + * present finite element object with the DoFs of @p fe_other, which is a + * reference to a finite element object representing one of the other finite + * elements active on this particular vertex. The function computes which of + * the degrees of freedom of the two finite element objects are equivalent, + * both numbered between zero and the corresponding value of + * n_dofs_per_vertex() of the two finite elements. The first index of each + * pair denotes one of the vertex dofs of the present element, whereas the + * second is the corresponding index of the other finite element. + * + * This being a discontinuous element, the set of such constraints is of + * course empty. + */ + virtual std::vector> + hp_vertex_dof_identities( + const FiniteElement &fe_other) const override; + + /** + * Same as hp_vertex_dof_indices(), except that the function treats degrees + * of freedom on lines. + * + * This being a discontinuous element, the set of such constraints is of + * course empty. + */ + virtual std::vector> + hp_line_dof_identities( + const FiniteElement &fe_other) const override; + + /** + * Same as hp_vertex_dof_indices(), except that the function treats degrees + * of freedom on quads. + * + * This being a discontinuous element, the set of such constraints is of + * course empty. + */ + virtual std::vector> + hp_quad_dof_identities(const FiniteElement &fe_other, + const unsigned int face_no = 0) const override; + + /** + * Return whether this element implements its hanging node constraints in + * the new way, which has to be used to make elements "hp-compatible". + * + * For the FE_AggloDGP class the result is always true (independent of the + * degree of the element), as it has no hanging nodes (being a discontinuous + * element). + */ + virtual bool + hp_constraints_are_implemented() const override; + + /** + * @copydoc FiniteElement::compare_for_domination() + */ + virtual FiniteElementDomination::Domination + compare_for_domination(const FiniteElement &fe_other, + const unsigned int codim = 0) const override final; + + /** + * @} + */ + + /** + * Return the matrix interpolating from a face of one element to the face + * of the neighboring element. The size of the matrix is then + * source.dofs_per_face times this->dofs_per_face. + * + * Derived elements will have to implement this function. They may only + * provide interpolation matrices for certain source finite elements, for + * example those from the same family. If they don't implement interpolation + * from a given element, then they must throw an exception of type + * FiniteElement::ExcInterpolationNotImplemented. + */ + virtual void + get_face_interpolation_matrix(const FiniteElement &source, + FullMatrix &matrix, + const unsigned int face_no = 0) const override; + + /** + * Return the matrix interpolating from a face of one element to the face + * of the neighboring element. The size of the matrix is then + * source.dofs_per_face times this->dofs_per_face. + * + * Derived elements will have to implement this function. They may only + * provide interpolation matrices for certain source finite elements, for + * example those from the same family. If they don't implement interpolation + * from a given element, then they must throw an exception of type + * FiniteElement::ExcInterpolationNotImplemented. + */ + virtual void + get_subface_interpolation_matrix( + const FiniteElement &source, + const unsigned int subface, + FullMatrix &matrix, + const unsigned int face_no = 0) const override; + + /** + * This function returns @p true, if the shape function @p shape_index has + * non-zero function values somewhere on the face @p face_index. + */ + virtual bool + has_support_on_face(const unsigned int shape_index, + const unsigned int face_index) const override; + + /** + * Determine an estimate for the memory consumption (in bytes) of this + * object. + * + * This function is made virtual, since finite element objects are usually + * accessed through pointers to their base class, rather than the class + * itself. + */ + virtual std::size_t + memory_consumption() const override; + + + /** + * Return a list of constant modes of the element. For this element, the + * first entry is true, all other are false. + */ + virtual std::pair, std::vector> + get_constant_modes() const override; + + virtual std::unique_ptr> + clone() const override; + +private: + /** + * Only for internal use. Its full name is @p get_dofs_per_object_vector + * function and it creates the @p dofs_per_object vector that is needed + * within the constructor to be passed to the constructor of @p + * FiniteElementData. + */ + static std::vector + get_dpo_vector(const unsigned int degree); +}; + +/** @} */ + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/source/CMakeLists.txt b/source/CMakeLists.txt index cac6a7e4..7e023d72 100644 --- a/source/CMakeLists.txt +++ b/source/CMakeLists.txt @@ -3,7 +3,7 @@ SET(_d2_build_types "Release;Debug") SET(Release_postfix "") SET(Debug_postfix ".g") -set(_files agglomeration_handler.cc multigrid_amg.cc mapping_box.cc) +set(_files agglomeration_handler.cc multigrid_amg.cc mapping_box.cc fe_agglodgp.cc) FOREACH(_build_type ${_d2_build_types}) # Postfix to use everywhere diff --git a/source/agglomeration_handler.cc b/source/agglomeration_handler.cc index cfa56b01..3c19892b 100644 --- a/source/agglomeration_handler.cc +++ b/source/agglomeration_handler.cc @@ -15,6 +15,7 @@ #include #include +#include template AgglomerationHandler::AgglomerationHandler( @@ -236,8 +237,8 @@ AgglomerationHandler::distribute_agglomerated_dofs( { if (dynamic_cast *>(&fe_space)) fe = std::make_unique>(fe_space.degree); - else if (dynamic_cast *>(&fe_space)) - fe = std::make_unique>(fe_space.degree); + else if (dynamic_cast *>(&fe_space)) + fe = std::make_unique>(fe_space.degree); else if (dynamic_cast *>(&fe_space)) fe = std::make_unique>(fe_space.degree); else diff --git a/source/fe_agglodgp.cc b/source/fe_agglodgp.cc new file mode 100644 index 00000000..6ebd75a0 --- /dev/null +++ b/source/fe_agglodgp.cc @@ -0,0 +1,298 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2002 - 2024 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + + +#include +#include + +#include + +#include +#include + + +DEAL_II_NAMESPACE_OPEN + +template +FE_AggloDGP::FE_AggloDGP(const unsigned int degree) + : FE_Poly( + PolynomialSpace( + Polynomials::Legendre::generate_complete_basis(degree)), + FiniteElementData(get_dpo_vector(degree), + 1, + degree, + FiniteElementData::L2), + std::vector( + FiniteElementData(get_dpo_vector(degree), 1, degree) + .n_dofs_per_cell(), + true), + std::vector( + FiniteElementData(get_dpo_vector(degree), 1, degree) + .n_dofs_per_cell(), + ComponentMask(std::vector(1, true)))) +{ + // Reinit the vectors of restriction and prolongation matrices to the right + // sizes + this->reinit_restriction_and_prolongation_matrices(); + // Fill prolongation matrices with embedding operators + if (dim == spacedim) + { + FETools::compute_embedding_matrices(*this, this->prolongation); + // Fill restriction matrices with L2-projection + FETools::compute_projection_matrices(*this, this->restriction); + } +} + + +template +std::string +FE_AggloDGP::get_name() const +{ + // note that the FETools::get_fe_by_name function depends on the + // particular format of the string this function returns, so they have to be + // kept in sync + + std::ostringstream namebuf; + namebuf << "FE_AggloDGP<" << Utilities::dim_string(dim, spacedim) << ">(" + << this->degree << ")"; + + return namebuf.str(); +} + + + +template +std::unique_ptr> +FE_AggloDGP::clone() const +{ + return std::make_unique>(*this); +} + + + +//--------------------------------------------------------------------------- +// Auxiliary functions +//--------------------------------------------------------------------------- + + +template +std::vector +FE_AggloDGP::get_dpo_vector(const unsigned int deg) +{ + std::vector dpo(dim + 1, 0U); + dpo[dim] = deg + 1; + for (unsigned int i = 1; i < dim; ++i) + { + dpo[dim] *= deg + 1 + i; + dpo[dim] /= i + 1; + } + return dpo; +} + + + +template +void +FE_AggloDGP::get_face_interpolation_matrix( + const FiniteElement &x_source_fe, + FullMatrix &interpolation_matrix, + const unsigned int) const +{ + // this is only implemented, if the source FE is also a DGP element. in that + // case, both elements have no dofs on their faces and the face + // interpolation matrix is necessarily empty -- i.e. there isn't much we + // need to do here. + (void)interpolation_matrix; + using FE = FiniteElement; + using FEDGP = FE_AggloDGP; + AssertThrow((x_source_fe.get_name().find("FE_AggloDGP<") == 0) || + (dynamic_cast(&x_source_fe) != nullptr), + typename FE::ExcInterpolationNotImplemented()); + + Assert(interpolation_matrix.m() == 0, + ExcDimensionMismatch(interpolation_matrix.m(), 0)); + Assert(interpolation_matrix.n() == 0, + ExcDimensionMismatch(interpolation_matrix.n(), 0)); +} + + + +template +void +FE_AggloDGP::get_subface_interpolation_matrix( + const FiniteElement &x_source_fe, + const unsigned int, + FullMatrix &interpolation_matrix, + const unsigned int) const +{ + // this is only implemented, if the source FE is also a DGP element. in that + // case, both elements have no dofs on their faces and the face + // interpolation matrix is necessarily empty -- i.e. there isn't much we + // need to do here. + (void)interpolation_matrix; + using FE = FiniteElement; + using FEDGP = FE_AggloDGP; + AssertThrow((x_source_fe.get_name().find("FE_AggloDGP<") == 0) || + (dynamic_cast(&x_source_fe) != nullptr), + typename FE::ExcInterpolationNotImplemented()); + + Assert(interpolation_matrix.m() == 0, + ExcDimensionMismatch(interpolation_matrix.m(), 0)); + Assert(interpolation_matrix.n() == 0, + ExcDimensionMismatch(interpolation_matrix.n(), 0)); +} + + + +template +bool +FE_AggloDGP::hp_constraints_are_implemented() const +{ + return true; +} + + + +template +std::vector> +FE_AggloDGP::hp_vertex_dof_identities( + const FiniteElement &fe_other) const +{ + (void)fe_other; + return std::vector>(); +} + + + +template +std::vector> +FE_AggloDGP::hp_line_dof_identities( + const FiniteElement &fe_other) const +{ + // there are no such constraints for DGP elements at all + if (dynamic_cast *>(&fe_other) != nullptr) + return std::vector>(); + else + { + DEAL_II_NOT_IMPLEMENTED(); + return std::vector>(); + } +} + + + +template +std::vector> +FE_AggloDGP::hp_quad_dof_identities( + const FiniteElement &fe_other, + const unsigned int) const +{ + // there are no such constraints for DGP elements at all + if (dynamic_cast *>(&fe_other) != nullptr) + return std::vector>(); + else + { + DEAL_II_NOT_IMPLEMENTED(); + return std::vector>(); + } +} + + + +template +FiniteElementDomination::Domination +FE_AggloDGP::compare_for_domination( + const FiniteElement &fe_other, + const unsigned int codim) const +{ + Assert(codim <= dim, ExcImpossibleInDim(dim)); + + // vertex/line/face domination + // --------------------------- + if (codim > 0) + // this is a discontinuous element, so by definition there will + // be no constraints wherever this element comes together with + // any other kind of element + return FiniteElementDomination::no_requirements; + + // cell domination + // --------------- + if (const FE_AggloDGP *fe_dgp_other = + dynamic_cast *>(&fe_other)) + { + if (this->degree < fe_dgp_other->degree) + return FiniteElementDomination::this_element_dominates; + else if (this->degree == fe_dgp_other->degree) + return FiniteElementDomination::either_element_can_dominate; + else + return FiniteElementDomination::other_element_dominates; + } + else if (const FE_Nothing *fe_nothing = + dynamic_cast *>(&fe_other)) + { + if (fe_nothing->is_dominating()) + return FiniteElementDomination::other_element_dominates; + else + // the FE_Nothing has no degrees of freedom and it is typically used + // in a context where we don't require any continuity along the + // interface + return FiniteElementDomination::no_requirements; + } + + DEAL_II_NOT_IMPLEMENTED(); + return FiniteElementDomination::neither_element_dominates; +} + + + +template +bool +FE_AggloDGP::has_support_on_face(const unsigned int, + const unsigned int) const +{ + // all shape functions have support on all faces + return true; +} + + + +template +std::pair, std::vector> +FE_AggloDGP::get_constant_modes() const +{ + Table<2, bool> constant_modes(1, this->n_dofs_per_cell()); + constant_modes(0, 0) = true; + return std::pair, std::vector>( + constant_modes, std::vector(1, 0)); +} + + + +template +std::size_t +FE_AggloDGP::memory_consumption() const +{ + DEAL_II_NOT_IMPLEMENTED(); + return 0; +} + + + +// explicit instantiations +template class FE_AggloDGP<1>; +template class FE_AggloDGP<2>; +template class FE_AggloDGP<3>; + + +DEAL_II_NAMESPACE_CLOSE diff --git a/test/polydeal/exact_solutions_dgp.cc b/test/polydeal/exact_solutions_dgp.cc index 11c82372..5b554eba 100644 --- a/test/polydeal/exact_solutions_dgp.cc +++ b/test/polydeal/exact_solutions_dgp.cc @@ -1,5 +1,3 @@ -#include - #include #include #include @@ -16,6 +14,7 @@ #include #include +#include #include #include @@ -241,7 +240,7 @@ class Poisson Triangulation tria; MappingQ mapping; - FE_DGP *dg_fe; + FE_AggloDGP *dg_fe; std::unique_ptr> ah; // no hanging node in DG discretization, we define an AffineConstraints // object @@ -276,12 +275,12 @@ Poisson::Poisson(const SolutionType &solution_type) { if (sol_type == SolutionType::LinearSolution) { - dg_fe = new FE_DGP{1}; + dg_fe = new FE_AggloDGP{1}; analytical_solution = new SolutionLinear(); } else { - dg_fe = new FE_DGP{2}; + dg_fe = new FE_AggloDGP{2}; analytical_solution = new SolutionQuadratic(); } }