]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove IndexSet overloads in DataOut
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Wed, 20 Mar 2019 15:20:44 +0000 (16:20 +0100)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Wed, 20 Mar 2019 15:56:45 +0000 (16:56 +0100)
doc/news/changes/incompatibilities/20190320DanielArndt [new file with mode: 0644]
examples/step-41/step-41.cc
examples/step-42/step-42.cc
include/deal.II/fe/fe_values.h
include/deal.II/numerics/data_out_dof_data.h
source/fe/fe_values.inst.in
source/numerics/data_out_dof_data.inst.in
source/numerics/data_out_dof_data_codim.inst.in

diff --git a/doc/news/changes/incompatibilities/20190320DanielArndt b/doc/news/changes/incompatibilities/20190320DanielArndt
new file mode 100644 (file)
index 0000000..e4cc414
--- /dev/null
@@ -0,0 +1,4 @@
+Removed: The IndexSet instantiations for FEValues::get_function_values()
+and DataOut_DoFData::add_data_vector() have been removed.
+<br>
+(Daniel Arndt, 2019/03/20)
index 41d2e233aeda5d63e692d7870bb697b51a0f02ec..662829c1507b96210810d02b250cdc4792f8f37d 100644 (file)
@@ -101,6 +101,7 @@ namespace Step41
     TrilinosWrappers::SparseMatrix complete_system_matrix;
 
     TrilinosWrappers::MPI::Vector solution;
+    TrilinosWrappers::MPI::Vector active_set_vector;
     TrilinosWrappers::MPI::Vector system_rhs;
     TrilinosWrappers::MPI::Vector complete_system_rhs;
     TrilinosWrappers::MPI::Vector diagonal_of_mass_matrix;
@@ -261,6 +262,7 @@ namespace Step41
 
     IndexSet solution_index_set = dof_handler.locally_owned_dofs();
     solution.reinit(solution_index_set, MPI_COMM_WORLD);
+    active_set_vector.reinit(solution_index_set, MPI_COMM_WORLD);
     system_rhs.reinit(solution_index_set, MPI_COMM_WORLD);
     complete_system_rhs.reinit(solution_index_set, MPI_COMM_WORLD);
     contact_force.reinit(solution_index_set, MPI_COMM_WORLD);
@@ -542,6 +544,10 @@ namespace Step41
                                              BoundaryValues<dim>(),
                                              constraints);
     constraints.close();
+
+    active_set_vector = 0.;
+    for (const auto index : active_set)
+      active_set_vector[index] = 1.;
   }
 
   // @sect4{ObstacleProblem::solve}
@@ -576,11 +582,7 @@ namespace Step41
   // @sect4{ObstacleProblem::output_results}
 
   // We use the vtk-format for the output.  The file contains the displacement
-  // and a numerical representation of the active set. The function looks
-  // standard but note that we can add an IndexSet object to the DataOut
-  // object in exactly the same way as a regular solution vector: it is simply
-  // interpreted as a function that is either zero (when a degree of freedom
-  // is not part of the IndexSet) or one (if it is).
+  // and a numerical representation of the active set.
   template <int dim>
   void ObstacleProblem<dim>::output_results(const unsigned int iteration) const
   {
@@ -590,7 +592,7 @@ namespace Step41
 
     data_out.attach_dof_handler(dof_handler);
     data_out.add_data_vector(solution, "displacement");
-    data_out.add_data_vector(active_set, "active_set");
+    data_out.add_data_vector(active_set_vector, "active_set");
     data_out.add_data_vector(contact_force, "lambda");
 
     data_out.build_patches();
index 94cc9e445aa904608b80d47bdaf878ecc92ff694..ec2b52df794492ca0b62c2f1b37bc5461dda8af9 100644 (file)
@@ -701,6 +701,7 @@ namespace Step42
     TrilinosWrappers::SparseMatrix newton_matrix;
 
     TrilinosWrappers::MPI::Vector solution;
+    TrilinosWrappers::MPI::Vector active_set_vector;
     TrilinosWrappers::MPI::Vector newton_rhs;
     TrilinosWrappers::MPI::Vector newton_rhs_uncondensed;
     TrilinosWrappers::MPI::Vector diag_mass_matrix_vector;
@@ -1022,6 +1023,7 @@ namespace Step42
     {
       TimerOutput::Scope t(computing_timer, "Setup: vectors");
       solution.reinit(locally_relevant_dofs, mpi_communicator);
+      active_set_vector.reinit(locally_relevant_dofs, MPI_COMM_WORLD);
       newton_rhs.reinit(locally_owned_dofs, mpi_communicator);
       newton_rhs_uncondensed.reinit(locally_owned_dofs, mpi_communicator);
       diag_mass_matrix_vector.reinit(locally_owned_dofs, mpi_communicator);
@@ -1339,6 +1341,14 @@ namespace Step42
     all_constraints.close();
     all_constraints.merge(constraints_dirichlet_and_hanging_nodes);
 
+    // We misuse distributed_solution for updating active_set_vector.
+    // We can do that because distributed_solution has a compatible underlying
+    // IndexSet and we don't need the vector anymore.
+    distributed_solution = 0.;
+    for (const auto index : active_set)
+      distributed_solution[index] = 1.;
+    active_set_vector = distributed_solution;
+
     pcout << "         Size of active set: "
           << Utilities::MPI::sum((active_set & locally_owned_dofs).n_elements(),
                                  mpi_communicator)
@@ -1388,13 +1398,9 @@ namespace Step42
 
     std::vector<types::global_dof_index> local_dof_indices(dofs_per_cell);
 
-    typename DoFHandler<dim>::active_cell_iterator cell =
-                                                     dof_handler.begin_active(),
-                                                   endc = dof_handler.end();
-
     const FEValuesExtractors::Vector displacement(0);
 
-    for (; cell != endc; ++cell)
+    for (const auto &cell : dof_handler.active_cell_iterators())
       if (cell->is_locally_owned())
         {
           fe_values.reinit(cell);
@@ -1553,10 +1559,9 @@ namespace Step42
 
     fraction_of_plastic_q_points_per_cell = 0;
 
-    typename DoFHandler<dim>::active_cell_iterator cell =
-                                                     dof_handler.begin_active(),
-                                                   endc = dof_handler.end();
-    unsigned int cell_number                            = 0;
+    auto         cell        = dof_handler.begin_active();
+    const auto   endc        = dof_handler.end();
+    unsigned int cell_number = 0;
     for (; cell != endc; ++cell, ++cell_number)
       if (cell->is_locally_owned())
         {
@@ -2048,7 +2053,7 @@ namespace Step42
                              std::vector<std::string>(dim, "contact_force"),
                              DataOut<dim>::type_dof_data,
                              data_component_interpretation);
-    data_out.add_data_vector(active_set,
+    data_out.add_data_vector(active_set_vector,
                              std::vector<std::string>(dim, "active_set"),
                              DataOut<dim>::type_dof_data,
                              data_component_interpretation);
index 713f1c64adad43fc7ffa67b2a500dc21f0f0d5fa..d810d3b6fb8f4a3aae7b1b530469ddba8ccada75 100644 (file)
@@ -2290,9 +2290,7 @@ public:
    * Vector&lt;T&gt;, BlockVector&lt;T&gt;, or one of the PETSc or Trilinos
    * vector wrapper classes. It represents a global vector of DoF values
    * associated with the DoFHandler object with which this FEValues object was
-   * last initialized. Alternatively, if the vector argument is of type
-   * IndexSet, then the function is represented as one that is either zero or
-   * one, depending on whether a DoF index is in the set or not.
+   * last initialized.
    *
    * @dealiiRequiresUpdateFlags{update_values}
    */
@@ -2449,9 +2447,7 @@ public:
    * Vector&lt;T&gt;, BlockVector&lt;T&gt;, or one of the PETSc or Trilinos
    * vector wrapper classes. It represents a global vector of DoF values
    * associated with the DoFHandler object with which this FEValues object was
-   * last initialized. Alternatively, if the vector argument is of type
-   * IndexSet, then the function is represented as one that is either zero or
-   * one, depending on whether a DoF index is in the set or not.
+   * last initialized.
    *
    * @dealiiRequiresUpdateFlags{update_gradients}
    */
@@ -2554,9 +2550,7 @@ public:
    * Vector&lt;T&gt;, BlockVector&lt;T&gt;, or one of the PETSc or Trilinos
    * vector wrapper classes. It represents a global vector of DoF values
    * associated with the DoFHandler object with which this FEValues object was
-   * last initialized. Alternatively, if the vector argument is of type
-   * IndexSet, then the function is represented as one that is either zero or
-   * one, depending on whether a DoF index is in the set or not.
+   * last initialized.
    *
    * @dealiiRequiresUpdateFlags{update_hessians}
    */
@@ -2657,9 +2651,7 @@ public:
    * Vector&lt;T&gt;, BlockVector&lt;T&gt;, or one of the PETSc or Trilinos
    * vector wrapper classes. It represents a global vector of DoF values
    * associated with the DoFHandler object with which this FEValues object was
-   * last initialized. Alternatively, if the vector argument is of type
-   * IndexSet, then the function is represented as one that is either zero or
-   * one, depending on whether a DoF index is in the set or not.
+   * last initialized.
    *
    * @dealiiRequiresUpdateFlags{update_hessians}
    */
@@ -2772,9 +2764,7 @@ public:
    * Vector&lt;T&gt;, BlockVector&lt;T&gt;, or one of the PETSc or Trilinos
    * vector wrapper classes. It represents a global vector of DoF values
    * associated with the DoFHandler object with which this FEValues object was
-   * last initialized. Alternatively, if the vector argument is of type
-   * IndexSet, then the function is represented as one that is either zero or
-   * one, depending on whether a DoF index is in the set or not.
+   * last initialized.
    *
    * @dealiiRequiresUpdateFlags{update_3rd_derivatives}
    */
index ea11144f7e3f3d5ecdced89e2fc3ecf47133ba7f..cdf385c1ccc93ad9aafe45300dfa8b4c3ffdfa93 100644 (file)
@@ -730,9 +730,7 @@ public:
    *
    * @note The actual type for the vector argument may be any vector type from
    * which FEValues can extract values on a cell using the
-   * FEValuesBase::get_function_values() function. In particular, this
-   * includes all of the usual vector types, but also IndexSet (see step-41
-   * for a use of this).
+   * FEValuesBase::get_function_values() function.
    */
   template <class VectorType>
   void
@@ -824,9 +822,7 @@ public:
    *
    * @note The actual type for the vector argument may be any vector type from
    * which FEValues can extract values on a cell using the
-   * FEValuesBase::get_function_values() function. In particular, this
-   * includes all of the usual vector types, but also IndexSet (see step-41
-   * for a use of this).
+   * FEValuesBase::get_function_values() function.
    *
    * @note The DataPostprocessor object (i.e., in reality the object of your
    * derived class) has to live until the DataOut object is destroyed as the
index 863a55f233b343bc4d2218522995432614457905..5073c4d06b22962724ac3b18c50ec68ade2242a7 100644 (file)
@@ -485,295 +485,3 @@ for (VEC : VECTOR_TYPES; deal_II_dimension : DIMENSIONS;
 
 #endif
   }
-
-
-// instantiations for VEC=IndexSet
-
-for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : DIMENSIONS)
-  {
-#if deal_II_dimension <= deal_II_space_dimension
-    template void
-    FEValuesViews::Scalar<deal_II_dimension, deal_II_space_dimension>::
-      get_function_values<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<ProductType<IndexSet::value_type, double>::type> &) const;
-    template void
-    FEValuesViews::Scalar<deal_II_dimension, deal_II_space_dimension>::
-      get_function_gradients<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<
-          ProductType<IndexSet::value_type,
-                      dealii::Tensor<1, deal_II_space_dimension>>::type> &)
-        const;
-    template void
-    FEValuesViews::Scalar<deal_II_dimension, deal_II_space_dimension>::
-      get_function_hessians<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<
-          ProductType<IndexSet::value_type,
-                      dealii::Tensor<2, deal_II_space_dimension>>::type> &)
-        const;
-    template void
-    FEValuesViews::Scalar<deal_II_dimension, deal_II_space_dimension>::
-      get_function_laplacians<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<ProductType<IndexSet::value_type, double>::type> &) const;
-    template void
-    FEValuesViews::Scalar<deal_II_dimension, deal_II_space_dimension>::
-      get_function_third_derivatives<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<
-          ProductType<IndexSet::value_type,
-                      dealii::Tensor<3, deal_II_space_dimension>>::type> &)
-        const;
-
-    template void
-    FEValuesViews::Vector<deal_II_dimension, deal_II_space_dimension>::
-      get_function_values<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<
-          ProductType<IndexSet::value_type,
-                      dealii::Tensor<1, deal_II_space_dimension>>::type> &)
-        const;
-    template void
-    FEValuesViews::Vector<deal_II_dimension, deal_II_space_dimension>::
-      get_function_gradients<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<
-          ProductType<IndexSet::value_type,
-                      dealii::Tensor<2, deal_II_space_dimension>>::type> &)
-        const;
-    template void
-    FEValuesViews::Vector<deal_II_dimension, deal_II_space_dimension>::
-      get_function_symmetric_gradients<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<ProductType<
-          IndexSet::value_type,
-          dealii::SymmetricTensor<2, deal_II_space_dimension>>::type> &) const;
-    template void
-    FEValuesViews::Vector<deal_II_dimension, deal_II_space_dimension>::
-      get_function_curls<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<ProductType<IndexSet::value_type, curl_type>::type> &)
-        const;
-    template void
-    FEValuesViews::Vector<deal_II_dimension, deal_II_space_dimension>::
-      get_function_divergences<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<ProductType<IndexSet::value_type, divergence_type>::type> &)
-        const;
-    template void
-    FEValuesViews::Vector<deal_II_dimension, deal_II_space_dimension>::
-      get_function_hessians<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<
-          ProductType<IndexSet::value_type,
-                      dealii::Tensor<3, deal_II_space_dimension>>::type> &)
-        const;
-    template void
-    FEValuesViews::Vector<deal_II_dimension, deal_II_space_dimension>::
-      get_function_laplacians<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<
-          ProductType<IndexSet::value_type,
-                      dealii::Tensor<1, deal_II_space_dimension>>::type> &)
-        const;
-    template void
-    FEValuesViews::Vector<deal_II_dimension, deal_II_space_dimension>::
-      get_function_third_derivatives<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<
-          ProductType<IndexSet::value_type,
-                      dealii::Tensor<4, deal_II_space_dimension>>::type> &)
-        const;
-
-    template void FEValuesViews::SymmetricTensor<2,
-                                                 deal_II_dimension,
-                                                 deal_II_space_dimension>::
-      get_function_values<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<ProductType<
-          IndexSet::value_type,
-          dealii::SymmetricTensor<2, deal_II_space_dimension>>::type> &) const;
-    template void FEValuesViews::
-      SymmetricTensor<2, deal_II_dimension, deal_II_space_dimension>::
-        get_function_divergences<dealii::IndexSet>(
-          const dealii::IndexSet &,
-          std::vector<
-            ProductType<IndexSet::value_type,
-                        dealii::Tensor<1, deal_II_space_dimension>>::type> &)
-          const;
-
-    template void
-    FEValuesViews::Tensor<2, deal_II_dimension, deal_II_space_dimension>::
-      get_function_values<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<
-          ProductType<IndexSet::value_type,
-                      dealii::Tensor<2, deal_II_space_dimension>>::type> &)
-        const;
-    template void
-    FEValuesViews::Tensor<2, deal_II_dimension, deal_II_space_dimension>::
-      get_function_divergences<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<
-          ProductType<IndexSet::value_type,
-                      dealii::Tensor<1, deal_II_space_dimension>>::type> &)
-        const;
-    template void
-    FEValuesViews::Tensor<2, deal_II_dimension, deal_II_space_dimension>::
-      get_function_gradients<dealii::IndexSet>(
-        const dealii::IndexSet &,
-        std::vector<
-          ProductType<IndexSet::value_type,
-                      dealii::Tensor<3, deal_II_space_dimension>>::type> &)
-        const;
-#endif
-  }
-
-
-
-// Instantiations of functions in FEValuesBase and IndexSet=IndexSet
-
-for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : DIMENSIONS)
-  {
-#if deal_II_dimension <= deal_II_space_dimension
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_values<IndexSet>(const IndexSet &,
-                                    std::vector<IndexSet::value_type> &) const;
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_values<IndexSet>(
-        const IndexSet &,
-        const ArrayView<const types::global_dof_index> &,
-        std::vector<IndexSet::value_type> &) const;
-
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_values<IndexSet>(const IndexSet &,
-                                    std::vector<Vector<IndexSet::value_type>> &)
-        const;
-
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_values<IndexSet>(
-        const IndexSet &,
-        const ArrayView<const types::global_dof_index> &,
-        std::vector<Vector<IndexSet::value_type>> &) const;
-
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_values<IndexSet>(
-        const IndexSet &,
-        const ArrayView<const types::global_dof_index> &,
-        ArrayView<std::vector<IndexSet::value_type>>,
-        bool) const;
-
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_gradients<IndexSet>(
-        const IndexSet &,
-        std::vector<
-          dealii::Tensor<1, deal_II_space_dimension, IndexSet::value_type>> &)
-        const;
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_gradients<IndexSet>(
-        const IndexSet &,
-        const ArrayView<const types::global_dof_index> &,
-        std::vector<
-          dealii::Tensor<1, deal_II_space_dimension, IndexSet::value_type>> &)
-        const;
-
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_gradients<IndexSet>(
-        const IndexSet &,
-        std::vector<std::vector<
-          dealii::Tensor<1, deal_II_space_dimension, IndexSet::value_type>>> &)
-        const;
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_gradients<IndexSet>(
-        const IndexSet &,
-        const ArrayView<const types::global_dof_index> &,
-        ArrayView<std::vector<
-          dealii::Tensor<1, deal_II_space_dimension, IndexSet::value_type>>>,
-        bool) const;
-
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_hessians<IndexSet>(
-        const IndexSet &,
-        std::vector<
-          dealii::Tensor<2, deal_II_space_dimension, IndexSet::value_type>> &)
-        const;
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_hessians<IndexSet>(
-        const IndexSet &,
-        const ArrayView<const types::global_dof_index> &,
-        std::vector<
-          dealii::Tensor<2, deal_II_space_dimension, IndexSet::value_type>> &)
-        const;
-
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_hessians<IndexSet>(
-        const IndexSet &,
-        std::vector<std::vector<
-          dealii::Tensor<2, deal_II_space_dimension, IndexSet::value_type>>> &,
-        bool) const;
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_hessians<IndexSet>(
-        const IndexSet &,
-        const ArrayView<const types::global_dof_index> &,
-        ArrayView<std::vector<
-          dealii::Tensor<2, deal_II_space_dimension, IndexSet::value_type>>>,
-        bool) const;
-
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_laplacians<IndexSet>(const IndexSet &,
-                                        std::vector<IndexSet::value_type> &)
-        const;
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_laplacians<IndexSet>(
-        const IndexSet &,
-        const ArrayView<const types::global_dof_index> &,
-        std::vector<IndexSet::value_type> &) const;
-
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_laplacians<IndexSet>(
-        const IndexSet &, std::vector<Vector<IndexSet::value_type>> &) const;
-
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_laplacians<IndexSet>(
-        const IndexSet &,
-        const ArrayView<const types::global_dof_index> &,
-        std::vector<Vector<IndexSet::value_type>> &) const;
-
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_laplacians<IndexSet>(
-        const IndexSet &,
-        const ArrayView<const types::global_dof_index> &,
-        std::vector<std::vector<IndexSet::value_type>> &,
-        bool) const;
-
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_third_derivatives<IndexSet>(
-        const IndexSet &,
-        std::vector<
-          dealii::Tensor<3, deal_II_space_dimension, IndexSet::value_type>> &)
-        const;
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_third_derivatives<IndexSet>(
-        const IndexSet &,
-        const ArrayView<const types::global_dof_index> &,
-        std::vector<
-          dealii::Tensor<3, deal_II_space_dimension, IndexSet::value_type>> &)
-        const;
-
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_third_derivatives<IndexSet>(
-        const IndexSet &,
-        std::vector<std::vector<
-          dealii::Tensor<3, deal_II_space_dimension, IndexSet::value_type>>> &,
-        bool) const;
-    template void FEValuesBase<deal_II_dimension, deal_II_space_dimension>::
-      get_function_third_derivatives<IndexSet>(
-        const IndexSet &,
-        const ArrayView<const types::global_dof_index> &,
-        ArrayView<std::vector<
-          dealii::Tensor<3, deal_II_space_dimension, IndexSet::value_type>>>,
-        bool) const;
-#endif
-  }
index 1d2cf2612df3191db355f83bc49051c1a63b6080..ecc15cb3602f3f7be9bc204495db27153cb05cb8 100644 (file)
@@ -93,82 +93,6 @@ for (VEC : VECTOR_TYPES; DH : DOFHANDLER_TEMPLATES;
 
 
 
-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_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 DataPostprocessor<
-          DH<deal_II_dimension, deal_II_dimension>::space_dimension> &);
-
-
-
-    // stuff needed for face data
-
-    template void DataOut_DoFData<DH<deal_II_dimension, deal_II_dimension>,
-                                  deal_II_dimension - 1,
-                                  deal_II_dimension>::
-      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 DH<deal_II_dimension, deal_II_dimension> &,
-        const IndexSet &,
-        const DataPostprocessor<
-          DH<deal_II_dimension, deal_II_dimension>::space_dimension> &);
-
-    // 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_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 DH<deal_II_dimension, deal_II_dimension> &,
-        const IndexSet &,
-        const DataPostprocessor<
-          DH<deal_II_dimension, deal_II_dimension>::space_dimension> &);
-#endif
-  }
-
-
-
 for (DH : DOFHANDLER_TEMPLATES; deal_II_dimension : DIMENSIONS)
   {
     template class DataOut_DoFData<DH<deal_II_dimension>, deal_II_dimension>;
index f5ae7d5b7cc224f57d1c88f50367284a40919d0b..005f78bccd763abf41c38c623e370a0920c6a6aa 100644 (file)
@@ -55,36 +55,6 @@ for (VEC : REAL_VECTOR_TYPES; DH : DOFHANDLER_TEMPLATES;
 
 
 
-for (DH : DOFHANDLER_TEMPLATES; deal_II_dimension : DIMENSIONS;
-     deal_II_space_dimension : DIMENSIONS)
-  {
-#if deal_II_dimension < deal_II_space_dimension
-    template void
-    DataOut_DoFData<DH<deal_II_dimension, deal_II_space_dimension>,
-                    deal_II_dimension,
-                    deal_II_space_dimension>::
-      add_data_vector_internal<IndexSet>(
-        const DH<deal_II_dimension, deal_II_space_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_space_dimension>,
-                    deal_II_dimension,
-                    deal_II_space_dimension>::
-      add_data_vector<IndexSet>(
-        const IndexSet &,
-        const DataPostprocessor<
-          DH<deal_II_dimension, deal_II_space_dimension>::space_dimension> &);
-#endif
-  }
-
-
-
 for (DH : DOFHANDLER_TEMPLATES; deal_II_dimension : DIMENSIONS)
   {
 #if deal_II_dimension < 3

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.