From: bangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Date: Fri, 17 Jan 2014 22:02:02 +0000 (+0000)
Subject: Disallow calling some functions in the hp context on non-active cells.
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=e8a86a7f12c372e4f45419c1a982751d2b5b4426;p=dealii-svn.git

Disallow calling some functions in the hp context on non-active cells.

git-svn-id: https://svn.dealii.org/trunk@32239 0785d39b-7218-0410-832d-ea1e28bc413d
---

diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h
index 54e3ae9172..aa8a9151e8 100644
--- a/deal.II/doc/news/changes.h
+++ b/deal.II/doc/news/changes.h
@@ -85,6 +85,15 @@ inconvenience this causes.
 
 <ol>
 
+  <li> Fixed: DoFCellAccessor::set_dof_values_by_interpolation and
+  DoFCellAccessor::get_interpolated_dof_values could previously be
+  called for hp::DoFHandler objects on cells that are non-active. This
+  makes no sense since these cells have no associated finite element
+  space. Doing so now raises an exception.
+  <br>
+  (Wolfgang Bangerth, 2014/01/17)
+  </li>
+
   <li> Fixed: The methods IndexSet::do_compress() and
   IndexSet::add_indices(IndexSet&) had quadratic complexity in the number of
   ranges. The algorithms have been changed into linear complexity ones.
diff --git a/deal.II/source/dofs/dof_accessor_get.cc b/deal.II/source/dofs/dof_accessor_get.cc
index 9b9ca8b7b8..33639cf5bc 100644
--- a/deal.II/source/dofs/dof_accessor_get.cc
+++ b/deal.II/source/dofs/dof_accessor_get.cc
@@ -1,7 +1,7 @@
 // ---------------------------------------------------------------------
 // $Id$
 //
-// Copyright (C) 1998 - 2013 by the deal.II authors
+// Copyright (C) 1998 - 2014 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -45,25 +45,41 @@ DoFCellAccessor<DH,lda>::
 get_interpolated_dof_values (const InputVector &values,
                              Vector<number>    &interpolated_values) const
 {
-  const FiniteElement<dim,spacedim> &fe            = this->get_fe();
-  const unsigned int        dofs_per_cell = fe.dofs_per_cell;
-
-  Assert (this->dof_handler != 0,
-          typename BaseClass::ExcInvalidObject());
-  Assert (&fe != 0,
-          typename BaseClass::ExcInvalidObject());
-  Assert (interpolated_values.size() == dofs_per_cell,
-          typename BaseClass::ExcVectorDoesNotMatch());
-  Assert (values.size() == this->dof_handler->n_dofs(),
-          typename BaseClass::ExcVectorDoesNotMatch());
-
   if (!this->has_children())
     // if this cell has no children: simply return the exact values on this
     // cell
     this->get_dof_values (values, interpolated_values);
   else
-    // otherwise clobber them from the children
+    // otherwise obtain them from the children
     {
+      // we are on a non-active cell. these do not have any finite
+      // element associated with them in the hp context (in the non-hp
+      // context, we can simply assume that the FE space to which we
+      // want to interpolate is the same as for all elements in the
+      // mesh). consequently, we cannot interpolate from children's FE
+      // space to this cell's (unknown) FE space
+      Assert ((dynamic_cast<DoFHandler<DH::dimension,DH::space_dimension>*>
+	       (this->dof_handler)
+	       != 0),
+	      ExcMessage ("You cannot call this function on non-active cells "
+			  "of hp::DoFHandler objects because they do not have "
+			  "associated finite element spaces associated: degrees "
+			  "of freedom are only distributed on active cells for which "
+			  "the active_fe_index has been set."));
+
+      const FiniteElement<dim,spacedim> &fe            = this->get_fe();
+      const unsigned int                 dofs_per_cell = fe.dofs_per_cell;
+
+      Assert (this->dof_handler != 0,
+	      typename BaseClass::ExcInvalidObject());
+      Assert (&fe != 0,
+	      typename BaseClass::ExcInvalidObject());
+      Assert (interpolated_values.size() == dofs_per_cell,
+	      typename BaseClass::ExcVectorDoesNotMatch());
+      Assert (values.size() == this->dof_handler->n_dofs(),
+	      typename BaseClass::ExcVectorDoesNotMatch());
+
+
       Vector<number> tmp1(dofs_per_cell);
       Vector<number> tmp2(dofs_per_cell);
 
diff --git a/deal.II/source/dofs/dof_accessor_set.cc b/deal.II/source/dofs/dof_accessor_set.cc
index f11fc845a0..4438c6b5b5 100644
--- a/deal.II/source/dofs/dof_accessor_set.cc
+++ b/deal.II/source/dofs/dof_accessor_set.cc
@@ -1,7 +1,7 @@
 // ---------------------------------------------------------------------
 // $Id$
 //
-// Copyright (C) 1998 - 2013 by the deal.II authors
+// Copyright (C) 1998 - 2014 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -46,23 +46,38 @@ DoFCellAccessor<DH,lda>::
 set_dof_values_by_interpolation (const Vector<number> &local_values,
                                  OutputVector         &values) const
 {
-  const unsigned int dofs_per_cell = this->get_fe().dofs_per_cell;
-
-  Assert (this->dof_handler != 0,
-          typename BaseClass::ExcInvalidObject());
-  Assert (&this->get_fe() != 0,
-          typename BaseClass::ExcInvalidObject());
-  Assert (local_values.size() == dofs_per_cell,
-          typename BaseClass::ExcVectorDoesNotMatch());
-  Assert (values.size() == this->dof_handler->n_dofs(),
-          typename BaseClass::ExcVectorDoesNotMatch());
-
   if (!this->has_children())
     // if this cell has no children: simply set the values on this cell
     this->set_dof_values (local_values, values);
   else
     // otherwise distribute them to the children
     {
+      // we are on a non-active cell. these do not have any finite
+      // element associated with them in the hp context (in the non-hp
+      // context, we can simply assume that the FE space to from which we
+      // want to interpolate is the same as for all elements in the
+      // mesh). consequently, we cannot interpolate to children's FE
+      // space from this cell's (unknown) FE space
+      Assert ((dynamic_cast<DoFHandler<DH::dimension,DH::space_dimension>*>
+	       (this->dof_handler)
+	       != 0),
+	      ExcMessage ("You cannot call this function on non-active cells "
+			  "of hp::DoFHandler objects because they do not have "
+			  "associated finite element spaces associated: degrees "
+			  "of freedom are only distributed on active cells for which "
+			  "the active_fe_index has been set."));
+
+      const unsigned int dofs_per_cell = this->get_fe().dofs_per_cell;
+
+      Assert (this->dof_handler != 0,
+	      typename BaseClass::ExcInvalidObject());
+      Assert (&this->get_fe() != 0,
+	      typename BaseClass::ExcInvalidObject());
+      Assert (local_values.size() == dofs_per_cell,
+	      typename BaseClass::ExcVectorDoesNotMatch());
+      Assert (values.size() == this->dof_handler->n_dofs(),
+	      typename BaseClass::ExcVectorDoesNotMatch());
+
       Vector<number> tmp(dofs_per_cell);
 
       for (unsigned int child=0; child<this->n_children(); ++child)