]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Also assert in read_dof_values_plain. 5212/head
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Mon, 9 Oct 2017 10:40:20 +0000 (12:40 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Mon, 9 Oct 2017 10:40:20 +0000 (12:40 +0200)
include/deal.II/matrix_free/fe_evaluation.h

index dafb5161844a9a88571cd24779de896cf363fdf0..1170a3075db8cb32b984b0334c63c52dd637910b 100644 (file)
@@ -3455,6 +3455,16 @@ FEEvaluationBase<dim,n_components_,Number>
   // and sit on a different vector each)
   if (n_fe_components == 1)
     {
+      for (unsigned int c=0; c<n_components; ++c)
+        Assert(src[c] != nullptr,
+               ExcMessage("The finite element underlying this FEEvaluation "
+                          "object is scalar, but you requested " +
+                          std::to_string(n_components) +
+                          " components via the template argument in "
+                          "FEEvaluation. In that case, you must pass an "
+                          "std::vector<VectorType> or a BlockVector to " +
+                          "read_dof_values_plain."));
+
       const unsigned int n_local_dofs =
         VectorizedArray<Number>::n_array_elements * dofs_per_cell;
       for (unsigned int comp=0; comp<n_components; ++comp)

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.