]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix several bugs that stood in the way of making the things with FE_Nothing that...
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Sat, 5 Mar 2011 19:32:49 +0000 (19:32 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Sat, 5 Mar 2011 19:32:49 +0000 (19:32 +0000)
git-svn-id: https://svn.dealii.org/trunk@23461 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/fe/fe_base.h
deal.II/include/deal.II/fe/fe_q_hierarchical.h
deal.II/source/dofs/dof_tools.cc
deal.II/source/fe/fe_dgp.cc
deal.II/source/fe/fe_dgp_monomial.cc
deal.II/source/fe/fe_dgp_nonparametric.cc
deal.II/source/fe/fe_dgq.cc
deal.II/source/fe/fe_q_hierarchical.cc

index ec46db31f740e82657d60cb366cc46d1ef133799..d4b6b485525130022d83eda2ac19a989cea775fa 100644 (file)
@@ -682,14 +682,16 @@ namespace FiniteElementDomination
       {
        case this_element_dominates:
              if ((d2 == this_element_dominates) ||
-                 (d2 == either_element_can_dominate))
+                 (d2 == either_element_can_dominate) ||
+                 (d2 == no_requirements))
                return this_element_dominates;
              else
                return neither_element_dominates;
 
        case other_element_dominates:
              if ((d2 == other_element_dominates) ||
-                 (d2 == either_element_can_dominate))
+                 (d2 == either_element_can_dominate) ||
+                 (d2 == no_requirements))
                return other_element_dominates;
              else
                return neither_element_dominates;
index a4daa1f28f7ed57f80b64249756e78990e55a6f2..bff542ffe9c3be80ed6f1bf3120143a4bad2a454 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 2002, 2003, 2004, 2005, 2006, 2010 by the deal.II authors
+//    Copyright (C) 2002, 2003, 2004, 2005, 2006, 2010, 2011 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -352,6 +352,22 @@ class FE_Q_Hierarchical : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                      */
     virtual void get_subface_interpolation_matrix (const FiniteElement<dim>& source, const unsigned int subface, FullMatrix<double>& matrix) const;
 
+                                    /**
+                                     * Return whether this element dominates
+                                     * the one given as argument when they
+                                     * meet at a common face,
+                                     * whether it is the other way around,
+                                     * whether neither dominates, or if
+                                     * either could dominate.
+                                     *
+                                     * For a definition of domination, see
+                                     * FiniteElementBase::Domination and in
+                                     * particular the @ref hp_paper "hp paper".
+                                     */
+    virtual
+    FiniteElementDomination::Domination
+    compare_for_face_domination (const FiniteElement<dim> &fe_other) const;
+
                                     /**
                                      * Determine an estimate for the
                                      * memory consumption (in bytes)
index c5e41814567c8bc93a0b7e8ef44f08cbdd34712f..c2c6df733a51c97a129ee353cfe64c7449609162 100644 (file)
@@ -2266,34 +2266,6 @@ namespace internal
                Assert(cell->face(face)->refinement_case()==RefinementCase<dim-1>::isotropic_refinement,
                       ExcNotImplemented());
 
-                                                // check to see if the child elements
-                                                // have no dofs on the
-                                                // shared face.  if none of them do, we
-                                                // we can simply continue to the next face.
-                                                // if only some of them have dofs, while
-                                                // others do not, then we don't know
-                                                // what to do and throw an exception.
-                                                //
-                                                // ignore all
-                                                // interfaces with
-                                                // artificial cells
-               bool any_are_zero = false;
-               bool all_are_zero = true;
-
-               for (unsigned int c=0; c<cell->face(face)->n_children(); ++c)
-                 if (!cell->neighbor_child_on_subface(face, c)->is_artificial())
-                   {
-                     if (cell->neighbor_child_on_subface(face, c)->get_fe().dofs_per_face == 0)
-                       any_are_zero = true;
-                     else
-                       all_are_zero = false;
-                   }
-
-               if(all_are_zero)
-                 continue;
-
-               Assert( all_are_zero || !any_are_zero, ExcNotImplemented() );
-
                                                 // so now we've found a
                                                 // face of an active
                                                 // cell that has
@@ -2998,6 +2970,12 @@ namespace internal
                          break;
                        }
 
+                       case FiniteElementDomination::no_requirements:
+                       {
+                                                          // nothing to do here
+                         break;
+                       }
+                       
                        default:
                                                               // we shouldn't get
                                                               // here
index 3e6fe2c836a6170226f3ac8a2f11416d11e4b66d..1ca1b52780979bcf0d2ab41ed274b76f1311e497 100644 (file)
@@ -223,11 +223,10 @@ FiniteElementDomination::Domination
 FE_DGP<dim,spacedim>::compare_for_face_domination (const FiniteElement<dim,spacedim> &fe_other) const
 {
                                   // check whether both are discontinuous
-                                  // elements and both could dominate, see
-                                  // the description of
+                                  // elements, see the description of
                                   // FiniteElementDomination::Domination
   if (dynamic_cast<const FE_DGP<dim,spacedim>*>(&fe_other) != 0)
-    return FiniteElementDomination::either_element_can_dominate;
+    return FiniteElementDomination::no_requirements;
 
   Assert (false, ExcNotImplemented());
   return FiniteElementDomination::neither_element_dominates;
index ad2f656cea2139ca1aebc9c6f7c4010030a52ea1..83969bc58608f50c0070b55eb54e83e54217f51d 100644 (file)
@@ -367,11 +367,11 @@ FE_DGPMonomial<dim>::
 compare_for_face_domination (const FiniteElement<dim> &fe_other) const
 {
                                   // check whether both are discontinuous
-                                  // elements and both could dominate, see
+                                  // elements, see
                                   // the description of
                                   // FiniteElementDomination::Domination
   if (dynamic_cast<const FE_DGPMonomial<dim>*>(&fe_other) != 0)
-    return FiniteElementDomination::either_element_can_dominate;
+    return FiniteElementDomination::no_requirements;
 
   Assert (false, ExcNotImplemented());
   return FiniteElementDomination::neither_element_dominates;
index 36fbeb0c26495fe690d660da61ff1b9e4d768a75..ba5ff93f1cdd3f664feb56459a0ab6ace2c0c782 100644 (file)
@@ -573,11 +573,10 @@ FE_DGPNonparametric<dim,spacedim>::
 compare_for_face_domination (const FiniteElement<dim,spacedim> &fe_other) const
 {
                                   // check whether both are discontinuous
-                                  // elements and both could dominate, see
-                                  // the description of
+                                  // elements, see the description of
                                   // FiniteElementDomination::Domination
   if (dynamic_cast<const FE_DGPNonparametric<dim,spacedim>*>(&fe_other) != 0)
-    return FiniteElementDomination::either_element_can_dominate;
+    return FiniteElementDomination::no_requirements;
 
   Assert (false, ExcNotImplemented());
   return FiniteElementDomination::neither_element_dominates;
index 11d3216cca11ad3c85d968f7d17054c015871d2a..d0667da04f5f8713232fa91126924ccde350d937 100644 (file)
@@ -576,11 +576,10 @@ FiniteElementDomination::Domination
 FE_DGQ<dim, spacedim>::compare_for_face_domination (const FiniteElement<dim, spacedim> &fe_other) const
 {
                                   // check whether both are discontinuous
-                                  // elements and both could dominate, see
-                                  // the description of
+                                  // elements, see the description of
                                   // FiniteElementDomination::Domination
   if (dynamic_cast<const FE_DGQ<dim, spacedim>*>(&fe_other) != 0)
-    return FiniteElementDomination::either_element_can_dominate;
+    return FiniteElementDomination::no_requirements;
 
   Assert (false, ExcNotImplemented());
   return FiniteElementDomination::neither_element_dominates;
index 6f6a888ac635fc700c6aea381dbac4b2bcb8cc60..348c7cfa56817d8dbc65d572242e24bb4118fc18 100644 (file)
@@ -164,6 +164,38 @@ hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const
     }
 }
 
+
+template <int dim>
+FiniteElementDomination::Domination
+FE_Q_Hierarchical<dim>::
+compare_for_face_domination (const FiniteElement<dim> &fe_other) const
+{
+  if (const FE_Q_Hierarchical<dim> *fe_q_other
+      = dynamic_cast<const FE_Q_Hierarchical<dim>*>(&fe_other))
+    {
+      if (this->degree < fe_q_other->degree)
+       return FiniteElementDomination::this_element_dominates;
+      else if (this->degree == fe_q_other->degree)
+       return FiniteElementDomination::either_element_can_dominate;
+      else
+       return FiniteElementDomination::other_element_dominates;
+    }
+  else if (dynamic_cast<const FE_Nothing<dim>*>(&fe_other) != 0)
+    {
+                                      // 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;
+    }
+
+  Assert (false, ExcNotImplemented());
+  return FiniteElementDomination::neither_element_dominates;
+}
+
+
 //---------------------------------------------------------------------------
 // Auxiliary functions
 //---------------------------------------------------------------------------

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.