]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add FE_Nothing implemented by Joshua White.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 25 Sep 2009 02:47:25 +0000 (02:47 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 25 Sep 2009 02:47:25 +0000 (02:47 +0000)
git-svn-id: https://svn.dealii.org/trunk@19537 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_nothing.h [new file with mode: 0644]
deal.II/deal.II/source/fe/fe_nothing.cc [new file with mode: 0644]

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 (file)
index 0000000..2949e45
--- /dev/null
@@ -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 <base/config.h>
+#include <fe/fe.h>
+
+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 <int dim>
+class FE_Nothing : public FiniteElement<dim>
+{
+  public:
+
+    FE_Nothing ();
+
+    virtual
+    FiniteElement<dim> *
+    clone() const;
+
+    virtual
+    std::string
+    get_name() const;
+
+    virtual
+    unsigned int
+    n_base_elements () const;
+
+    virtual
+    const FiniteElement<dim> &
+    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<dim> &p) const;
+
+    virtual
+    void
+    fill_fe_values (const Mapping<dim> & mapping,
+                   const typename Triangulation<dim>::cell_iterator & cell,
+                   const Quadrature<dim> & quadrature,
+                   typename Mapping<dim>::InternalDataBase & mapping_data,
+                   typename Mapping<dim>::InternalDataBase & fedata,
+                   FEValuesData<dim,dim> & data,
+                   CellSimilarity::Similarity & cell_similarity) const;
+
+    virtual
+    void
+    fill_fe_face_values (const Mapping<dim> & mapping,
+                        const typename Triangulation<dim> :: cell_iterator & cell,
+                        const unsigned int face,
+                        const Quadrature<dim-1> & quadrature,
+                        typename Mapping<dim> :: InternalDataBase & mapping_data,
+                        typename Mapping<dim> :: InternalDataBase & fedata,
+                        FEValuesData<dim,dim> & data) const;
+
+    virtual
+    void
+    fill_fe_subface_values (const Mapping<dim> & mapping,
+                           const typename Triangulation<dim>::cell_iterator & cell,
+                           const unsigned int face,
+                           const unsigned int subface,
+                           const Quadrature<dim-1> & quadrature,
+                           typename Mapping<dim>::InternalDataBase & mapping_data,
+                           typename Mapping<dim>::InternalDataBase & fedata,
+                           FEValuesData<dim,dim> & data) const;
+
+    virtual
+    typename Mapping<dim>::InternalDataBase *
+    get_data (const UpdateFlags     update_flags,
+             const Mapping<dim>    & mapping,
+             const Quadrature<dim> & 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 (file)
index 0000000..edf5626
--- /dev/null
@@ -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 <fe/fe_nothing.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace
+{
+  const char*
+  zero_dof_message = "This element has no shape functions.";
+}
+
+
+
+
+template <int dim>
+FE_Nothing<dim>::FE_Nothing ()
+                :
+                FiniteElement<dim>
+               (FiniteElementData<dim>(std::vector<unsigned>(dim+1,0),
+                                       1, 0,
+                                       FiniteElementData<dim>::unknown),
+                std::vector<bool>(),
+                std::vector<std::vector<bool> >() )
+{}
+
+
+template <int dim>
+FiniteElement<dim> *
+FE_Nothing<dim>::clone() const
+{
+  return new FE_Nothing<dim>(*this);
+}
+
+
+
+template <int dim>
+std::string
+FE_Nothing<dim>::get_name () const
+{
+  std::ostringstream namebuf;
+  namebuf << "FE_Nothing<" << dim << ">()";
+  return namebuf.str();
+}
+
+
+
+template <int dim>
+unsigned int
+FE_Nothing<dim>::n_base_elements () const
+{
+  return 1;
+}
+
+
+
+template <int dim>
+const FiniteElement<dim> &
+FE_Nothing<dim>::base_element (const unsigned int index) const
+{
+  Assert (index==0, ExcIndexRange(index, 0, 1));
+  return *this;
+}
+
+
+
+template <int dim>
+unsigned int
+FE_Nothing<dim>::element_multiplicity (const unsigned int index) const
+{
+  Assert (index==0, ExcIndexRange(index, 0, 1));
+  return 1;
+}
+
+
+
+template <int dim>
+UpdateFlags
+FE_Nothing<dim>::update_once (const UpdateFlags /*flags*/) const
+{
+  return update_default;
+}
+
+
+
+template <int dim>
+UpdateFlags
+FE_Nothing<dim>::update_each (const UpdateFlags /*flags*/) const
+{
+  return update_default;
+}
+
+
+
+template <int dim>
+double
+FE_Nothing<dim>::shape_value (const unsigned int /*i*/,
+                             const Point<dim> & /*p*/) const
+{
+  Assert(false,ExcMessage(zero_dof_message));
+  return 0;
+}
+
+
+
+template <int dim>
+typename Mapping<dim>::InternalDataBase *
+FE_Nothing<dim>::get_data (const UpdateFlags  /*flags*/,
+                          const Mapping<dim> & /*mapping*/,
+                          const Quadrature<dim> & /*quadrature*/) const
+{
+  Assert(false,ExcMessage(zero_dof_message));
+  return NULL;
+}
+
+
+
+template <int dim>
+void
+FE_Nothing<dim>::
+fill_fe_values (const Mapping<dim>                      & /*mapping*/,
+               const typename Triangulation<dim>::cell_iterator & /*cell*/,
+               const Quadrature<dim> & /*quadrature*/,
+               typename Mapping<dim>::InternalDataBase & /*mapping_data*/,
+               typename Mapping<dim>::InternalDataBase & /*fedata*/,
+               FEValuesData<dim,dim> & /*data*/,
+               CellSimilarity::Similarity & /*cell_similarity*/) const
+{
+  Assert(false,ExcMessage(zero_dof_message));
+}
+
+
+
+template <int dim>
+void
+FE_Nothing<dim>::
+fill_fe_face_values (const Mapping<dim> & /*mapping*/,
+                    const typename Triangulation<dim>::cell_iterator & /*cell*/,
+                    const unsigned int /*face*/,
+                    const Quadrature<dim-1> & /*quadrature*/,
+                    typename Mapping<dim>::InternalDataBase & /*mapping_data*/,
+                    typename Mapping<dim>::InternalDataBase & /*fedata*/,
+                    FEValuesData<dim,dim> & /*data*/) const
+{
+  Assert(false,ExcMessage(zero_dof_message));
+}
+
+
+
+template <int dim>
+void
+FE_Nothing<dim>::
+fill_fe_subface_values (const Mapping<dim> & /*mapping*/,
+                       const typename Triangulation<dim>::cell_iterator & /*cell*/,
+                       const unsigned int /*face*/,
+                       const unsigned int /*subface*/,
+                       const Quadrature<dim-1> & /*quadrature*/,
+                       typename Mapping<dim>::InternalDataBase & /*mapping_data*/,
+                       typename Mapping<dim>::InternalDataBase & /*fedata*/,
+                       FEValuesData<dim,dim> & /*data*/) const
+{
+  Assert(false,ExcMessage(zero_dof_message));
+}
+
+
+
+template class FE_Nothing<deal_II_dimension>;
+
+DEAL_II_NAMESPACE_CLOSE
+

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.