From: Wolfgang Bangerth Date: Tue, 28 Mar 2006 22:54:31 +0000 (+0000) Subject: Add one assertion. X-Git-Tag: v8.0.0~11980 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=4a5a7330d9fd204e7d67247510bf7f42ed365900;p=dealii.git Add one assertion. git-svn-id: https://svn.dealii.org/trunk@12730 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/numerics/vectors.templates.h b/deal.II/deal.II/include/numerics/vectors.templates.h index 30daa10fd6..315069c2f2 100644 --- a/deal.II/deal.II/include/numerics/vectors.templates.h +++ b/deal.II/deal.II/include/numerics/vectors.templates.h @@ -887,7 +887,31 @@ interpolate_boundary_values (const Mapping &mapping, ++face_no) { const FiniteElement &fe = cell->get_fe(); - + + // we can presently deal only with + // primitive elements for boundary + // values. this does not preclude + // us using non-primitive elements + // in components that we aren't + // interested in, however. make + // sure that all shape functions + // that are non-zero for the + // components we are interested in, + // are in fact primitive + for (unsigned int i=0; iget_fe().dofs_per_cell; ++i) + { + const std::vector &nonzero_component_array + = cell->get_fe().get_nonzero_components (i); + for (unsigned int c=0; cget_fe().is_primitive (i), + ExcMessage ("This function can only deal with requested boundary " + "values that correspond to primitive (scalar) base " + "elements")); + } + typename DH::face_iterator face = cell->face(face_no); const unsigned char boundary_component = face->boundary_indicator(); if (function_map.find(boundary_component) != function_map.end()) @@ -928,7 +952,7 @@ interpolate_boundary_values (const Mapping &mapping, for (unsigned int i=0; i aux_quad (unit_support_points);