From: David Wells <wellsd2@rpi.edu> Date: Sat, 13 May 2017 17:47:04 +0000 (-0400) Subject: Combine four DataOut_DoFData::attach_data_vector functions. X-Git-Tag: v9.0.0-rc1~1516^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=5e18902c784e54c179700fa1b4cb607c3aba6f8a;p=dealii.git Combine four DataOut_DoFData::attach_data_vector functions. --- diff --git a/include/deal.II/numerics/data_out_dof_data.h b/include/deal.II/numerics/data_out_dof_data.h index 34afa0be14..c23b710b04 100644 --- a/include/deal.II/numerics/data_out_dof_data.h +++ b/include/deal.II/numerics/data_out_dof_data.h @@ -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> diff --git a/source/numerics/data_out_dof_data.cc b/source/numerics/data_out_dof_data.cc index 3dba0c37c2..00d89adee2 100644 --- a/source/numerics/data_out_dof_data.cc +++ b/source/numerics/data_out_dof_data.cc @@ -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" diff --git a/source/numerics/data_out_dof_data.inst.in b/source/numerics/data_out_dof_data.inst.in index 123cf995d3..32bf5f57fd 100644 --- a/source/numerics/data_out_dof_data.inst.in +++ b/source/numerics/data_out_dof_data.inst.in @@ -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 diff --git a/tests/numerics/data_out_11.debug.output b/tests/numerics/data_out_11.debug.output index 471bff544c..ff48e342be 100644 --- a/tests/numerics/data_out_11.debug.output +++ b/tests/numerics/data_out_11.debug.output @@ -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