]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Combine four DataOut_DoFData::attach_data_vector functions.
authorDavid Wells <wellsd2@rpi.edu>
Sat, 13 May 2017 17:47:04 +0000 (13:47 -0400)
committerDavid Wells <wellsd2@rpi.edu>
Fri, 2 Jun 2017 12:52:53 +0000 (08:52 -0400)
include/deal.II/numerics/data_out_dof_data.h
source/numerics/data_out_dof_data.cc
source/numerics/data_out_dof_data.inst.in
tests/numerics/data_out_11.debug.output

index 34afa0be14d237f981833b29a06f6ec9471a11ee..c23b710b047089b05437ff198001371f92fbab5d 100644 (file)
@@ -35,9 +35,6 @@
 
 DEAL_II_NAMESPACE_OPEN
 
-template <int, int> class FEValuesBase;
-
-
 namespace Exceptions
 {
   /**
@@ -875,11 +872,99 @@ protected:
    */
   template <class, int, int>
   friend class DataOut_DoFData;
+private:
+  /**
+   * Common function called by the four public add_data_vector methods.
+   */
+  template <class VectorType>
+  void
+  add_data_vector_internal
+  (const DoFHandlerType           *dof_handler,
+   const VectorType               &data,
+   const std::vector<std::string> &names,
+   const DataVectorType            type,
+   const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation,
+   const bool                      deduce_output_names);
 };
 
 
 
 // -------------------- template and inline functions ------------------------
+template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
+template <typename VectorType>
+void
+DataOut_DoFData<DoFHandlerType,patch_dim,patch_space_dim>::
+add_data_vector
+(const VectorType                         &vec,
+ const std::string                        &name,
+ const DataVectorType                      type,
+ const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation)
+{
+  Assert (triangulation != nullptr, Exceptions::DataOut::ExcNoTriangulationSelected ());
+  std::vector<std::string> names(1, name);
+  add_data_vector_internal (dofs, vec, names, type, data_component_interpretation, true);
+}
+
+
+
+template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
+template <typename VectorType>
+void
+DataOut_DoFData<DoFHandlerType,patch_dim,patch_space_dim>::
+add_data_vector
+(const VectorType                         &vec,
+ const std::vector<std::string>           &names,
+ const DataVectorType                      type,
+ const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation)
+{
+  Assert (triangulation != nullptr, Exceptions::DataOut::ExcNoTriangulationSelected ());
+  add_data_vector_internal(dofs, vec, names, type, data_component_interpretation, false);
+}
+
+
+
+template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
+template <typename VectorType>
+void
+DataOut_DoFData<DoFHandlerType,patch_dim,patch_space_dim>::
+add_data_vector
+(const DoFHandlerType           &dof_handler,
+ const VectorType               &data,
+ const std::string              &name,
+ const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation)
+{
+  std::vector<std::string> names(1, name);
+  add_data_vector_internal (&dof_handler, data, names, type_dof_data, data_component_interpretation, true);
+}
+
+
+
+template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
+template <typename VectorType>
+void
+DataOut_DoFData<DoFHandlerType,patch_dim,patch_space_dim>::
+add_data_vector
+(const DoFHandlerType           &dof_handler,
+ const VectorType               &data,
+ const std::vector<std::string> &names,
+ const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation)
+{
+  add_data_vector_internal(&dof_handler, data, names, type_dof_data, data_component_interpretation, false);
+}
+
+
+
+template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
+template <typename VectorType>
+void
+DataOut_DoFData<DoFHandlerType,patch_dim,patch_space_dim>::
+add_data_vector (const VectorType                       &vec,
+                 const DataPostprocessor<DoFHandlerType::space_dimension> &data_postprocessor)
+{
+  Assert (dofs != nullptr, Exceptions::DataOut::ExcNoDoFHandlerSelected ());
+  add_data_vector(*dofs, vec, data_postprocessor);
+}
+
 
 
 template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
index 3dba0c37c2cb87f477a36b0ee3e0820f9332bd54..00d89adee2734a69c9a592f97a6d3b8d6808616a 100644 (file)
@@ -883,246 +883,157 @@ attach_triangulation (const Triangulation<DoFHandlerType::dimension,DoFHandlerTy
 
 
 
-
-template <typename DoFHandlerType,
-          int patch_dim, int patch_space_dim>
+template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
 template <typename VectorType>
 void
 DataOut_DoFData<DoFHandlerType,patch_dim,patch_space_dim>::
-add_data_vector (const VectorType                         &vec,
-                 const std::string                        &name,
-                 const DataVectorType                      type,
-                 const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation)
+add_data_vector (const DoFHandlerType                   &dof_handler,
+                 const VectorType                       &vec,
+                 const DataPostprocessor<DoFHandlerType::space_dimension> &data_postprocessor)
 {
-  Assert (triangulation != nullptr,
-          Exceptions::DataOut::ExcNoTriangulationSelected ());
-  const unsigned int n_components =
-    dofs != nullptr ? dofs->get_fe().n_components () : 1;
-
-  std::vector<std::string> names;
-  // if only one component or vector is cell vector: we only need one name
-  if ((n_components == 1) ||
-      (vec.size() == triangulation->n_active_cells()))
+  // this is a specialized version of the other function where we have a
+  // postprocessor. if we do, we know that we have type_dof_data, which makes
+  // things a bit simpler, we also don't need to deal with some of the other
+  // stuff and use a different constructor of DataEntry
+  if (triangulation != nullptr)
     {
-      names.resize (1, name);
+      Assert (&dof_handler.get_triangulation() == triangulation,
+              ExcMessage("The triangulation attached to the DoFHandler does not "
+                         "match with the one set previously"));
     }
   else
-    // otherwise append _i to the given name
     {
-      names.resize (n_components);
-      for (unsigned int i=0; i<n_components; ++i)
-        {
-          std::ostringstream namebuf;
-          namebuf << '_' << i;
-          names[i] = name + namebuf.str();
-        }
+      triangulation = SmartPointer<const Triangulation<DoFHandlerType::dimension,DoFHandlerType::space_dimension>>
+        (&dof_handler.get_triangulation(), typeid(*this).name());
     }
 
-  add_data_vector (vec, names, type, data_component_interpretation);
+  Assert (vec.size() == dof_handler.n_dofs(),
+          Exceptions::DataOut::ExcInvalidVectorSize (vec.size(),
+                                                     dof_handler.n_dofs(),
+                                                     dof_handler.get_triangulation().n_active_cells()));
+
+
+  auto new_entry = std_cxx14::make_unique<internal::DataOut::DataEntry<DoFHandlerType,VectorType>>
+                   (&dof_handler, &vec, &data_postprocessor);
+  dof_data.emplace_back (std::move(new_entry));
 }
 
 
 
-template <typename DoFHandlerType,
-          int patch_dim, int patch_space_dim>
+template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
 template <typename VectorType>
 void
 DataOut_DoFData<DoFHandlerType,patch_dim,patch_space_dim>::
-add_data_vector (const VectorType                         &vec,
-                 const std::vector<std::string>           &names,
-                 const DataVectorType                      type,
-                 const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation_)
+add_data_vector_internal
+(const DoFHandlerType           *dof_handler,
+ const VectorType               &data_vector,
+ const std::vector<std::string> &names,
+ const DataVectorType            type,
+ const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation_,
+ const bool                      deduce_output_names)
 {
-  Assert (triangulation != nullptr,
-          Exceptions::DataOut::ExcNoTriangulationSelected ());
+  // Check available mesh information:
+  if (triangulation == nullptr)
+    {
+      Assert(dof_handler != nullptr, ExcInternalError());
+      triangulation = SmartPointer<const Triangulation<DoFHandlerType::dimension,DoFHandlerType::space_dimension>>
+                      (&dof_handler->get_triangulation(), typeid(*this).name());
+    }
 
-  const std::vector<DataComponentInterpretation::DataComponentInterpretation> &
-  data_component_interpretation
-    = (data_component_interpretation_.size() != 0
-       ?
-       data_component_interpretation_
-       :
-       std::vector<DataComponentInterpretation::DataComponentInterpretation>
-       (names.size(), DataComponentInterpretation::component_is_scalar));
+  if (dof_handler != nullptr)
+    {
+      Assert (&dof_handler->get_triangulation() == triangulation,
+              ExcMessage("The triangulation attached to the DoFHandler does not "
+                         "match with the one set previously"));
+    }
 
-  // either cell data and one name,
-  // or dof data and n_components names
+  // Figure out the data type:
   DataVectorType actual_type = type;
   if (type == type_automatic)
     {
-      // in the rare case that someone has a DGP(0) attached, we can not decide what she wants here:
-      Assert((dofs == nullptr) || (triangulation->n_active_cells() != dofs->n_dofs()),
+      Assert((dof_handler == nullptr) || (triangulation->n_active_cells() != dof_handler->n_dofs()),
              ExcMessage("Unable to determine the type of vector automatically because the number of DoFs "
                         "is equal to the number of cells. Please specify DataVectorType."));
 
-      if (vec.size() == triangulation->n_active_cells())
+      if (data_vector.size() == triangulation->n_active_cells())
         actual_type = type_cell_data;
       else
         actual_type = type_dof_data;
     }
+  Assert(actual_type == type_cell_data || actual_type == type_dof_data,
+         ExcInternalError());
 
+  // If necessary, append '_1', '_2', etc. to component names:
+  std::vector<std::string> deduced_names;
+  if (deduce_output_names && actual_type == type_dof_data)
+    {
+      Assert(names.size() == 1, ExcInternalError());
+      Assert(dof_handler != nullptr, ExcInternalError());
+      Assert(dof_handler->n_dofs() > 0,
+             ExcMessage("The DoF handler attached to the current output vector "
+                        "does not have any degrees of freedom, so it is not "
+                        "possible to output DoF data in this context."));
+      const std::string name = names[0];
+      const unsigned int n_components = dof_handler->get_fe().n_components();
+      deduced_names.resize (n_components);
+      if (n_components > 1)
+        {
+          for (unsigned int i=0; i<n_components; ++i)
+            {
+              deduced_names[i] = name + '_' + Utilities::to_string(i);
+            }
+        }
+      else
+        {
+          deduced_names[0] = name;
+        }
+    }
+  else
+    {
+      deduced_names = names;
+    }
+
+  // Check that things have the right sizes for the data type:
   switch (actual_type)
     {
     case type_cell_data:
-      Assert (vec.size() == triangulation->n_active_cells(),
-              ExcDimensionMismatch (vec.size(),
+      Assert (data_vector.size() == triangulation->n_active_cells(),
+              ExcDimensionMismatch (data_vector.size(),
                                     triangulation->n_active_cells()));
-      Assert (names.size() == 1,
-              Exceptions::DataOut::ExcInvalidNumberOfNames (names.size(), 1));
+      Assert (deduced_names.size() == 1,
+              Exceptions::DataOut::ExcInvalidNumberOfNames (deduced_names.size(), 1));
       break;
-
     case type_dof_data:
-      Assert (dofs != nullptr,
+      Assert (dof_handler != nullptr,
               Exceptions::DataOut::ExcNoDoFHandlerSelected ());
-      Assert (vec.size() == dofs->n_dofs(),
-              Exceptions::DataOut::ExcInvalidVectorSize (vec.size(),
-                                                         dofs->n_dofs(),
+      Assert (data_vector.size() == dof_handler->n_dofs(),
+              Exceptions::DataOut::ExcInvalidVectorSize (data_vector.size(),
+                                                         dof_handler->n_dofs(),
                                                          triangulation->n_active_cells()));
-      Assert (names.size() == dofs->get_fe().n_components(),
-              Exceptions::DataOut::ExcInvalidNumberOfNames (names.size(),
-                                                            dofs->get_fe().n_components()));
+      Assert (deduced_names.size() == dof_handler->get_fe().n_components(),
+              Exceptions::DataOut::ExcInvalidNumberOfNames (deduced_names.size(),
+                                                            dof_handler->get_fe().n_components()));
       break;
-
-    case type_automatic:
-      // this case should have been handled above...
+    default:
       Assert (false, ExcInternalError());
     }
 
-  auto new_entry = std_cxx14::make_unique<internal::DataOut::DataEntry<DoFHandlerType,VectorType>>
-                   (dofs, &vec, names, data_component_interpretation);
-
-  if (actual_type == type_dof_data)
-    dof_data.emplace_back (std::move(new_entry));
-  else
-    cell_data.emplace_back (std::move(new_entry));
-}
-
-
-
-template <typename DoFHandlerType,
-          int patch_dim, int patch_space_dim>
-template <typename VectorType>
-void
-DataOut_DoFData<DoFHandlerType,patch_dim,patch_space_dim>::
-add_data_vector (const VectorType                       &vec,
-                 const DataPostprocessor<DoFHandlerType::space_dimension> &data_postprocessor)
-{
-  // this is a specialized version of the other function where we have a
-  // postprocessor. if we do, we know that we have type_dof_data, which makes
-  // things a bit simpler, we also don't need to deal with some of the other
-  // stuff and use a different constructor of DataEntry
-
-  Assert (dofs != nullptr,
-          Exceptions::DataOut::ExcNoDoFHandlerSelected ());
-
-  Assert (vec.size() == dofs->n_dofs(),
-          Exceptions::DataOut::ExcInvalidVectorSize (vec.size(),
-                                                     dofs->n_dofs(),
-                                                     dofs->get_triangulation().n_active_cells()));
-
-  auto new_entry = std_cxx14::make_unique<internal::DataOut::DataEntry<DoFHandlerType,VectorType>>
-                   (dofs, &vec, &data_postprocessor);
-  dof_data.emplace_back (std::move(new_entry));
-}
-
-
-
-template <typename DoFHandlerType,
-          int patch_dim, int patch_space_dim>
-template <typename VectorType>
-void
-DataOut_DoFData<DoFHandlerType,patch_dim,patch_space_dim>::
-add_data_vector (const DoFHandlerType                   &dof_handler,
-                 const VectorType                       &vec,
-                 const DataPostprocessor<DoFHandlerType::space_dimension> &data_postprocessor)
-{
-  // this is a specialized version of the other function where we have a
-  // postprocessor. if we do, we know that we have type_dof_data, which makes
-  // things a bit simpler, we also don't need to deal with some of the other
-  // stuff and use a different constructor of DataEntry
-
-  AssertDimension (vec.size(), dof_handler.n_dofs());
-
-  auto new_entry = std_cxx14::make_unique<internal::DataOut::DataEntry<DoFHandlerType,VectorType>>
-                   (&dof_handler, &vec, &data_postprocessor);
-  dof_data.emplace_back (std::move(new_entry));
-}
-
-
-
-template <typename DoFHandlerType,
-          int patch_dim, int patch_space_dim>
-template <typename VectorType>
-void
-DataOut_DoFData<DoFHandlerType,patch_dim,patch_space_dim>::
-add_data_vector
-(const DoFHandlerType           &dof_handler,
- const VectorType               &data,
- const std::string              &name,
- const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation)
-{
-  const unsigned int n_components = dof_handler.get_fe().n_components ();
-
-  std::vector<std::string> names;
-  // if only one component: we only need one name
-  if (n_components == 1)
-    names.resize (1, name);
-  else
-    // otherwise append _i to the given name
-    {
-      names.resize (n_components);
-      for (unsigned int i=0; i<n_components; ++i)
-        {
-          std::ostringstream namebuf;
-          namebuf << '_' << i;
-          names[i] = name + namebuf.str();
-        }
-    }
-
-  add_data_vector (dof_handler, data, names, data_component_interpretation);
-}
-
-
-
-template <typename DoFHandlerType,
-          int patch_dim, int patch_space_dim>
-template <typename VectorType>
-void
-DataOut_DoFData<DoFHandlerType,patch_dim,patch_space_dim>::
-add_data_vector
-(const DoFHandlerType           &dof_handler,
- const VectorType               &data,
- const std::vector<std::string> &names,
- const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation_)
-{
-  // this is an extended version of the other functions where we pass a vector
-  // together with its DoFHandler. if we do, we know that we have
-  // type_dof_data, which makes things a bit simpler
-  if (triangulation == nullptr)
-    triangulation = SmartPointer<const Triangulation<DoFHandlerType::dimension,DoFHandlerType::space_dimension> >(&dof_handler.get_triangulation(), typeid(*this).name());
-
-  Assert (&dof_handler.get_triangulation() == triangulation,
-          ExcMessage("The triangulation attached to the DoFHandler does not "
-                     "match with the one set previously"));
-
-  Assert (data.size() == dof_handler.n_dofs(),
-          ExcDimensionMismatch (data.size(), dof_handler.n_dofs()));
-  Assert (names.size() == dof_handler.get_fe().n_components(),
-          Exceptions::DataOut::ExcInvalidNumberOfNames (names.size(),
-                                                        dof_handler.get_fe().n_components()));
-
-  const std::vector<DataComponentInterpretation::DataComponentInterpretation> &
-  data_component_interpretation
+  const auto &data_component_interpretation
     = (data_component_interpretation_.size() != 0
        ?
        data_component_interpretation_
        :
        std::vector<DataComponentInterpretation::DataComponentInterpretation>
-       (names.size(), DataComponentInterpretation::component_is_scalar));
+       (deduced_names.size(), DataComponentInterpretation::component_is_scalar));
 
+  // finally, add the data vector:
   auto new_entry = std_cxx14::make_unique<internal::DataOut::DataEntry<DoFHandlerType,VectorType>>
-                   (&dof_handler, &data, names, data_component_interpretation);
-  dof_data.emplace_back (std::move(new_entry));
+                   (dof_handler, &data_vector, deduced_names, data_component_interpretation);
+
+  if (actual_type == type_dof_data)
+    dof_data.emplace_back (std::move(new_entry));
+  else
+    cell_data.emplace_back (std::move(new_entry));
 }
 
 
@@ -1363,7 +1274,6 @@ DataOut_DoFData<DoFHandlerType,patch_dim,patch_space_dim>::memory_consumption ()
 }
 
 
-
 // explicit instantiations
 #include "data_out_dof_data.inst"
 
index 123cf995d38b7effd6d7360f697542cb43a0e7d9..32bf5f57fd6fc9e45251c04d66d372b70b4641b4 100644 (file)
@@ -20,36 +20,12 @@ for (VEC : REAL_SERIAL_VECTORS; DH : DOFHANDLER_TEMPLATES; deal_II_dimension : D
 
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension,deal_II_dimension>::
-    add_data_vector<VEC> (const VEC            &,
-                          const std::string   &,
-                          const DataVectorType,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension,deal_II_dimension>::
-    add_data_vector<VEC> (const VEC                       &,
-                          const std::vector<std::string> &,
-                          const DataVectorType,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension,deal_II_dimension>::
-    add_data_vector<VEC> (const VEC                 &,
-                          const DataPostprocessor<DH<deal_II_dimension,deal_II_dimension>::space_dimension> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension,deal_II_dimension>::
-    add_data_vector<VEC> (const DH<deal_II_dimension,deal_II_dimension> &,
-                          const VEC            &,
-                          const std::string   &,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension,deal_II_dimension>::
-    add_data_vector<VEC> (const DH<deal_II_dimension,deal_II_dimension> &,
-                          const VEC                       &,
-                          const std::vector<std::string> &,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
+    add_data_vector_internal<VEC> (const DH<deal_II_dimension,deal_II_dimension> *,
+                                   const VEC                       &,
+                                   const std::vector<std::string> &,
+                                   const DataVectorType,
+                                   const std::vector<DataComponentInterpretation::DataComponentInterpretation> &,
+                                   const bool);
 
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension,deal_II_dimension>::
@@ -63,36 +39,12 @@ for (VEC : REAL_SERIAL_VECTORS; DH : DOFHANDLER_TEMPLATES; deal_II_dimension : D
 
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension-1,deal_II_dimension>::
-    add_data_vector<VEC> (const VEC            &,
-                          const std::string   &,
-                          const DataVectorType,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension-1,deal_II_dimension>::
-    add_data_vector<VEC> (const VEC                       &,
-                          const std::vector<std::string> &,
-                          const DataVectorType,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension-1,deal_II_dimension>::
-    add_data_vector<VEC> (const VEC                 &,
-                          const DataPostprocessor<DH<deal_II_dimension,deal_II_dimension>::space_dimension> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension-1,deal_II_dimension>::
-    add_data_vector<VEC> (const DH<deal_II_dimension,deal_II_dimension> &,
-                          const VEC            &,
-                          const std::string   &,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension-1,deal_II_dimension>::
-    add_data_vector<VEC> (const DH<deal_II_dimension,deal_II_dimension> &,
-                          const VEC                       &,
-                          const std::vector<std::string> &,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
+    add_data_vector_internal<VEC> (const DH<deal_II_dimension,deal_II_dimension> *,
+                                   const VEC                       &,
+                                   const std::vector<std::string> &,
+                                   const DataVectorType,
+                                   const std::vector<DataComponentInterpretation::DataComponentInterpretation> &,
+                                   const bool);
 
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension-1,deal_II_dimension>::
@@ -107,36 +59,12 @@ for (VEC : REAL_SERIAL_VECTORS; DH : DOFHANDLER_TEMPLATES; deal_II_dimension : D
 #if deal_II_dimension < 3
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension+1,deal_II_dimension+1>::
-    add_data_vector<VEC> (const VEC            &,
-                          const std::string   &,
-                          const DataVectorType,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension+1,deal_II_dimension+1>::
-    add_data_vector<VEC> (const VEC                       &,
-                          const std::vector<std::string> &,
-                          const DataVectorType,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension+1,deal_II_dimension+1>::
-    add_data_vector<VEC> (const VEC                 &,
-                          const DataPostprocessor<DH<deal_II_dimension,deal_II_dimension>::space_dimension> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension+1,deal_II_dimension+1>::
-    add_data_vector<VEC> (const DH<deal_II_dimension,deal_II_dimension> &,
-                          const VEC            &,
-                          const std::string   &,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension+1,deal_II_dimension+1>::
-    add_data_vector<VEC> (const DH<deal_II_dimension,deal_II_dimension> &,
-                          const VEC                       &,
-                          const std::vector<std::string> &,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
+    add_data_vector_internal<VEC> (const DH<deal_II_dimension,deal_II_dimension> *,
+                                   const VEC                       &,
+                                   const std::vector<std::string> &,
+                                   const DataVectorType,
+                                   const std::vector<DataComponentInterpretation::DataComponentInterpretation> &,
+                                   const bool);
 
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension+1,deal_II_dimension+1>::
@@ -151,36 +79,12 @@ for (VEC : REAL_SERIAL_VECTORS; DH : DOFHANDLER_TEMPLATES; deal_II_dimension : D
 #if deal_II_dimension < 3
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension+1>,deal_II_dimension,deal_II_dimension+1>::
-    add_data_vector<VEC> (const VEC            &,
-                          const std::string   &,
-                          const DataVectorType,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension+1>,deal_II_dimension,deal_II_dimension+1>::
-    add_data_vector<VEC> (const VEC                       &,
-                          const std::vector<std::string> &,
-                          const DataVectorType,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension+1>,deal_II_dimension,deal_II_dimension+1>::
-    add_data_vector<VEC> (const VEC                 &,
-                          const DataPostprocessor<DH<deal_II_dimension,deal_II_dimension+1>::space_dimension> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension+1>,deal_II_dimension,deal_II_dimension+1>::
-    add_data_vector<VEC> (const DH<deal_II_dimension,deal_II_dimension+1> &,
-                          const VEC            &,
-                          const std::string   &,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension+1>,deal_II_dimension,deal_II_dimension+1>::
-    add_data_vector<VEC> (const DH<deal_II_dimension,deal_II_dimension+1> &,
-                          const VEC                       &,
-                          const std::vector<std::string> &,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
+    add_data_vector_internal<VEC> (const DH<deal_II_dimension,deal_II_dimension+1> *,
+                                   const VEC                       &,
+                                   const std::vector<std::string> &,
+                                   const DataVectorType,
+                                   const std::vector<DataComponentInterpretation::DataComponentInterpretation> &,
+                                   const bool);
 
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension+1>,deal_II_dimension,deal_II_dimension+1>::
@@ -196,36 +100,12 @@ for (VEC : REAL_SERIAL_VECTORS; DH : DOFHANDLER_TEMPLATES; deal_II_dimension : D
 #if deal_II_dimension == 3
     template void
     DataOut_DoFData<DH<1,3>,1,3>::
-    add_data_vector<VEC> (const VEC            &,
-                          const std::string   &,
-                          const DataVectorType,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<1,3>,1,3>::
-    add_data_vector<VEC> (const VEC                       &,
-                          const std::vector<std::string> &,
-                          const DataVectorType,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<1,3>,1,3>::
-    add_data_vector<VEC> (const VEC                 &,
-                          const DataPostprocessor<DH<1,3>::space_dimension> &);
-
-    template void
-    DataOut_DoFData<DH<1,3>,1,3>::
-    add_data_vector<VEC> (const DH<1,3> &,
-                          const VEC            &,
-                          const std::string   &,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<1,3>,1,3>::
-    add_data_vector<VEC> (const DH<1,3> &,
-                          const VEC                       &,
-                          const std::vector<std::string> &,
-                          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
+      add_data_vector_internal<VEC> (const DH<1, 3> *,
+                                     const VEC                       &,
+                                     const std::vector<std::string> &,
+                                     const DataVectorType,
+                                     const std::vector<DataComponentInterpretation::DataComponentInterpretation> &,
+                                     const bool);
 
     template void
     DataOut_DoFData<DH<1,3>,1,3>::
@@ -244,36 +124,18 @@ for (DH : DOFHANDLER_TEMPLATES; deal_II_dimension : DIMENSIONS)
 
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension,deal_II_dimension>::
-    add_data_vector<IndexSet> (const IndexSet            &,
-                               const std::string   &,
-                               const DataVectorType,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension,deal_II_dimension>::
-    add_data_vector<IndexSet> (const IndexSet                       &,
-                               const std::vector<std::string> &,
-                               const DataVectorType,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension,deal_II_dimension>::
-    add_data_vector<IndexSet> (const DH<deal_II_dimension,deal_II_dimension> &,
-                               const IndexSet            &,
-                               const std::string   &,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
+    add_data_vector_internal<IndexSet> (const DH<deal_II_dimension,deal_II_dimension> *,
+                                        const IndexSet                       &,
+                                        const std::vector<std::string> &,
+                                        const DataVectorType,
+                                        const std::vector<DataComponentInterpretation::DataComponentInterpretation> &,
+                                        const bool);
 
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension,deal_II_dimension>::
-    add_data_vector<IndexSet> (const DH<deal_II_dimension,deal_II_dimension> &,
-                               const IndexSet                       &,
-                               const std::vector<std::string> &,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension,deal_II_dimension>::
-    add_data_vector<IndexSet> (const IndexSet                 &,
-                               const DataPostprocessor<DH<deal_II_dimension,deal_II_dimension>::space_dimension> &);
+      add_data_vector<IndexSet> (const DH<deal_II_dimension, deal_II_dimension> &,
+                                 const IndexSet                 &,
+                                 const DataPostprocessor<DH<deal_II_dimension,deal_II_dimension>::space_dimension> &);
 
 
 
@@ -281,74 +143,36 @@ for (DH : DOFHANDLER_TEMPLATES; deal_II_dimension : DIMENSIONS)
 
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension-1,deal_II_dimension>::
-    add_data_vector<IndexSet> (const IndexSet            &,
-                               const std::string   &,
-                               const DataVectorType,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension-1,deal_II_dimension>::
-    add_data_vector<IndexSet> (const IndexSet                       &,
-                               const std::vector<std::string> &,
-                               const DataVectorType,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
+    add_data_vector_internal<IndexSet> (const DH<deal_II_dimension,deal_II_dimension> *,
+                                        const IndexSet                       &,
+                                        const std::vector<std::string> &,
+                                        const DataVectorType,
+                                        const std::vector<DataComponentInterpretation::DataComponentInterpretation> &,
+                                        const bool);
 
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension-1,deal_II_dimension>::
-    add_data_vector<IndexSet> (const IndexSet                 &,
+    add_data_vector<IndexSet> (const DH<deal_II_dimension, deal_II_dimension> &,
+                               const IndexSet                 &,
                                const DataPostprocessor<DH<deal_II_dimension,deal_II_dimension>::space_dimension> &);
 
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension-1,deal_II_dimension>::
-    add_data_vector<IndexSet> (const DH<deal_II_dimension,deal_II_dimension> &,
-                               const IndexSet            &,
-                               const std::string   &,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension-1,deal_II_dimension>::
-    add_data_vector<IndexSet> (const DH<deal_II_dimension,deal_II_dimension> &,
-                               const IndexSet                       &,
-                               const std::vector<std::string> &,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
 // things for DataOutRotation
 
 #if deal_II_dimension < 3
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension+1,deal_II_dimension+1>::
-    add_data_vector<IndexSet> (const IndexSet            &,
-                               const std::string   &,
-                               const DataVectorType,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
+    add_data_vector_internal<IndexSet> (const DH<deal_II_dimension,deal_II_dimension> *,
+                                        const IndexSet                       &,
+                                        const std::vector<std::string> &,
+                                        const DataVectorType,
+                                        const std::vector<DataComponentInterpretation::DataComponentInterpretation> &,
+                                        const bool);
 
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension+1,deal_II_dimension+1>::
-    add_data_vector<IndexSet> (const IndexSet                       &,
-                               const std::vector<std::string> &,
-                               const DataVectorType,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension+1,deal_II_dimension+1>::
-    add_data_vector<IndexSet> (const IndexSet                 &,
+    add_data_vector<IndexSet> (const DH<deal_II_dimension, deal_II_dimension> &,
+                               const IndexSet                 &,
                                const DataPostprocessor<DH<deal_II_dimension,deal_II_dimension>::space_dimension> &);
-
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension+1,deal_II_dimension+1>::
-    add_data_vector<IndexSet> (const DH<deal_II_dimension,deal_II_dimension> &,
-                               const IndexSet            &,
-                               const std::string   &,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension>,deal_II_dimension+1,deal_II_dimension+1>::
-    add_data_vector<IndexSet> (const DH<deal_II_dimension,deal_II_dimension> &,
-                               const IndexSet                       &,
-                               const std::vector<std::string> &,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
 #endif
 
 // codim 1
@@ -356,21 +180,17 @@ for (DH : DOFHANDLER_TEMPLATES; deal_II_dimension : DIMENSIONS)
 #if deal_II_dimension < 3
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension+1>,deal_II_dimension,deal_II_dimension+1>::
-    add_data_vector<IndexSet> (const IndexSet            &,
-                               const std::string   &,
-                               const DataVectorType,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
-
-    template void
-    DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension+1>,deal_II_dimension,deal_II_dimension+1>::
-    add_data_vector<IndexSet> (const IndexSet                       &,
-                               const std::vector<std::string> &,
-                               const DataVectorType,
-                               const std::vector<DataComponentInterpretation::DataComponentInterpretation> &);
+    add_data_vector_internal<IndexSet> (const DH<deal_II_dimension,deal_II_dimension+1> *,
+                                        const IndexSet                       &,
+                                        const std::vector<std::string> &,
+                                        const DataVectorType,
+                                        const std::vector<DataComponentInterpretation::DataComponentInterpretation> &,
+                                        const bool);
 
     template void
     DataOut_DoFData<DH<deal_II_dimension,deal_II_dimension+1>,deal_II_dimension,deal_II_dimension+1>::
-    add_data_vector<IndexSet> (const IndexSet                 &,
+    add_data_vector<IndexSet> (const DH<deal_II_dimension, deal_II_dimension+1> &,
+                               const IndexSet                 &,
                                const DataPostprocessor<DH<deal_II_dimension,deal_II_dimension+1>::space_dimension> &);
 #endif
 
index 471bff544c29dfa3a3f91bdab75b79e551ba89d2..ff48e342be70382bfe2a2d6c7e879175a1de0afe 100644 (file)
@@ -1,4 +1,4 @@
 
-DEAL::Exceptions::DataOut::ExcInvalidNumberOfNames (names.size(), dof_handler.get_fe().n_components())
-DEAL::Exceptions::DataOut::ExcInvalidNumberOfNames (names.size(), dofs->get_fe().n_components())
+DEAL::Exceptions::DataOut::ExcInvalidNumberOfNames (deduced_names.size(), dof_handler->get_fe().n_components())
+DEAL::Exceptions::DataOut::ExcInvalidNumberOfNames (deduced_names.size(), dof_handler->get_fe().n_components())
 DEAL::OK

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.