From: bangerth Date: Fri, 25 Sep 2009 02:47:25 +0000 (+0000) Subject: Add FE_Nothing implemented by Joshua White. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=08af457849ff80ce6629d235b0e1088be742567a;p=dealii-svn.git Add FE_Nothing implemented by Joshua White. git-svn-id: https://svn.dealii.org/trunk@19537 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_nothing.h b/deal.II/deal.II/include/fe/fe_nothing.h new file mode 100644 index 0000000000..2949e45c06 --- /dev/null +++ b/deal.II/deal.II/include/fe/fe_nothing.h @@ -0,0 +1,120 @@ +//--------------------------------------------------------------------------- +// $Id: fe_q.h 19241 2009-08-12 14:27:43Z bangerth $ +// Version: $Name$ +// +// Copyright (C) 2009 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//--------------------------------------------------------------------------- +#ifndef __deal2__fe_nothing_h +#define __deal2__fe_nothing_h + +#include +#include + +DEAL_II_NAMESPACE_OPEN + + +/*!@addtogroup fe */ +/*@{*/ + +/** + * Definition of a finite element with zero degrees of freedom. This + * class is useful (in the context of an hp method) to represent empty + * cells in the triangulation on which no degrees of freedom should + * be allocated. Thus a triangulation may be divided into two regions: + * an 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. + * + * @author Joshua White + */ +template +class FE_Nothing : public FiniteElement +{ + public: + + FE_Nothing (); + + virtual + FiniteElement * + clone() const; + + virtual + std::string + get_name() const; + + virtual + unsigned int + n_base_elements () const; + + virtual + const FiniteElement & + base_element (const unsigned int index) const; + + virtual + unsigned int + element_multiplicity (const unsigned int index) const; + + virtual + UpdateFlags + update_once (const UpdateFlags flags) const; + + virtual + UpdateFlags + update_each (const UpdateFlags flags) const; + + virtual + double + shape_value (const unsigned int i, const Point &p) const; + + virtual + void + fill_fe_values (const Mapping & mapping, + const typename Triangulation::cell_iterator & cell, + const Quadrature & quadrature, + typename Mapping::InternalDataBase & mapping_data, + typename Mapping::InternalDataBase & fedata, + FEValuesData & data, + CellSimilarity::Similarity & cell_similarity) const; + + virtual + void + fill_fe_face_values (const Mapping & mapping, + const typename Triangulation :: cell_iterator & cell, + const unsigned int face, + const Quadrature & quadrature, + typename Mapping :: InternalDataBase & mapping_data, + typename Mapping :: InternalDataBase & fedata, + FEValuesData & data) const; + + virtual + void + fill_fe_subface_values (const Mapping & mapping, + const typename Triangulation::cell_iterator & cell, + const unsigned int face, + const unsigned int subface, + const Quadrature & quadrature, + typename Mapping::InternalDataBase & mapping_data, + typename Mapping::InternalDataBase & fedata, + FEValuesData & data) const; + + virtual + typename Mapping::InternalDataBase * + get_data (const UpdateFlags update_flags, + const Mapping & mapping, + const Quadrature & quadrature) const; +}; + + +/*@}*/ + +DEAL_II_NAMESPACE_CLOSE + +#endif + diff --git a/deal.II/deal.II/source/fe/fe_nothing.cc b/deal.II/deal.II/source/fe/fe_nothing.cc new file mode 100644 index 0000000000..edf56265aa --- /dev/null +++ b/deal.II/deal.II/source/fe/fe_nothing.cc @@ -0,0 +1,181 @@ +//--------------------------------------------------------------------------- +// $Id: fe_q.cc 18913 2009-06-05 08:02:16Z kronbichler $ +// Version: $Name$ +// +// Copyright (C) 2009 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//--------------------------------------------------------------------------- + +#include + +DEAL_II_NAMESPACE_OPEN + +namespace +{ + const char* + zero_dof_message = "This element has no shape functions."; +} + + + + +template +FE_Nothing::FE_Nothing () + : + FiniteElement + (FiniteElementData(std::vector(dim+1,0), + 1, 0, + FiniteElementData::unknown), + std::vector(), + std::vector >() ) +{} + + +template +FiniteElement * +FE_Nothing::clone() const +{ + return new FE_Nothing(*this); +} + + + +template +std::string +FE_Nothing::get_name () const +{ + std::ostringstream namebuf; + namebuf << "FE_Nothing<" << dim << ">()"; + return namebuf.str(); +} + + + +template +unsigned int +FE_Nothing::n_base_elements () const +{ + return 1; +} + + + +template +const FiniteElement & +FE_Nothing::base_element (const unsigned int index) const +{ + Assert (index==0, ExcIndexRange(index, 0, 1)); + return *this; +} + + + +template +unsigned int +FE_Nothing::element_multiplicity (const unsigned int index) const +{ + Assert (index==0, ExcIndexRange(index, 0, 1)); + return 1; +} + + + +template +UpdateFlags +FE_Nothing::update_once (const UpdateFlags /*flags*/) const +{ + return update_default; +} + + + +template +UpdateFlags +FE_Nothing::update_each (const UpdateFlags /*flags*/) const +{ + return update_default; +} + + + +template +double +FE_Nothing::shape_value (const unsigned int /*i*/, + const Point & /*p*/) const +{ + Assert(false,ExcMessage(zero_dof_message)); + return 0; +} + + + +template +typename Mapping::InternalDataBase * +FE_Nothing::get_data (const UpdateFlags /*flags*/, + const Mapping & /*mapping*/, + const Quadrature & /*quadrature*/) const +{ + Assert(false,ExcMessage(zero_dof_message)); + return NULL; +} + + + +template +void +FE_Nothing:: +fill_fe_values (const Mapping & /*mapping*/, + const typename Triangulation::cell_iterator & /*cell*/, + const Quadrature & /*quadrature*/, + typename Mapping::InternalDataBase & /*mapping_data*/, + typename Mapping::InternalDataBase & /*fedata*/, + FEValuesData & /*data*/, + CellSimilarity::Similarity & /*cell_similarity*/) const +{ + Assert(false,ExcMessage(zero_dof_message)); +} + + + +template +void +FE_Nothing:: +fill_fe_face_values (const Mapping & /*mapping*/, + const typename Triangulation::cell_iterator & /*cell*/, + const unsigned int /*face*/, + const Quadrature & /*quadrature*/, + typename Mapping::InternalDataBase & /*mapping_data*/, + typename Mapping::InternalDataBase & /*fedata*/, + FEValuesData & /*data*/) const +{ + Assert(false,ExcMessage(zero_dof_message)); +} + + + +template +void +FE_Nothing:: +fill_fe_subface_values (const Mapping & /*mapping*/, + const typename Triangulation::cell_iterator & /*cell*/, + const unsigned int /*face*/, + const unsigned int /*subface*/, + const Quadrature & /*quadrature*/, + typename Mapping::InternalDataBase & /*mapping_data*/, + typename Mapping::InternalDataBase & /*fedata*/, + FEValuesData & /*data*/) const +{ + Assert(false,ExcMessage(zero_dof_message)); +} + + + +template class FE_Nothing; + +DEAL_II_NAMESPACE_CLOSE +