From 0ef6e96485b45c8cb11882eab87ee4749e22b8dc Mon Sep 17 00:00:00 2001 From: David Wells Date: Sun, 22 Apr 2018 19:46:05 -0400 Subject: [PATCH] Change step-10 to use FE_Nothing. --- examples/step-10/step-10.cc | 24 +++++++++++------------- include/deal.II/fe/fe_nothing.h | 6 +++--- 2 files changed, 14 insertions(+), 16 deletions(-) diff --git a/examples/step-10/step-10.cc b/examples/step-10/step-10.cc index 88469c70f8..cd1c5ac496 100644 --- a/examples/step-10/step-10.cc +++ b/examples/step-10/step-10.cc @@ -30,7 +30,7 @@ #include #include #include -#include +#include #include // This is the only new one: in it, we declare the MappingQ class @@ -207,15 +207,14 @@ namespace Step10 const MappingQ mapping (degree); - // We now create a dummy finite element. Here we could choose any - // finite element, as we are only interested in the `JxW' values - // provided by the FEValues object below. Nevertheless, we have to - // provide a finite element since in this example we abuse the - // FEValues class a little in that we only ask it to provide us with - // the weights of certain quadrature points, in contrast to the usual - // purpose (and name) of the FEValues class which is to provide the - // values of finite elements at these points. - const FE_Q dummy_fe (1); + // We now create a finite element. Unlike the rest of the example + // programs, we do not actually need to do any computations with shape + // functions; we only need the `JxW' values from an FEValues + // object. Hence we use the special finite element class FE_Nothing + // which has exactly zero degrees of freedom per cell (as the name + // implies, the local basis on each cell is the empty set). A more + // typical usage of FE_Nothing is shown in step-46. + const FE_Nothing fe; // Likewise, we need to create a DoFHandler object. We do not actually // use it, but it will provide us with `active_cell_iterators' that @@ -236,8 +235,7 @@ namespace Step10 // computation of the mapping from unit to real cell. In previous // examples, this argument was omitted, resulting in the implicit use // of an object of type MappingQ1. - FEValues fe_values (mapping, dummy_fe, quadrature, - update_JxW_values); + FEValues fe_values (mapping, fe, quadrature, update_JxW_values); // We employ an object of the ConvergenceTable class to store all // important data like the approximated values for $\pi$ and the error @@ -261,7 +259,7 @@ namespace Step10 // our special case but we call it to make the DoFHandler happy -- // otherwise it would throw an assertion in the FEValues::reinit // function below. - dof_handler.distribute_dofs (dummy_fe); + dof_handler.distribute_dofs (fe); // We define the variable area as `long double' like we did for // the pi variable before. diff --git a/include/deal.II/fe/fe_nothing.h b/include/deal.II/fe/fe_nothing.h index 88a0ad6f47..9e4c8f60ca 100644 --- a/include/deal.II/fe/fe_nothing.h +++ b/include/deal.II/fe/fe_nothing.h @@ -34,9 +34,9 @@ DEAL_II_NAMESPACE_OPEN * active region where normal elements are used, and an inactive region where * FE_Nothing elements are used. The hp::DoFHandler will therefore assign no * degrees of freedom to the FE_Nothing cells, and this subregion is therefore - * implicitly deleted from the computation. step-46 shows a use case for this - * element. An interesting application for this element is also presented in - * the paper A. Cangiani, J. Chapman, E. Georgoulis, M. Jensen: + * implicitly deleted from the computation. step-10 and step-46 show use cases + * for this element. An interesting application for this element is also + * presented in the paper A. Cangiani, J. Chapman, E. Georgoulis, M. Jensen: * Implementation of the Continuous-Discontinuous Galerkin Finite Element * Method, arXiv:1201.2878v1 [math.NA], 2012 (see * http://arxiv.org/abs/1201.2878). -- 2.39.5