From 6a6a25b047787d30b9dce4515c4a0fcd157c477b Mon Sep 17 00:00:00 2001 From: bangerth Date: Fri, 11 Aug 2006 15:14:28 +0000 Subject: [PATCH] Implement face and subface interpolation matrices for the FESystem element git-svn-id: https://svn.dealii.org/trunk@13665 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe_system.h | 54 ++++ deal.II/deal.II/source/fe/fe_system.cc | 256 ++++++++++++++- tests/bits/face_interpolation/cmp/generic | 108 +++++++ tests/bits/subface_interpolation/cmp/generic | 324 +++++++++++++++++++ 4 files changed, 740 insertions(+), 2 deletions(-) diff --git a/deal.II/deal.II/include/fe/fe_system.h b/deal.II/deal.II/include/fe/fe_system.h index f09bd41d9f..7a6636a8fd 100644 --- a/deal.II/deal.II/include/fe/fe_system.h +++ b/deal.II/deal.II/include/fe/fe_system.h @@ -409,6 +409,60 @@ class FESystem : public FiniteElement */ virtual bool hp_constraints_are_implemented () const; + /** + * Return the matrix + * interpolating from a face of + * of one element to the face of + * the neighboring element. + * The size of the matrix is + * then @p dofs_per_face times + * source.dofs_per_face. + * + * Base elements of this element 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, + * which will get propagated out from + * this element. + */ + virtual void + get_face_interpolation_matrix (const FiniteElement &source, + FullMatrix &matrix) const; + + + /** + * Return the matrix + * interpolating from a face of + * of one element to the subface of + * the neighboring element. + * The size of the matrix is + * then @p dofs_per_face times + * source.dofs_per_face. + * + * Base elements of this element 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, + * which will get propagated out from + * this element. + */ + virtual void + get_subface_interpolation_matrix (const FiniteElement &source, + const unsigned int subface, + FullMatrix &matrix) const; + /** * If, on a vertex, several * finite elements are active, diff --git a/deal.II/deal.II/source/fe/fe_system.cc b/deal.II/deal.II/source/fe/fe_system.cc index 6aeb731c2f..cd4f6982d3 100644 --- a/deal.II/deal.II/source/fe/fe_system.cc +++ b/deal.II/deal.II/source/fe/fe_system.cc @@ -1851,8 +1851,6 @@ bool FESystem:: hp_constraints_are_implemented () const { - return false; - for (unsigned int b=0; b +void +FESystem:: +get_face_interpolation_matrix (const FiniteElement &x_source_fe, + FullMatrix &interpolation_matrix) const +{ + AssertThrow ((x_source_fe.get_name().find ("FE_System<") == 0) + || + (dynamic_cast*>(&x_source_fe) != 0), + typename FiniteElement:: + ExcInterpolationNotImplemented()); + + Assert (interpolation_matrix.m() == this->dofs_per_face, + ExcDimensionMismatch (interpolation_matrix.m(), + this->dofs_per_face)); + Assert (interpolation_matrix.n() == x_source_fe.dofs_per_face, + ExcDimensionMismatch (interpolation_matrix.m(), + x_source_fe.dofs_per_face)); + + // since dofs for each base are + // independent, we only have to stack + // things up from base element to base + // element + // + // the problem is that we have to work with + // two FEs (this and fe_other). only deal + // with the case that both are FESystems + // and that they both have the same number + // of bases (counting multiplicity) each of + // which match in their number of + // components. this covers + // FESystem(FE_Q(p),1,FE_Q(q),2) vs + // FESystem(FE_Q(r),2,FE_Q(s),1), but not + // FESystem(FE_Q(p),1,FE_Q(q),2) vs + // FESystem(FESystem(FE_Q(r),2),1,FE_Q(s),1) + const FESystem *fe_other_system + = dynamic_cast*>(&x_source_fe); + + // clear matrix, since we will not get to + // set all elements + interpolation_matrix = 0; + + // loop over all the base elements of + // this and the other element, counting + // their multiplicities + unsigned int base_index = 0, + base_index_other = 0; + unsigned int multiplicity = 0, + multiplicity_other = 0; + + FullMatrix base_to_base_interpolation; + + while (true) + { + const FiniteElement + &base = base_element(base_index), + &base_other = fe_other_system->base_element(base_index_other); + + Assert (base.n_components() == base_other.n_components(), + ExcNotImplemented()); + + // get the interpolation from the bases + base_to_base_interpolation.reinit (base.dofs_per_face, + base_other.dofs_per_face); + base.get_face_interpolation_matrix (base_other, + base_to_base_interpolation); + + // now translate entries. we'd like to + // have something like + // face_base_to_system_index, but that + // doesn't exist. rather, all we have + // is the reverse. well, use that then + for (unsigned int i=0; idofs_per_face; ++i) + if (this->face_system_to_base_index(i).first + == + std::make_pair (base_index, multiplicity)) + for (unsigned int j=0; jdofs_per_face; ++j) + if (fe_other_system->face_system_to_base_index(j).first + == + std::make_pair (base_index_other, multiplicity_other)) + interpolation_matrix(i, j) + = base_to_base_interpolation(this->face_system_to_base_index(i).second, + fe_other_system->face_system_to_base_index(j).second); + + // advance to the next base element + // for this and the other + // fe_system; see if we can simply + // advance the multiplicity by one, + // or if have to move on to the + // next base element + ++multiplicity; + if (multiplicity == element_multiplicity(base_index)) + { + multiplicity = 0; + ++base_index; + } + ++multiplicity_other; + if (multiplicity_other == + fe_other_system->element_multiplicity(base_index_other)) + { + multiplicity_other = 0; + ++base_index_other; + } + + // see if we have reached the end + // of the present element. if so, + // we should have reached the end + // of the other one as well + if (base_index == n_base_elements()) + { + Assert (base_index_other == fe_other_system->n_base_elements(), + ExcInternalError()); + break; + } + + // if we haven't reached the end of + // this element, we shouldn't have + // reached the end of the other one + // either + Assert (base_index_other != fe_other_system->n_base_elements(), + ExcInternalError()); + } +} + + + +template +void +FESystem:: +get_subface_interpolation_matrix (const FiniteElement &x_source_fe, + const unsigned int subface, + FullMatrix &interpolation_matrix) const +{ + AssertThrow ((x_source_fe.get_name().find ("FE_System<") == 0) + || + (dynamic_cast*>(&x_source_fe) != 0), + typename FiniteElement:: + ExcInterpolationNotImplemented()); + + Assert (interpolation_matrix.m() == this->dofs_per_face, + ExcDimensionMismatch (interpolation_matrix.m(), + this->dofs_per_face)); + Assert (interpolation_matrix.n() == x_source_fe.dofs_per_face, + ExcDimensionMismatch (interpolation_matrix.m(), + x_source_fe.dofs_per_face)); + + // since dofs for each base are + // independent, we only have to stack + // things up from base element to base + // element + // + // the problem is that we have to work with + // two FEs (this and fe_other). only deal + // with the case that both are FESystems + // and that they both have the same number + // of bases (counting multiplicity) each of + // which match in their number of + // components. this covers + // FESystem(FE_Q(p),1,FE_Q(q),2) vs + // FESystem(FE_Q(r),2,FE_Q(s),1), but not + // FESystem(FE_Q(p),1,FE_Q(q),2) vs + // FESystem(FESystem(FE_Q(r),2),1,FE_Q(s),1) + const FESystem *fe_other_system + = dynamic_cast*>(&x_source_fe); + + // clear matrix, since we will not get to + // set all elements + interpolation_matrix = 0; + + // loop over all the base elements of + // this and the other element, counting + // their multiplicities + unsigned int base_index = 0, + base_index_other = 0; + unsigned int multiplicity = 0, + multiplicity_other = 0; + + FullMatrix base_to_base_interpolation; + + while (true) + { + const FiniteElement + &base = base_element(base_index), + &base_other = fe_other_system->base_element(base_index_other); + + Assert (base.n_components() == base_other.n_components(), + ExcNotImplemented()); + + // get the interpolation from the bases + base_to_base_interpolation.reinit (base.dofs_per_face, + base_other.dofs_per_face); + base.get_subface_interpolation_matrix (base_other, + subface, + base_to_base_interpolation); + + // now translate entries. we'd like to + // have something like + // face_base_to_system_index, but that + // doesn't exist. rather, all we have + // is the reverse. well, use that then + for (unsigned int i=0; idofs_per_face; ++i) + if (this->face_system_to_base_index(i).first + == + std::make_pair (base_index, multiplicity)) + for (unsigned int j=0; jdofs_per_face; ++j) + if (fe_other_system->face_system_to_base_index(j).first + == + std::make_pair (base_index_other, multiplicity_other)) + interpolation_matrix(i, j) + = base_to_base_interpolation(this->face_system_to_base_index(i).second, + fe_other_system->face_system_to_base_index(j).second); + + // advance to the next base element + // for this and the other + // fe_system; see if we can simply + // advance the multiplicity by one, + // or if have to move on to the + // next base element + ++multiplicity; + if (multiplicity == element_multiplicity(base_index)) + { + multiplicity = 0; + ++base_index; + } + ++multiplicity_other; + if (multiplicity_other == + fe_other_system->element_multiplicity(base_index_other)) + { + multiplicity_other = 0; + ++base_index_other; + } + + // see if we have reached the end + // of the present element. if so, + // we should have reached the end + // of the other one as well + if (base_index == n_base_elements()) + { + Assert (base_index_other == fe_other_system->n_base_elements(), + ExcInternalError()); + break; + } + + // if we haven't reached the end of + // this element, we shouldn't have + // reached the end of the other one + // either + Assert (base_index_other != fe_other_system->n_base_elements(), + ExcInternalError()); + } +} + + + template template std::vector > diff --git a/tests/bits/face_interpolation/cmp/generic b/tests/bits/face_interpolation/cmp/generic index 695c78f8ed..b682f4faa3 100644 --- a/tests/bits/face_interpolation/cmp/generic +++ b/tests/bits/face_interpolation/cmp/generic @@ -1130,14 +1130,122 @@ DEAL::Checking FE_DGQ<1>(3)2 against FE_DGQ<1>(2)2 in 1d: DEAL::Checking FE_DGP<1>(3)1 against FE_DGP<1>(1)1 in 1d: DEAL::Checking FE_DGP<1>(1)1 against FE_DGP<1>(3)1 in 1d: DEAL::Checking FE_Q<2>(1)3 against FE_Q<2>(2)3 in 2d: +DEAL::FESystem<2>[FE_Q<2>(1)^3] vs. FESystem<2>[FE_Q<2>(1)^3] +DEAL::1.0 1.0 +DEAL::2.4 +DEAL::1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 +DEAL::FESystem<2>[FE_Q<2>(2)^3] vs. FESystem<2>[FE_Q<2>(2)^3] +DEAL::1.0 1.0 +DEAL::3.0 +DEAL::1.0 0 1.0 0 1.0 0 1.0 0 1.0 1.0 1.0 0 1.0 0 1.0 0 1.0 0 +DEAL::FESystem<2>[FE_Q<2>(1)^3] vs. FESystem<2>[FE_Q<2>(2)^3] +DEAL::1.0 1.5 +DEAL::1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 DEAL::Checking FE_Q<2>(2)3 against FE_Q<2>(1)3 in 2d: +DEAL::FESystem<2>[FE_Q<2>(2)^3] vs. FESystem<2>[FE_Q<2>(2)^3] +DEAL::1.0 1.0 +DEAL::3.0 +DEAL::1.0 0 1.0 0 1.0 0 1.0 0 1.0 1.0 1.0 0 1.0 0 1.0 0 1.0 0 +DEAL::FESystem<2>[FE_Q<2>(1)^3] vs. FESystem<2>[FE_Q<2>(1)^3] +DEAL::1.0 1.0 +DEAL::2.4 +DEAL::1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 +DEAL::FESystem<2>[FE_Q<2>(1)^3] vs. FESystem<2>[FE_Q<2>(2)^3] +DEAL::1.0 1.5 +DEAL::1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 DEAL::Checking FE_DGQ<2>(2)2 against FE_DGQ<2>(3)2 in 2d: +DEAL::FESystem<2>[FE_DGQ<2>(2)^2] vs. FESystem<2>[FE_DGQ<2>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(3)^2] vs. FESystem<2>[FE_DGQ<2>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(2)^2] vs. FESystem<2>[FE_DGQ<2>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(3)^2] vs. FESystem<2>[FE_DGQ<2>(2)^2] +DEAL::(Empty matrix) DEAL::Checking FE_DGQ<2>(3)2 against FE_DGQ<2>(2)2 in 2d: +DEAL::FESystem<2>[FE_DGQ<2>(3)^2] vs. FESystem<2>[FE_DGQ<2>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(2)^2] vs. FESystem<2>[FE_DGQ<2>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(3)^2] vs. FESystem<2>[FE_DGQ<2>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(2)^2] vs. FESystem<2>[FE_DGQ<2>(3)^2] +DEAL::(Empty matrix) DEAL::Checking FE_DGP<2>(3)1 against FE_DGP<2>(1)1 in 2d: +DEAL::FESystem<2>[FE_DGP<2>(3)] vs. FESystem<2>[FE_DGP<2>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(1)] vs. FESystem<2>[FE_DGP<2>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(3)] vs. FESystem<2>[FE_DGP<2>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(1)] vs. FESystem<2>[FE_DGP<2>(3)] +DEAL::(Empty matrix) DEAL::Checking FE_DGP<2>(1)1 against FE_DGP<2>(3)1 in 2d: +DEAL::FESystem<2>[FE_DGP<2>(1)] vs. FESystem<2>[FE_DGP<2>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(3)] vs. FESystem<2>[FE_DGP<2>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(1)] vs. FESystem<2>[FE_DGP<2>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(3)] vs. FESystem<2>[FE_DGP<2>(1)] +DEAL::(Empty matrix) DEAL::Checking FE_Q<3>(1)3 against FE_Q<3>(2)3 in 3d: +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(1)^3] +DEAL::1.0 1.0 +DEAL::3.5 +DEAL::1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 +DEAL::FESystem<3>[FE_Q<3>(2)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.0 1.0 +DEAL::5.2 +DEAL::1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 1.0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.0 2.2 +DEAL::1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 DEAL::Checking FE_Q<3>(2)3 against FE_Q<3>(1)3 in 3d: +DEAL::FESystem<3>[FE_Q<3>(2)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.0 1.0 +DEAL::5.2 +DEAL::1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 1.0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(1)^3] +DEAL::1.0 1.0 +DEAL::3.5 +DEAL::1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.0 2.2 +DEAL::1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 1.0 0 DEAL::Checking FE_DGQ<3>(2)2 against FE_DGQ<3>(3)2 in 3d: +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) DEAL::Checking FE_DGQ<3>(3)2 against FE_DGQ<3>(2)2 in 3d: +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) DEAL::Checking FE_DGP<3>(3)1 against FE_DGP<3>(1)1 in 3d: +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) DEAL::Checking FE_DGP<3>(1)1 against FE_DGP<3>(3)1 in 3d: +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) diff --git a/tests/bits/subface_interpolation/cmp/generic b/tests/bits/subface_interpolation/cmp/generic index c5bdcf8d2b..29e0febc9e 100644 --- a/tests/bits/subface_interpolation/cmp/generic +++ b/tests/bits/subface_interpolation/cmp/generic @@ -3058,14 +3058,338 @@ DEAL::Checking FE_DGQ<1>(3)2 against FE_DGQ<1>(2)2 in 1d: DEAL::Checking FE_DGP<1>(3)1 against FE_DGP<1>(1)1 in 1d: DEAL::Checking FE_DGP<1>(1)1 against FE_DGP<1>(3)1 in 1d: DEAL::Checking FE_Q<2>(1)3 against FE_Q<2>(2)3 in 2d: +DEAL::FESystem<2>[FE_Q<2>(1)^3] vs. FESystem<2>[FE_Q<2>(1)^3] +DEAL::1.0 1.5 +DEAL::2.1 +DEAL::1.0 0 1.0 0.50 1.0 0 0.50 0 0.50 0 0.50 0 +DEAL::FESystem<2>[FE_Q<2>(2)^3] vs. FESystem<2>[FE_Q<2>(2)^3] +DEAL::1.2 1.8 +DEAL::2.9 +DEAL::1.0 0 1.0 0.38 1.0 0 0 0 0 0 0 0 0.75 0 0.75 0 0.75 0 +DEAL::FESystem<2>[FE_Q<2>(1)^3] vs. FESystem<2>[FE_Q<2>(2)^3] +DEAL::1.0 2.2 +DEAL::1.0 0 1.0 0.50 1.0 0 0.50 0 0.50 0 0.50 0 +DEAL::FESystem<2>[FE_Q<2>(1)^3] vs. FESystem<2>[FE_Q<2>(1)^3] +DEAL::1.0 1.5 +DEAL::2.1 +DEAL::0.50 0 0.50 0 0.50 0 1.0 0 1.0 0.50 1.0 0 +DEAL::FESystem<2>[FE_Q<2>(2)^3] vs. FESystem<2>[FE_Q<2>(2)^3] +DEAL::1.2 1.8 +DEAL::2.9 +DEAL::0 0 0 -0.12 0 0 1.0 0 1.0 1.0 1.0 0 0.75 0 0.75 1.0 0.75 0 +DEAL::FESystem<2>[FE_Q<2>(1)^3] vs. FESystem<2>[FE_Q<2>(2)^3] +DEAL::1.0 2.2 +DEAL::0.50 0 0.50 0 0.50 0 1.0 0 1.0 0.50 1.0 0 DEAL::Checking FE_Q<2>(2)3 against FE_Q<2>(1)3 in 2d: +DEAL::FESystem<2>[FE_Q<2>(2)^3] vs. FESystem<2>[FE_Q<2>(2)^3] +DEAL::1.2 1.8 +DEAL::2.9 +DEAL::1.0 0 1.0 0.38 1.0 0 0 0 0 0 0 0 0.75 0 0.75 0 0.75 0 +DEAL::FESystem<2>[FE_Q<2>(1)^3] vs. FESystem<2>[FE_Q<2>(1)^3] +DEAL::1.0 1.5 +DEAL::2.1 +DEAL::1.0 0 1.0 0.50 1.0 0 0.50 0 0.50 0 0.50 0 +DEAL::FESystem<2>[FE_Q<2>(1)^3] vs. FESystem<2>[FE_Q<2>(2)^3] +DEAL::1.0 2.2 +DEAL::1.0 0 1.0 0.50 1.0 0 0.50 0 0.50 0 0.50 0 +DEAL::FESystem<2>[FE_Q<2>(2)^3] vs. FESystem<2>[FE_Q<2>(2)^3] +DEAL::1.2 1.8 +DEAL::2.9 +DEAL::0 0 0 -0.12 0 0 1.0 0 1.0 1.0 1.0 0 0.75 0 0.75 1.0 0.75 0 +DEAL::FESystem<2>[FE_Q<2>(1)^3] vs. FESystem<2>[FE_Q<2>(1)^3] +DEAL::1.0 1.5 +DEAL::2.1 +DEAL::0.50 0 0.50 0 0.50 0 1.0 0 1.0 0.50 1.0 0 +DEAL::FESystem<2>[FE_Q<2>(1)^3] vs. FESystem<2>[FE_Q<2>(2)^3] +DEAL::1.0 2.2 +DEAL::0.50 0 0.50 0 0.50 0 1.0 0 1.0 0.50 1.0 0 DEAL::Checking FE_DGQ<2>(2)2 against FE_DGQ<2>(3)2 in 2d: +DEAL::FESystem<2>[FE_DGQ<2>(2)^2] vs. FESystem<2>[FE_DGQ<2>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(3)^2] vs. FESystem<2>[FE_DGQ<2>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(2)^2] vs. FESystem<2>[FE_DGQ<2>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(3)^2] vs. FESystem<2>[FE_DGQ<2>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(2)^2] vs. FESystem<2>[FE_DGQ<2>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(3)^2] vs. FESystem<2>[FE_DGQ<2>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(2)^2] vs. FESystem<2>[FE_DGQ<2>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(3)^2] vs. FESystem<2>[FE_DGQ<2>(2)^2] +DEAL::(Empty matrix) DEAL::Checking FE_DGQ<2>(3)2 against FE_DGQ<2>(2)2 in 2d: +DEAL::FESystem<2>[FE_DGQ<2>(3)^2] vs. FESystem<2>[FE_DGQ<2>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(2)^2] vs. FESystem<2>[FE_DGQ<2>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(3)^2] vs. FESystem<2>[FE_DGQ<2>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(2)^2] vs. FESystem<2>[FE_DGQ<2>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(3)^2] vs. FESystem<2>[FE_DGQ<2>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(2)^2] vs. FESystem<2>[FE_DGQ<2>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(3)^2] vs. FESystem<2>[FE_DGQ<2>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGQ<2>(2)^2] vs. FESystem<2>[FE_DGQ<2>(3)^2] +DEAL::(Empty matrix) DEAL::Checking FE_DGP<2>(3)1 against FE_DGP<2>(1)1 in 2d: +DEAL::FESystem<2>[FE_DGP<2>(3)] vs. FESystem<2>[FE_DGP<2>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(1)] vs. FESystem<2>[FE_DGP<2>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(3)] vs. FESystem<2>[FE_DGP<2>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(1)] vs. FESystem<2>[FE_DGP<2>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(3)] vs. FESystem<2>[FE_DGP<2>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(1)] vs. FESystem<2>[FE_DGP<2>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(3)] vs. FESystem<2>[FE_DGP<2>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(1)] vs. FESystem<2>[FE_DGP<2>(3)] +DEAL::(Empty matrix) DEAL::Checking FE_DGP<2>(1)1 against FE_DGP<2>(3)1 in 2d: +DEAL::FESystem<2>[FE_DGP<2>(1)] vs. FESystem<2>[FE_DGP<2>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(3)] vs. FESystem<2>[FE_DGP<2>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(1)] vs. FESystem<2>[FE_DGP<2>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(3)] vs. FESystem<2>[FE_DGP<2>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(1)] vs. FESystem<2>[FE_DGP<2>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(3)] vs. FESystem<2>[FE_DGP<2>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(1)] vs. FESystem<2>[FE_DGP<2>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<2>[FE_DGP<2>(3)] vs. FESystem<2>[FE_DGP<2>(1)] +DEAL::(Empty matrix) DEAL::Checking FE_Q<3>(1)3 against FE_Q<3>(2)3 in 3d: +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(1)^3] +DEAL::1.0 2.2 +DEAL::2.6 +DEAL::1.0 0 1.0 0.25 1.0 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 0.25 0 0.25 0 0.25 0 +DEAL::FESystem<3>[FE_Q<3>(2)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.6 3.1 +DEAL::4.7 +DEAL::1.0 0 1.0 0.14 1.0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.75 0 0.75 0.75 0.75 0 0 0 0 0 0 0 0.75 0 0.75 0 0.75 0 0 0 0 0 0 0 0.56 0 0.56 0 0.56 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.0 5.1 +DEAL::1.0 0 1.0 0.25 1.0 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 0.25 0 0.25 0 0.25 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(1)^3] +DEAL::1.0 2.2 +DEAL::2.6 +DEAL::0.50 0 0.50 0 0.50 0 1.0 0 1.0 0.25 1.0 0 0.25 0 0.25 0 0.25 0 0.50 0 0.50 0 0.50 0 +DEAL::FESystem<3>[FE_Q<3>(2)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.6 3.1 +DEAL::4.7 +DEAL::0 0 0 -0.047 0 0 1.0 0 1.0 0 1.0 0 0 0 0 0 0 0 0 0 0 -0.12 0 0 0 0 0 0 0 0 0.75 0 0.75 1.0 0.75 0 0.75 0 0.75 0 0.75 0 0 0 0 0 0 0 0.56 0 0.56 0 0.56 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.0 5.1 +DEAL::0.50 0 0.50 0 0.50 0 1.0 0 1.0 0.25 1.0 0 0.25 0 0.25 0 0.25 0 0.50 0 0.50 0 0.50 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(1)^3] +DEAL::1.0 2.2 +DEAL::2.6 +DEAL::0.50 0 0.50 0 0.50 0 0.25 0 0.25 0 0.25 0 1.0 0 1.0 0.25 1.0 0 0.50 0 0.50 0 0.50 0 +DEAL::FESystem<3>[FE_Q<3>(2)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.6 3.1 +DEAL::4.7 +DEAL::0 0 0 -0.047 0 0 0 0 0 0 0 0 1.0 0 1.0 0 1.0 0 0 0 0 0 0 0 0.75 0 0.75 0.75 0.75 0 0 0 0 0 0 0 0 0 0 0 0 0 0.75 0 0.75 0 0.75 0 0.56 0 0.56 0 0.56 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.0 5.1 +DEAL::0.50 0 0.50 0 0.50 0 0.25 0 0.25 0 0.25 0 1.0 0 1.0 0.25 1.0 0 0.50 0 0.50 0 0.50 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(1)^3] +DEAL::1.0 2.2 +DEAL::2.6 +DEAL::0.25 0 0.25 0 0.25 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 1.0 0 1.0 0.25 1.0 0 +DEAL::FESystem<3>[FE_Q<3>(2)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.6 3.1 +DEAL::4.7 +DEAL::0 0 0 0.016 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1.0 0 1.0 0.38 1.0 0 0 0 0 0 0 0 0.75 0 0.75 0 0.75 0 0 0 0 0 0 0 0.75 0 0.75 0 0.75 0 0.56 0 0.56 1.0 0.56 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.0 5.1 +DEAL::0.25 0 0.25 0 0.25 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 1.0 0 1.0 0.25 1.0 0 DEAL::Checking FE_Q<3>(2)3 against FE_Q<3>(1)3 in 3d: +DEAL::FESystem<3>[FE_Q<3>(2)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.6 3.1 +DEAL::4.7 +DEAL::1.0 0 1.0 0.14 1.0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.75 0 0.75 0.75 0.75 0 0 0 0 0 0 0 0.75 0 0.75 0 0.75 0 0 0 0 0 0 0 0.56 0 0.56 0 0.56 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(1)^3] +DEAL::1.0 2.2 +DEAL::2.6 +DEAL::1.0 0 1.0 0.25 1.0 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 0.25 0 0.25 0 0.25 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.0 5.1 +DEAL::1.0 0 1.0 0.25 1.0 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 0.25 0 0.25 0 0.25 0 +DEAL::FESystem<3>[FE_Q<3>(2)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.6 3.1 +DEAL::4.7 +DEAL::0 0 0 -0.047 0 0 1.0 0 1.0 0 1.0 0 0 0 0 0 0 0 0 0 0 -0.12 0 0 0 0 0 0 0 0 0.75 0 0.75 1.0 0.75 0 0.75 0 0.75 0 0.75 0 0 0 0 0 0 0 0.56 0 0.56 0 0.56 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(1)^3] +DEAL::1.0 2.2 +DEAL::2.6 +DEAL::0.50 0 0.50 0 0.50 0 1.0 0 1.0 0.25 1.0 0 0.25 0 0.25 0 0.25 0 0.50 0 0.50 0 0.50 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.0 5.1 +DEAL::0.50 0 0.50 0 0.50 0 1.0 0 1.0 0.25 1.0 0 0.25 0 0.25 0 0.25 0 0.50 0 0.50 0 0.50 0 +DEAL::FESystem<3>[FE_Q<3>(2)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.6 3.1 +DEAL::4.7 +DEAL::0 0 0 -0.047 0 0 0 0 0 0 0 0 1.0 0 1.0 0 1.0 0 0 0 0 0 0 0 0.75 0 0.75 0.75 0.75 0 0 0 0 0 0 0 0 0 0 0 0 0 0.75 0 0.75 0 0.75 0 0.56 0 0.56 0 0.56 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(1)^3] +DEAL::1.0 2.2 +DEAL::2.6 +DEAL::0.50 0 0.50 0 0.50 0 0.25 0 0.25 0 0.25 0 1.0 0 1.0 0.25 1.0 0 0.50 0 0.50 0 0.50 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.0 5.1 +DEAL::0.50 0 0.50 0 0.50 0 0.25 0 0.25 0 0.25 0 1.0 0 1.0 0.25 1.0 0 0.50 0 0.50 0 0.50 0 +DEAL::FESystem<3>[FE_Q<3>(2)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.6 3.1 +DEAL::4.7 +DEAL::0 0 0 0.016 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1.0 0 1.0 0.38 1.0 0 0 0 0 0 0 0 0.75 0 0.75 0 0.75 0 0 0 0 0 0 0 0.75 0 0.75 0 0.75 0 0.56 0 0.56 1.0 0.56 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(1)^3] +DEAL::1.0 2.2 +DEAL::2.6 +DEAL::0.25 0 0.25 0 0.25 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 1.0 0 1.0 0.25 1.0 0 +DEAL::FESystem<3>[FE_Q<3>(1)^3] vs. FESystem<3>[FE_Q<3>(2)^3] +DEAL::1.0 5.1 +DEAL::0.25 0 0.25 0 0.25 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 0.50 0 1.0 0 1.0 0.25 1.0 0 DEAL::Checking FE_DGQ<3>(2)2 against FE_DGQ<3>(3)2 in 3d: +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) DEAL::Checking FE_DGQ<3>(3)2 against FE_DGQ<3>(2)2 in 3d: +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(3)^2] vs. FESystem<3>[FE_DGQ<3>(2)^2] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGQ<3>(2)^2] vs. FESystem<3>[FE_DGQ<3>(3)^2] +DEAL::(Empty matrix) DEAL::Checking FE_DGP<3>(3)1 against FE_DGP<3>(1)1 in 3d: +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) DEAL::Checking FE_DGP<3>(1)1 against FE_DGP<3>(3)1 in 3d: +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(1)] vs. FESystem<3>[FE_DGP<3>(3)] +DEAL::(Empty matrix) +DEAL::FESystem<3>[FE_DGP<3>(3)] vs. FESystem<3>[FE_DGP<3>(1)] +DEAL::(Empty matrix) -- 2.39.5