#include <deal.II/base/config.h>
#include <deal.II/base/data_out_base.h>
+#include <deal.II/base/mg_level_object.h>
#include <deal.II/base/smartpointer.h>
#include <deal.II/dofs/dof_handler.h>
const DataPostprocessor<DoFHandlerType::space_dimension>
&data_postprocessor);
+ /**
+ * Add a multilevel data vector.
+ *
+ * This function adds the vector-valued multilevel vector @p data in the
+ * form of a vector on each level that belongs to the DoFHandler @p
+ * dof_handler to the graphical output. This function is typically used in
+ * conjunction with a call to set_cell_selection() that selects cells on a
+ * specific level and not the active cells (the default).
+ *
+ * A vector @p data can be obtained in several ways, for example by using
+ * Multigrid::solution or Multigrid::defect during or after a multigrid
+ * cycle or by interpolating a solution via
+ * MGTransferMatrixFree::interpolate_to_mg().
+ *
+ * The handling of @p names and @p data_component_interpretation is identical
+ * to the add_data_vector() function.
+ */
+ template <class VectorType>
+ void
+ add_mg_data_vector(
+ const DoFHandlerType & dof_handler,
+ const MGLevelObject<VectorType> &data,
+ const std::vector<std::string> & names,
+ const std::vector<DataComponentInterpretation::DataComponentInterpretation>
+ &data_component_interpretation = std::vector<
+ DataComponentInterpretation::DataComponentInterpretation>());
+
+ /**
+ * Scalar version of the function above.
+ */
+ template <class VectorType>
+ void
+ add_mg_data_vector(const DoFHandlerType & dof_handler,
+ const MGLevelObject<VectorType> &data,
+ const std::string & name);
+
/**
* Release the pointers to the data vectors. This allows output of a new set
* of vectors without supplying the DoF handler again. Therefore, the
template <class, int, int>
friend class DataOut_DoFData;
+ /**
+ */
+ template <int, class>
+ friend class MGDataOut;
+
private:
/**
* Common function called by the four public add_data_vector methods.
if (duplicate == false)
{
- typename DoFHandlerType::active_cell_iterator dh_cell(
- &cell->get_triangulation(),
- cell->level(),
- cell->index(),
- dof_data[dataset]->dof_handler);
- if (x_fe_values.empty())
+ if (cell->active())
{
- AssertIndexRange(face, GeometryInfo<dim>::faces_per_cell);
- x_fe_face_values[dataset]->reinit(dh_cell, face);
+ typename DoFHandlerType::active_cell_iterator dh_cell(
+ &cell->get_triangulation(),
+ cell->level(),
+ cell->index(),
+ dof_data[dataset]->dof_handler);
+ if (x_fe_values.empty())
+ {
+ AssertIndexRange(face, GeometryInfo<dim>::faces_per_cell);
+ x_fe_face_values[dataset]->reinit(dh_cell, face);
+ }
+ else
+ x_fe_values[dataset]->reinit(dh_cell);
}
else
- x_fe_values[dataset]->reinit(dh_cell);
+ x_fe_values[dataset]->reinit(cell);
}
}
if (dof_data.empty())
}
+ /**
+ * Helper class templated on vector type to allow different implementations
+ * to extract information from a vector.
+ */
+ template <typename VectorType>
+ struct VectorHelper
+ {
+ /**
+ * extract the @p indices from @p vector and put them into @p values.
+ */
+ static void
+ extract(const VectorType & vector,
+ const std::vector<types::global_dof_index> &indices,
+ const ComponentExtractor extract_component,
+ std::vector<double> & values);
+ };
+
+
+
+ template <typename VectorType>
+ void
+ VectorHelper<VectorType>::extract(
+ const VectorType & vector,
+ const std::vector<types::global_dof_index> &indices,
+ const ComponentExtractor extract_component,
+ std::vector<double> & values)
+ {
+ for (unsigned int i = 0; i < values.size(); ++i)
+ values[i] = get_component(vector[indices[i]], extract_component);
+ }
+
+
+
+#ifdef DEAL_II_WITH_TRILINOS
+ template <>
+ inline void
+ VectorHelper<LinearAlgebra::EpetraWrappers::Vector>::extract(
+ const LinearAlgebra::EpetraWrappers::Vector & /*vector*/,
+ const std::vector<types::global_dof_index> & /*indices*/,
+ const ComponentExtractor /*extract_component*/,
+ std::vector<double> & /*values*/)
+ {
+ // TODO: we don't have element access
+ Assert(false, ExcNotImplemented());
+ }
+#endif
+
+
+
template <typename DoFHandlerType>
DataEntryBase<DoFHandlerType>::DataEntryBase(
const DoFHandlerType * dofs,
vector = nullptr;
this->dof_handler = nullptr;
}
+
+
+
+ /**
+ * Like DataEntry, but used to look up data from multigrid computations.
+ * Data will use level-DoF indices to look up in a
+ * MGLevelObject<VectorType> given on the specific level instead of
+ * interpolating data to coarser cells.
+ */
+ template <typename DoFHandlerType, typename VectorType>
+ class MGDataEntry : public DataEntryBase<DoFHandlerType>
+ {
+ public:
+ MGDataEntry(const DoFHandlerType * dofs,
+ const MGLevelObject<VectorType> *vectors,
+ const std::vector<std::string> & names,
+ const std::vector<
+ DataComponentInterpretation::DataComponentInterpretation>
+ &data_component_interpretation)
+ : DataEntryBase<DoFHandlerType>(dofs,
+ names,
+ data_component_interpretation)
+ , vectors(vectors)
+ {}
+
+ virtual double
+ get_cell_data_value(
+ const unsigned int cell_number,
+ const ComponentExtractor extract_component) const override;
+
+ virtual void
+ get_function_values(
+ const FEValuesBase<DoFHandlerType::dimension,
+ DoFHandlerType::space_dimension> &fe_patch_values,
+ const ComponentExtractor extract_component,
+ std::vector<double> &patch_values) const override;
+
+ /**
+ * Given a FEValuesBase object, extract the values on the present cell
+ * from the vector we actually store. This function does the same as the
+ * one above but for vector-valued finite elements.
+ */
+ virtual void
+ get_function_values(
+ const FEValuesBase<DoFHandlerType::dimension,
+ DoFHandlerType::space_dimension> &fe_patch_values,
+ const ComponentExtractor extract_component,
+ std::vector<dealii::Vector<double>> &patch_values_system)
+ const override;
+
+ /**
+ * Given a FEValuesBase object, extract the gradients on the present
+ * cell from the vector we actually store.
+ */
+ virtual void
+ get_function_gradients(
+ const FEValuesBase<DoFHandlerType::dimension,
+ DoFHandlerType::space_dimension>
+ & /*fe_patch_values*/,
+ const ComponentExtractor /*extract_component*/,
+ std::vector<Tensor<1, DoFHandlerType::space_dimension>>
+ & /*patch_gradients*/) const override
+ {
+ Assert(false, ExcNotImplemented());
+ }
+
+ /**
+ * Given a FEValuesBase object, extract the gradients on the present
+ * cell from the vector we actually store. This function does the same
+ * as the one above but for vector-valued finite elements.
+ */
+ virtual void
+ get_function_gradients(
+ const FEValuesBase<DoFHandlerType::dimension,
+ DoFHandlerType::space_dimension>
+ & /*fe_patch_values*/,
+ const ComponentExtractor /*extract_component*/,
+ std::vector<std::vector<Tensor<1, DoFHandlerType::space_dimension>>>
+ & /*patch_gradients_system*/) const override
+ {
+ Assert(false, ExcNotImplemented());
+ }
+
+
+ /**
+ * Given a FEValuesBase object, extract the second derivatives on the
+ * present cell from the vector we actually store.
+ */
+ virtual void
+ get_function_hessians(
+ const FEValuesBase<DoFHandlerType::dimension,
+ DoFHandlerType::space_dimension>
+ & /*fe_patch_values*/,
+ const ComponentExtractor /*extract_component*/,
+ std::vector<Tensor<2, DoFHandlerType::space_dimension>>
+ & /*patch_hessians*/) const override
+ {
+ Assert(false, ExcNotImplemented());
+ }
+
+ /**
+ * Given a FEValuesBase object, extract the second derivatives on the
+ * present cell from the vector we actually store. This function does
+ * the same as the one above but for vector-valued finite elements.
+ */
+ virtual void
+ get_function_hessians(
+ const FEValuesBase<DoFHandlerType::dimension,
+ DoFHandlerType::space_dimension>
+ & /*fe_patch_values*/,
+ const ComponentExtractor /*extract_component*/,
+ std::vector<std::vector<Tensor<2, DoFHandlerType::space_dimension>>>
+ & /*patch_hessians_system*/) const override
+ {
+ Assert(false, ExcNotImplemented());
+ }
+
+ /**
+ * Return whether the data represented by (a derived class of) this object
+ * represents a complex-valued (as opposed to real-valued) information.
+ */
+ virtual bool
+ is_complex_valued() const override
+ {
+ Assert(
+ numbers::NumberTraits<typename VectorType::value_type>::is_complex ==
+ false,
+ ExcNotImplemented());
+ return numbers::NumberTraits<
+ typename VectorType::value_type>::is_complex;
+ }
+
+ /**
+ * Clear all references to the vectors.
+ */
+ virtual void
+ clear() override
+ {
+ vectors = nullptr;
+ }
+
+ /**
+ * Determine an estimate for the memory consumption (in bytes) of this
+ * object.
+ */
+ virtual std::size_t
+ memory_consumption() const override
+ {
+ return sizeof(vectors);
+ }
+
+ private:
+ const MGLevelObject<VectorType> *vectors;
+ };
+
+
+
+ template <typename DoFHandlerType, typename VectorType>
+ double
+ MGDataEntry<DoFHandlerType, VectorType>::get_cell_data_value(
+ const unsigned int cell_number,
+ const ComponentExtractor extract_component) const
+ {
+ Assert(false, ExcNotImplemented());
+
+ (void)cell_number;
+ (void)extract_component;
+ return 0.0;
+ }
+
+
+
+ template <typename DoFHandlerType, typename VectorType>
+ void
+ MGDataEntry<DoFHandlerType, VectorType>::get_function_values(
+ const FEValuesBase<DoFHandlerType::dimension,
+ DoFHandlerType::space_dimension> &fe_patch_values,
+ const ComponentExtractor extract_component,
+ std::vector<double> & patch_values) const
+ {
+ Assert(extract_component == ComponentExtractor::real_part,
+ ExcNotImplemented());
+
+ const typename DoFHandlerType::level_cell_iterator dof_cell(
+ &fe_patch_values.get_cell()->get_triangulation(),
+ fe_patch_values.get_cell()->level(),
+ fe_patch_values.get_cell()->index(),
+ this->dof_handler);
+
+ const VectorType *vector = &((*vectors)[dof_cell->level()]);
+
+ const unsigned int dofs_per_cell =
+ this->dof_handler->get_fe()[0].dofs_per_cell;
+
+ std::vector<types::global_dof_index> dof_indices(dofs_per_cell);
+ dof_cell->get_mg_dof_indices(dof_indices);
+ std::vector<double> values(dofs_per_cell);
+ VectorHelper<VectorType>::extract(*vector,
+ dof_indices,
+ extract_component,
+ values);
+ std::fill(patch_values.begin(), patch_values.end(), 0.0);
+
+ for (unsigned int i = 0; i < dofs_per_cell; ++i)
+ for (unsigned int q_point = 0; q_point < patch_values.size(); ++q_point)
+ patch_values[q_point] +=
+ values[i] * fe_patch_values.shape_value(i, q_point);
+ }
+
+
+
+ template <typename DoFHandlerType, typename VectorType>
+ void
+ MGDataEntry<DoFHandlerType, VectorType>::get_function_values(
+ const FEValuesBase<DoFHandlerType::dimension,
+ DoFHandlerType::space_dimension> &fe_patch_values,
+ const ComponentExtractor extract_component,
+ std::vector<dealii::Vector<double>> &patch_values_system) const
+ {
+ Assert(extract_component == ComponentExtractor::real_part,
+ ExcNotImplemented());
+
+ typename DoFHandlerType::level_cell_iterator dof_cell(
+ &fe_patch_values.get_cell()->get_triangulation(),
+ fe_patch_values.get_cell()->level(),
+ fe_patch_values.get_cell()->index(),
+ this->dof_handler);
+
+ const VectorType *vector = &((*vectors)[dof_cell->level()]);
+
+ const unsigned int dofs_per_cell =
+ this->dof_handler->get_fe()[0].dofs_per_cell;
+
+ std::vector<types::global_dof_index> dof_indices(dofs_per_cell);
+ dof_cell->get_mg_dof_indices(dof_indices);
+ std::vector<double> values(dofs_per_cell);
+ VectorHelper<VectorType>::extract(*vector,
+ dof_indices,
+ extract_component,
+ values);
+
+ const unsigned int n_components = fe_patch_values.get_fe().n_components();
+ const unsigned int n_eval_points = fe_patch_values.n_quadrature_points;
+
+ AssertDimension(patch_values_system.size(), n_eval_points);
+ for (unsigned int q = 0; q < n_eval_points; q++)
+ {
+ AssertDimension(patch_values_system[q].size(), n_components);
+ patch_values_system[q] = 0.0;
+
+ for (unsigned int i = 0; i < dofs_per_cell; ++i)
+ for (unsigned int c = 0; c < n_components; ++c)
+ patch_values_system[q](c) +=
+ values[i] * fe_patch_values.shape_value_component(i, q, c);
+ }
+ }
+
} // namespace DataOutImplementation
} // namespace internal
+template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
+template <class VectorType>
+void
+DataOut_DoFData<DoFHandlerType, patch_dim, patch_space_dim>::add_mg_data_vector(
+ const DoFHandlerType & dof_handler,
+ const MGLevelObject<VectorType> &data,
+ const std::string & name)
+{
+ // forward the call to the vector version:
+ std::vector<std::string> names(1, name);
+ add_mg_data_vector(dof_handler, data, names);
+}
+
+
+
+template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
+template <class VectorType>
+void
+DataOut_DoFData<DoFHandlerType, patch_dim, patch_space_dim>::add_mg_data_vector(
+ const DoFHandlerType & dof_handler,
+ const MGLevelObject<VectorType> &data,
+ const std::vector<std::string> & names,
+ const std::vector<DataComponentInterpretation::DataComponentInterpretation>
+ &data_component_interpretation_)
+{
+ 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"));
+
+ const unsigned int n_components = dof_handler.get_fe(0).n_components();
+ std::vector<std::string> deduced_names = names;
+
+ if (names.size() == 1 && n_components > 1)
+ {
+ deduced_names.resize(n_components);
+ for (unsigned int i = 0; i < n_components; ++i)
+ {
+ deduced_names[i] = names[0] + '_' + std::to_string(i);
+ }
+ }
+
+ Assert(deduced_names.size() == n_components,
+ ExcMessage("Invalid number of names given."));
+
+ const std::vector<DataComponentInterpretation::DataComponentInterpretation>
+ &data_component_interpretation =
+ (data_component_interpretation_.size() != 0 ?
+ data_component_interpretation_ :
+ std::vector<DataComponentInterpretation::DataComponentInterpretation>(
+ n_components, DataComponentInterpretation::component_is_scalar));
+
+ Assert(data_component_interpretation.size() == n_components,
+ ExcMessage(
+ "Invalid number of entries in data_component_interpretation."));
+
+ auto new_entry = std_cxx14::make_unique<
+ internal::DataOutImplementation::MGDataEntry<DoFHandlerType, VectorType>>(
+ &dof_handler, &data, deduced_names, data_component_interpretation);
+ dof_data.emplace_back(std::move(new_entry));
+}
+
+
+
template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
void
DataOut_DoFData<DoFHandlerType, patch_dim, patch_space_dim>::
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_mg_data_vector<VEC>(const DH<deal_II_dimension, deal_II_dimension> &,
+ const MGLevelObject<VEC> &,
+ const std::string &);
+
+ template void DataOut_DoFData<DH<deal_II_dimension, deal_II_dimension>,
+ deal_II_dimension,
+ deal_II_dimension>::
+ add_mg_data_vector<VEC>(
+ const DH<deal_II_dimension, deal_II_dimension> &,
+ const MGLevelObject<VEC> &,
+ const std::vector<std::string> &,
+ const std::vector<
+ DataComponentInterpretation::DataComponentInterpretation> &);
+
+
// stuff needed for face data
deallog << "* owned on level " << level << std::endl;
data_out.set_cell_selection(
[level](const typename Triangulation<dim>::cell_iterator &cell) {
- return (cell->level() == level && cell->is_locally_owned_on_level());
+ return (cell->level() == static_cast<int>(level) &&
+ cell->is_locally_owned_on_level());
});
print(data_out, triangulation);
}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2006 - 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// Test DataOut::add_mg_data_vector for vector data in serial
+
+#include <deal.II/base/function_lib.h>
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/timer.h>
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+
+#include <deal.II/multigrid/mg_tools.h>
+#include <deal.II/multigrid/mg_transfer.h>
+#include <deal.II/multigrid/mg_transfer_matrix_free.h>
+
+#include <deal.II/numerics/data_out.h>
+#include <deal.II/numerics/vector_tools.h>
+
+#include "../tests.h"
+
+using namespace dealii;
+
+template <int dim>
+void
+do_test()
+{
+ FE_Q<dim> fe(1);
+ Triangulation<dim> triangulation(
+ Triangulation<dim>::limit_level_difference_at_vertices);
+ GridGenerator::hyper_cube(triangulation);
+ triangulation.refine_global(1);
+ triangulation.begin_active()->set_refine_flag();
+ triangulation.execute_coarsening_and_refinement();
+
+ DoFHandler<dim> dof_handler(triangulation);
+ dof_handler.distribute_dofs(fe);
+ dof_handler.distribute_mg_dofs(fe);
+
+ Vector<double> global_dof_vector(dof_handler.n_dofs());
+ VectorTools::interpolate(dof_handler,
+ Functions::SquareFunction<dim>(),
+ global_dof_vector);
+
+ {
+ DataOut<dim> data_out;
+ data_out.attach_dof_handler(dof_handler);
+ data_out.add_data_vector(dof_handler,
+ global_dof_vector,
+ std::vector<std::string>(1, "data"));
+ data_out.build_patches(0);
+ data_out.write_gnuplot(deallog.get_file_stream());
+ std::ofstream st("base.vtu");
+ data_out.write_vtu(st);
+ }
+
+ {
+ MGLevelObject<Vector<double>> dof_vector(0,
+ triangulation.n_global_levels() -
+ 1);
+
+ MGTransferPrebuilt<Vector<double>> transfer;
+
+ transfer.build_matrices(dof_handler);
+
+ // This is not quite the correct thing to do, but "interpolate_to_mg"
+ // doesn't exist for serial vectors. Here it means that the level 0 is 0
+ // instead of having correct values.
+ transfer.copy_to_mg(dof_handler, dof_vector, global_dof_vector);
+
+ for (unsigned int level = 0; level < triangulation.n_global_levels();
+ ++level)
+ {
+ deallog << "* level " << level << std::endl;
+
+ DataOut<dim> data_out;
+ data_out.set_cell_selection(
+ [level](const typename Triangulation<dim>::cell_iterator &cell) {
+ return cell->level() == static_cast<int>(level);
+ });
+
+ data_out.attach_triangulation(triangulation);
+ data_out.add_mg_data_vector(dof_handler, dof_vector, "data");
+ data_out.build_patches(0);
+ data_out.write_gnuplot(deallog.get_file_stream());
+ std::ofstream st(std::string("level") + Utilities::to_string(level) +
+ ".vtu");
+ data_out.write_vtu(st);
+ }
+ }
+}
+
+
+int
+main(int argc, char **argv)
+{
+ Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1);
+ initlog();
+
+ do_test<2>();
+ /*
+ do_test<3>();*/
+ return 0;
+}
--- /dev/null
+
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data>
+0.500000 0.00000 0.250000
+1.00000 0.00000 1.00000
+
+0.500000 0.500000 0.500000
+1.00000 0.500000 1.25000
+
+
+0.00000 0.500000 0.250000
+0.500000 0.500000 0.500000
+
+0.00000 1.00000 1.00000
+0.500000 1.00000 1.25000
+
+
+0.500000 0.500000 0.500000
+1.00000 0.500000 1.25000
+
+0.500000 1.00000 1.25000
+1.00000 1.00000 2.00000
+
+
+0.00000 0.00000 0.00000
+0.250000 0.00000 0.0625000
+
+0.00000 0.250000 0.0625000
+0.250000 0.250000 0.125000
+
+
+0.250000 0.00000 0.0625000
+0.500000 0.00000 0.250000
+
+0.250000 0.250000 0.125000
+0.500000 0.250000 0.312500
+
+
+0.00000 0.250000 0.0625000
+0.250000 0.250000 0.125000
+
+0.00000 0.500000 0.250000
+0.250000 0.500000 0.312500
+
+
+0.250000 0.250000 0.125000
+0.500000 0.250000 0.312500
+
+0.250000 0.500000 0.312500
+0.500000 0.500000 0.500000
+
+
+DEAL::* level 0
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data>
+0.00000 0.00000 0.00000
+1.00000 0.00000 0.00000
+
+0.00000 1.00000 0.00000
+1.00000 1.00000 0.00000
+
+
+DEAL::* level 1
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data>
+0.00000 0.00000 0.00000
+0.500000 0.00000 0.250000
+
+0.00000 0.500000 0.250000
+0.500000 0.500000 0.500000
+
+
+0.500000 0.00000 0.250000
+1.00000 0.00000 1.00000
+
+0.500000 0.500000 0.500000
+1.00000 0.500000 1.25000
+
+
+0.00000 0.500000 0.250000
+0.500000 0.500000 0.500000
+
+0.00000 1.00000 1.00000
+0.500000 1.00000 1.25000
+
+
+0.500000 0.500000 0.500000
+1.00000 0.500000 1.25000
+
+0.500000 1.00000 1.25000
+1.00000 1.00000 2.00000
+
+
+DEAL::* level 2
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data>
+0.00000 0.00000 0.00000
+0.250000 0.00000 0.0625000
+
+0.00000 0.250000 0.0625000
+0.250000 0.250000 0.125000
+
+
+0.250000 0.00000 0.0625000
+0.500000 0.00000 0.250000
+
+0.250000 0.250000 0.125000
+0.500000 0.250000 0.312500
+
+
+0.00000 0.250000 0.0625000
+0.250000 0.250000 0.125000
+
+0.00000 0.500000 0.250000
+0.250000 0.500000 0.312500
+
+
+0.250000 0.250000 0.125000
+0.500000 0.250000 0.312500
+
+0.250000 0.500000 0.312500
+0.500000 0.500000 0.500000
+
+
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2006 - 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// Test DataOut::add_mg_data_vector in parallel
+
+#include <deal.II/base/function_lib.h>
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+
+#include <deal.II/multigrid/mg_tools.h>
+#include <deal.II/multigrid/mg_transfer.h>
+#include <deal.II/multigrid/mg_transfer_matrix_free.h>
+
+#include <deal.II/numerics/data_out.h>
+#include <deal.II/numerics/vector_tools.h>
+
+#include "../tests.h"
+
+using namespace dealii;
+
+template <int dim>
+void
+do_test()
+{
+ parallel::distributed::Triangulation<dim> triangulation(
+ MPI_COMM_WORLD,
+ dealii::Triangulation<dim>::none,
+ parallel::distributed::Triangulation<dim>::construct_multigrid_hierarchy);
+
+ GridGenerator::hyper_cube(triangulation);
+ triangulation.refine_global(1);
+ triangulation.begin_active()->set_refine_flag();
+ triangulation.execute_coarsening_and_refinement();
+
+ FE_Q<dim> fe(1);
+ DoFHandler<dim> dof_handler(triangulation);
+ dof_handler.distribute_dofs(fe);
+ dof_handler.distribute_mg_dofs(fe);
+
+ // Make FE vector
+ IndexSet locally_relevant_dofs;
+ DoFTools::extract_locally_relevant_dofs(dof_handler, locally_relevant_dofs);
+
+ using VectorType = LinearAlgebra::distributed::Vector<double>;
+ VectorType global_dof_vector;
+ global_dof_vector.reinit(dof_handler.locally_owned_dofs(),
+ locally_relevant_dofs,
+ MPI_COMM_WORLD);
+
+ VectorTools::interpolate(dof_handler,
+ Functions::SquareFunction<dim>(),
+ global_dof_vector);
+ global_dof_vector.update_ghost_values();
+
+ {
+ DataOut<dim> data_out;
+ data_out.attach_dof_handler(dof_handler);
+ data_out.add_data_vector(dof_handler,
+ global_dof_vector,
+ std::vector<std::string>(1, "data"));
+ data_out.build_patches(0);
+ data_out.write_gnuplot(deallog.get_file_stream());
+ data_out.write_vtu_in_parallel("base.vtu", MPI_COMM_WORLD);
+ }
+
+ {
+ MGLevelObject<VectorType> dof_vector(0,
+ triangulation.n_global_levels() - 1);
+ MGTransferMatrixFree<dim, double> transfer;
+
+ transfer.build(dof_handler);
+ transfer.interpolate_to_mg(dof_handler, dof_vector, global_dof_vector);
+
+ for (unsigned int level = 0; level < triangulation.n_global_levels();
+ ++level)
+ {
+ deallog << "* level " << level << std::endl;
+ DataOut<dim> data_out;
+ data_out.set_cell_selection(
+ [level](const typename Triangulation<dim>::cell_iterator &cell) {
+ return (cell->level() == static_cast<int>(level) &&
+ cell->is_locally_owned_on_level());
+ });
+ data_out.attach_triangulation(triangulation);
+ data_out.add_mg_data_vector(dof_handler, dof_vector, "data");
+ data_out.build_patches(0);
+ data_out.write_gnuplot(deallog.get_file_stream());
+ data_out.write_vtu_in_parallel(std::string("level") +
+ Utilities::to_string(level) + ".vtu",
+ MPI_COMM_WORLD);
+ }
+ }
+}
+
+
+int
+main(int argc, char **argv)
+{
+ Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1);
+ MPILogInitAll log;
+
+ do_test<2>();
+ // do_test<3>();
+ return 0;
+}
--- /dev/null
+
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data>
+0.00000 0.00000 0.00000
+0.250000 0.00000 0.0625000
+
+0.00000 0.250000 0.0625000
+0.250000 0.250000 0.125000
+
+
+0.250000 0.00000 0.0625000
+0.500000 0.00000 0.250000
+
+0.250000 0.250000 0.125000
+0.500000 0.250000 0.312500
+
+
+0.00000 0.250000 0.0625000
+0.250000 0.250000 0.125000
+
+0.00000 0.500000 0.250000
+0.250000 0.500000 0.312500
+
+
+0.250000 0.250000 0.125000
+0.500000 0.250000 0.312500
+
+0.250000 0.500000 0.312500
+0.500000 0.500000 0.500000
+
+
+DEAL:0::* level 0
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data>
+0.00000 0.00000 0.00000
+1.00000 0.00000 1.00000
+
+0.00000 1.00000 1.00000
+1.00000 1.00000 2.00000
+
+
+DEAL:0::* level 1
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data>
+0.00000 0.00000 0.00000
+0.500000 0.00000 0.250000
+
+0.00000 0.500000 0.250000
+0.500000 0.500000 0.500000
+
+
+DEAL:0::* level 2
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data>
+0.00000 0.00000 0.00000
+0.250000 0.00000 0.0625000
+
+0.00000 0.250000 0.0625000
+0.250000 0.250000 0.125000
+
+
+0.250000 0.00000 0.0625000
+0.500000 0.00000 0.250000
+
+0.250000 0.250000 0.125000
+0.500000 0.250000 0.312500
+
+
+0.00000 0.250000 0.0625000
+0.250000 0.250000 0.125000
+
+0.00000 0.500000 0.250000
+0.250000 0.500000 0.312500
+
+
+0.250000 0.250000 0.125000
+0.500000 0.250000 0.312500
+
+0.250000 0.500000 0.312500
+0.500000 0.500000 0.500000
+
+
+
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data>
+0.500000 0.00000 0.250000
+1.00000 0.00000 1.00000
+
+0.500000 0.500000 0.500000
+1.00000 0.500000 1.25000
+
+
+0.00000 0.500000 0.250000
+0.500000 0.500000 0.500000
+
+0.00000 1.00000 1.00000
+0.500000 1.00000 1.25000
+
+
+0.500000 0.500000 0.500000
+1.00000 0.500000 1.25000
+
+0.500000 1.00000 1.25000
+1.00000 1.00000 2.00000
+
+
+DEAL:1::* level 0
+DEAL:1::* level 1
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data>
+0.500000 0.00000 0.250000
+1.00000 0.00000 1.00000
+
+0.500000 0.500000 0.500000
+1.00000 0.500000 1.25000
+
+
+0.00000 0.500000 0.250000
+0.500000 0.500000 0.500000
+
+0.00000 1.00000 1.00000
+0.500000 1.00000 1.25000
+
+
+0.500000 0.500000 0.500000
+1.00000 0.500000 1.25000
+
+0.500000 1.00000 1.25000
+1.00000 1.00000 2.00000
+
+
+DEAL:1::* level 2
+
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2006 - 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// Test DataOut::add_mg_data_vector in parallel for a vector valued dof function
+
+#include <deal.II/base/function_lib.h>
+#include <deal.II/base/function_parser.h>
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_system.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+
+#include <deal.II/multigrid/mg_tools.h>
+#include <deal.II/multigrid/mg_transfer.h>
+#include <deal.II/multigrid/mg_transfer_matrix_free.h>
+
+#include <deal.II/numerics/data_out.h>
+#include <deal.II/numerics/vector_tools.h>
+
+#include "../tests.h"
+
+using namespace dealii;
+
+template <int dim>
+void
+do_test()
+{
+ parallel::distributed::Triangulation<dim> triangulation(
+ MPI_COMM_WORLD,
+ dealii::Triangulation<dim>::none,
+ parallel::distributed::Triangulation<dim>::construct_multigrid_hierarchy);
+
+ GridGenerator::hyper_cube(triangulation);
+ triangulation.refine_global(1);
+ triangulation.begin_active()->set_refine_flag();
+ triangulation.execute_coarsening_and_refinement();
+
+ FESystem<dim> fe(FE_DGQ<dim>(1), dim + 1);
+
+ DoFHandler<dim> dof_handler(triangulation);
+ dof_handler.distribute_dofs(fe);
+ dof_handler.distribute_mg_dofs(fe);
+
+ // Make FE vector
+ IndexSet locally_relevant_dofs;
+ DoFTools::extract_locally_relevant_dofs(dof_handler, locally_relevant_dofs);
+
+ using VectorType = LinearAlgebra::distributed::Vector<double>;
+ VectorType global_dof_vector;
+ global_dof_vector.reinit(dof_handler.locally_owned_dofs(),
+ locally_relevant_dofs,
+ MPI_COMM_WORLD);
+
+ QGauss<dim> quadrature(3);
+
+ FunctionParser<dim> function(dim == 2 ? "y;-x;x*x+y*y" :
+ "y;-x;0;x*x+y*y+z*z");
+
+ AffineConstraints<double> constraints;
+ DoFTools::make_hanging_node_constraints(dof_handler, constraints);
+ constraints.close();
+ VectorTools::project(
+ dof_handler, constraints, quadrature, function, global_dof_vector);
+ /* VectorTools::interpolate(dof_handler,
+ Functions::SquareFunction<dim>(),
+ global_dof_vector);
+ */
+
+ global_dof_vector.update_ghost_values();
+
+ std::vector<std::string> names(dim, "vec");
+ names.emplace_back("my_scalar");
+ std::vector<DataComponentInterpretation::DataComponentInterpretation>
+ interpretation(dim,
+ DataComponentInterpretation::component_is_part_of_vector);
+ interpretation.emplace_back(DataComponentInterpretation::component_is_scalar);
+
+ {
+ DataOut<dim> data_out;
+ data_out.attach_dof_handler(dof_handler);
+ data_out.add_data_vector(dof_handler,
+ global_dof_vector,
+ names,
+ interpretation);
+ data_out.build_patches(0);
+ data_out.write_gnuplot(deallog.get_file_stream());
+ data_out.write_vtu_in_parallel("base.vtu", MPI_COMM_WORLD);
+ }
+
+ {
+ MGLevelObject<VectorType> dof_vector(0,
+ triangulation.n_global_levels() - 1);
+ MGTransferMatrixFree<dim, double> transfer;
+
+ transfer.build(dof_handler);
+ transfer.interpolate_to_mg(dof_handler, dof_vector, global_dof_vector);
+
+ for (unsigned int level = 0; level < triangulation.n_global_levels();
+ ++level)
+ {
+ deallog << "* level " << level << std::endl;
+ DataOut<dim> data_out;
+ data_out.set_cell_selection(
+ [level](const typename Triangulation<dim>::cell_iterator &cell) {
+ return (cell->level() == static_cast<int>(level) &&
+ cell->is_locally_owned_on_level());
+ });
+ data_out.attach_triangulation(triangulation);
+ data_out.add_mg_data_vector(dof_handler, dof_vector, "data");
+ data_out.add_mg_data_vector(dof_handler,
+ dof_vector,
+ names,
+ interpretation);
+ data_out.build_patches(0);
+ data_out.write_gnuplot(deallog.get_file_stream());
+ data_out.write_vtu_in_parallel(std::string("level") +
+ Utilities::to_string(level) + ".vtu",
+ MPI_COMM_WORLD);
+ }
+ }
+}
+
+
+int
+main(int argc, char **argv)
+{
+ Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1);
+ MPILogInitAll log;
+
+ do_test<2>();
+ // do_test<3>();
+ return 0;
+}
--- /dev/null
+
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <vec> <vec> <my_scalar>
+0.00000 0.00000 -6.93889e-17 2.77556e-17 -0.0208333
+0.250000 0.00000 -5.55112e-17 -0.250000 0.0416667
+
+0.00000 0.250000 0.250000 2.77556e-17 0.0416667
+0.250000 0.250000 0.250000 -0.250000 0.104167
+
+
+0.250000 0.00000 -6.93889e-17 -0.250000 0.0416667
+0.500000 0.00000 -5.55112e-17 -0.500000 0.229167
+
+0.250000 0.250000 0.250000 -0.250000 0.104167
+0.500000 0.250000 0.250000 -0.500000 0.291667
+
+
+0.00000 0.250000 0.250000 2.77556e-17 0.0416667
+0.250000 0.250000 0.250000 -0.250000 0.104167
+
+0.00000 0.500000 0.500000 2.77556e-17 0.229167
+0.250000 0.500000 0.500000 -0.250000 0.291667
+
+
+0.250000 0.250000 0.250000 -0.250000 0.104167
+0.500000 0.250000 0.250000 -0.500000 0.291667
+
+0.250000 0.500000 0.500000 -0.250000 0.291667
+0.500000 0.500000 0.500000 -0.500000 0.479167
+
+
+DEAL:0::* level 0
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data_0> <data_1> <data_2> <vec> <vec> <my_scalar>
+0.00000 0.00000 -2.35922e-16 1.24900e-16 -0.333333 -2.35922e-16 1.24900e-16 -0.333333
+1.00000 0.00000 5.55112e-17 -1.00000 0.666667 5.55112e-17 -1.00000 0.666667
+
+0.00000 1.00000 1.00000 1.66533e-16 0.666667 1.00000 1.66533e-16 0.666667
+1.00000 1.00000 1.00000 -1.00000 1.66667 1.00000 -1.00000 1.66667
+
+
+DEAL:0::* level 1
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data_0> <data_1> <data_2> <vec> <vec> <my_scalar>
+0.00000 0.00000 -3.46945e-17 6.93889e-18 -0.0833333 -3.46945e-17 6.93889e-18 -0.0833333
+0.500000 0.00000 0.00000 -0.500000 0.166667 0.00000 -0.500000 0.166667
+
+0.00000 0.500000 0.500000 1.11022e-16 0.166667 0.500000 1.11022e-16 0.166667
+0.500000 0.500000 0.500000 -0.500000 0.416667 0.500000 -0.500000 0.416667
+
+
+DEAL:0::* level 2
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data_0> <data_1> <data_2> <vec> <vec> <my_scalar>
+0.00000 0.00000 -6.93889e-17 2.77556e-17 -0.0208333 -6.93889e-17 2.77556e-17 -0.0208333
+0.250000 0.00000 -5.55112e-17 -0.250000 0.0416667 -5.55112e-17 -0.250000 0.0416667
+
+0.00000 0.250000 0.250000 2.77556e-17 0.0416667 0.250000 2.77556e-17 0.0416667
+0.250000 0.250000 0.250000 -0.250000 0.104167 0.250000 -0.250000 0.104167
+
+
+0.250000 0.00000 -6.93889e-17 -0.250000 0.0416667 -6.93889e-17 -0.250000 0.0416667
+0.500000 0.00000 -5.55112e-17 -0.500000 0.229167 -5.55112e-17 -0.500000 0.229167
+
+0.250000 0.250000 0.250000 -0.250000 0.104167 0.250000 -0.250000 0.104167
+0.500000 0.250000 0.250000 -0.500000 0.291667 0.250000 -0.500000 0.291667
+
+
+0.00000 0.250000 0.250000 2.77556e-17 0.0416667 0.250000 2.77556e-17 0.0416667
+0.250000 0.250000 0.250000 -0.250000 0.104167 0.250000 -0.250000 0.104167
+
+0.00000 0.500000 0.500000 2.77556e-17 0.229167 0.500000 2.77556e-17 0.229167
+0.250000 0.500000 0.500000 -0.250000 0.291667 0.500000 -0.250000 0.291667
+
+
+0.250000 0.250000 0.250000 -0.250000 0.104167 0.250000 -0.250000 0.104167
+0.500000 0.250000 0.250000 -0.500000 0.291667 0.250000 -0.500000 0.291667
+
+0.250000 0.500000 0.500000 -0.250000 0.291667 0.500000 -0.250000 0.291667
+0.500000 0.500000 0.500000 -0.500000 0.479167 0.500000 -0.500000 0.479167
+
+
+
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <vec> <vec> <my_scalar>
+0.500000 0.00000 -1.38778e-16 -0.500000 0.166667
+1.00000 0.00000 -1.11022e-16 -1.00000 0.916667
+
+0.500000 0.500000 0.500000 -0.500000 0.416667
+1.00000 0.500000 0.500000 -1.00000 1.16667
+
+
+0.00000 0.500000 0.500000 5.55112e-17 0.166667
+0.500000 0.500000 0.500000 -0.500000 0.416667
+
+0.00000 1.00000 1.00000 5.55112e-17 0.916667
+0.500000 1.00000 1.00000 -0.500000 1.16667
+
+
+0.500000 0.500000 0.500000 -0.500000 0.416667
+1.00000 0.500000 0.500000 -1.00000 1.16667
+
+0.500000 1.00000 1.00000 -0.500000 1.16667
+1.00000 1.00000 1.00000 -1.00000 1.91667
+
+
+DEAL:1::* level 0
+DEAL:1::* level 1
+# This file was generated by the deal.II library.
+
+
+#
+# For a description of the GNUPLOT format see the GNUPLOT manual.
+#
+# <x> <y> <data_0> <data_1> <data_2> <vec> <vec> <my_scalar>
+0.500000 0.00000 -1.38778e-16 -0.500000 0.166667 -1.38778e-16 -0.500000 0.166667
+1.00000 0.00000 -1.11022e-16 -1.00000 0.916667 -1.11022e-16 -1.00000 0.916667
+
+0.500000 0.500000 0.500000 -0.500000 0.416667 0.500000 -0.500000 0.416667
+1.00000 0.500000 0.500000 -1.00000 1.16667 0.500000 -1.00000 1.16667
+
+
+0.00000 0.500000 0.500000 5.55112e-17 0.166667 0.500000 5.55112e-17 0.166667
+0.500000 0.500000 0.500000 -0.500000 0.416667 0.500000 -0.500000 0.416667
+
+0.00000 1.00000 1.00000 5.55112e-17 0.916667 1.00000 5.55112e-17 0.916667
+0.500000 1.00000 1.00000 -0.500000 1.16667 1.00000 -0.500000 1.16667
+
+
+0.500000 0.500000 0.500000 -0.500000 0.416667 0.500000 -0.500000 0.416667
+1.00000 0.500000 0.500000 -1.00000 1.16667 0.500000 -1.00000 1.16667
+
+0.500000 1.00000 1.00000 -0.500000 1.16667 1.00000 -0.500000 1.16667
+1.00000 1.00000 1.00000 -1.00000 1.91667 1.00000 -1.00000 1.91667
+
+
+DEAL:1::* level 2
+
projects / dealii.git / commitdiff