point_value_history.h

// ---------------------------------------------------------------------
//
// Copyright (C) 2009 - 2018 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.md at
// the top level directory of deal.II.
//
// ---------------------------------------------------------------------


#ifndef dealii_point_value_history_h
#define dealii_point_value_history_h

#include <deal.II/base/exceptions.h>
#include <deal.II/base/point.h>
#include <deal.II/base/quadrature_lib.h>
#include <deal.II/base/smartpointer.h>
#include <deal.II/base/utilities.h>

#include <deal.II/dofs/dof_accessor.h>
#include <deal.II/dofs/dof_handler.h>

#include <deal.II/fe/component_mask.h>
#include <deal.II/fe/fe_q.h>
#include <deal.II/fe/fe_values.h>
#include <deal.II/fe/mapping.h>
#include <deal.II/fe/mapping_q1.h>

#include <deal.II/grid/grid_tools.h>

#include <deal.II/lac/vector.h>

#include <deal.II/numerics/data_postprocessor.h>

#include <fstream>
#include <iostream>
#include <map>
#include <sstream>
#include <string>
#include <vector>

DEAL_II_NAMESPACE_OPEN

namespace internal
{
  namespace PointValueHistoryImplementation
  {
    template <int dim>
    class PointGeometryData
    {
    public:
      PointGeometryData(
        const Point<dim> &                          new_requested_location,
        const std::vector<Point<dim>> &             new_locations,
        const std::vector<types::global_dof_index> &new_sol_indices);
      Point<dim>                           requested_location;
      std::vector<Point<dim>>              support_point_locations;
      std::vector<types::global_dof_index> solution_indices;
    };
  } // namespace PointValueHistoryImplementation
} // namespace internal



template <int dim>
class PointValueHistory
{
public:
  PointValueHistory(const unsigned int n_independent_variables = 0);

  PointValueHistory(const DoFHandler<dim> &dof_handler,
                    const unsigned int     n_independent_variables = 0);

  PointValueHistory(const PointValueHistory &point_value_history);

  PointValueHistory &
  operator=(const PointValueHistory &point_value_history);

  ~PointValueHistory();

  void
  add_point(const Point<dim> &location);

  void
  add_points(const std::vector<Point<dim>> &locations);



  void
  add_field_name(const std::string &  vector_name,
                 const ComponentMask &component_mask = ComponentMask());

  void
  add_field_name(const std::string &vector_name,
                 const unsigned int n_components);

  void
  add_component_names(const std::string &             vector_name,
                      const std::vector<std::string> &component_names);

  void
  add_independent_names(const std::vector<std::string> &independent_names);



  template <class VectorType>
  void
  evaluate_field(const std::string &name, const VectorType &solution);


  template <class VectorType>
  void
  evaluate_field(const std::vector<std::string> &names,
                 const VectorType &              solution,
                 const DataPostprocessor<dim> &  data_postprocessor,
                 const Quadrature<dim> &         quadrature);

  template <class VectorType>
  void
  evaluate_field(const std::string &           name,
                 const VectorType &            solution,
                 const DataPostprocessor<dim> &data_postprocessor,
                 const Quadrature<dim> &       quadrature);


  template <class VectorType>
  void
  evaluate_field_at_requested_location(const std::string &name,
                                       const VectorType & solution);


  void
  start_new_dataset(const double key);

  void
  push_back_independent(const std::vector<double> &independent_values);


  void
  write_gnuplot(const std::string &            base_name,
                const std::vector<Point<dim>> &postprocessor_locations =
                  std::vector<Point<dim>>());


  Vector<double>
  mark_support_locations();

  void
  get_support_locations(std::vector<std::vector<Point<dim>>> &locations);

  void
  get_points(std::vector<std::vector<Point<dim>>> &locations);

  void
  get_postprocessor_locations(const Quadrature<dim> &  quadrature,
                              std::vector<Point<dim>> &locations);

  void
  close();


  void
  clear();

  void
  status(std::ostream &out);


  bool
  deep_check(const bool strict);

  DeclExceptionMsg(ExcNoIndependent,
                   "A call has been made to push_back_independent() when "
                   "no independent values were requested.");

  DeclExceptionMsg(
    ExcDataLostSync,
    "This error is thrown to indicate that the data sets appear to be out of "
    "sync. The class requires that the number of dataset keys is the same as "
    "the number of independent values sets and mesh linked value sets. The "
    "number of each of these is allowed to differ by one to allow new values "
    "to be added with out restricting the order the user choses to do so. "
    "Special cases of no FHandler and no independent values should not "
    "trigger this error.");


  DeclExceptionMsg(
    ExcDoFHandlerRequired,
    "A method which requires access to a @p DoFHandler to be meaningful has "
    "been called when have_dof_handler is false (most likely due to default "
    "constructor being called). Only independent variables may be logged with "
    "no DoFHandler.");

  DeclExceptionMsg(
    ExcDoFHandlerChanged,
    "The triangulation has been refined or coarsened in some way. This "
    "suggests that the internal DoF indices stored by the current "
    "object are no longer meaningful.");

private:
  std::vector<double> dataset_key;

  std::vector<std::vector<double>> independent_values;

  std::vector<std::string> indep_names;

  std::map<std::string, std::vector<std::vector<double>>> data_store;

  std::map<std::string, ComponentMask> component_mask;


  std::map<std::string, std::vector<std::string>> component_names_map;

  std::vector<internal::PointValueHistoryImplementation::PointGeometryData<dim>>
    point_geometry_data;


  bool closed;

  bool cleared;


  SmartPointer<const DoFHandler<dim>, PointValueHistory<dim>> dof_handler;


  bool triangulation_changed;

  bool have_dof_handler;

  boost::signals2::connection tria_listener;

  unsigned int n_indep;


  void
  tria_change_listener();
};


DEAL_II_NAMESPACE_CLOSE
#endif /* dealii_point_value_history_h */