Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
point_value_history.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2009 - 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 
17 #ifndef dealii_point_value_history_h
18 #define dealii_point_value_history_h
19 
20 #include <deal.II/base/exceptions.h>
21 #include <deal.II/base/point.h>
22 #include <deal.II/base/quadrature_lib.h>
23 #include <deal.II/base/smartpointer.h>
24 #include <deal.II/base/utilities.h>
25 
26 #include <deal.II/dofs/dof_accessor.h>
27 #include <deal.II/dofs/dof_handler.h>
28 
29 #include <deal.II/fe/component_mask.h>
30 #include <deal.II/fe/fe_q.h>
31 #include <deal.II/fe/fe_values.h>
32 #include <deal.II/fe/mapping.h>
33 #include <deal.II/fe/mapping_q1.h>
34 
35 #include <deal.II/grid/grid_tools.h>
36 
37 #include <deal.II/lac/vector.h>
38 
39 #include <deal.II/numerics/data_postprocessor.h>
40 
41 #include <fstream>
42 #include <iostream>
43 #include <map>
44 #include <sstream>
45 #include <string>
46 #include <vector>
47 
48 DEAL_II_NAMESPACE_OPEN
49 
50 namespace internal
51 {
52  namespace PointValueHistoryImplementation
53  {
58  template <int dim>
60  {
61  public:
63  const Point<dim> & new_requested_location,
64  const std::vector<Point<dim>> & new_locations,
65  const std::vector<types::global_dof_index> &new_sol_indices);
66  Point<dim> requested_location;
67  std::vector<Point<dim>> support_point_locations;
68  std::vector<types::global_dof_index> solution_indices;
69  };
70  } // namespace PointValueHistoryImplementation
71 } // namespace internal
72 
73 
74 
211 template <int dim>
213 {
214 public:
220  PointValueHistory(const unsigned int n_independent_variables = 0);
221 
237  const unsigned int n_independent_variables = 0);
238 
244  PointValueHistory(const PointValueHistory &point_value_history);
245 
254  operator=(const PointValueHistory &point_value_history);
255 
260 
268  void
269  add_point(const Point<dim> &location);
270 
281  void
282  add_points(const std::vector<Point<dim>> &locations);
283 
284 
285 
293  void
294  add_field_name(const std::string & vector_name,
296 
305  void
306  add_field_name(const std::string &vector_name,
307  const unsigned int n_components);
308 
313  void
314  add_component_names(const std::string & vector_name,
315  const std::vector<std::string> &component_names);
316 
321  void
322  add_independent_names(const std::vector<std::string> &independent_names);
323 
324 
325 
334  template <class VectorType>
335  void
336  evaluate_field(const std::string &name, const VectorType &solution);
337 
338 
356  template <class VectorType>
357  void
358  evaluate_field(const std::vector<std::string> &names,
359  const VectorType & solution,
360  const DataPostprocessor<dim> & data_postprocessor,
361  const Quadrature<dim> & quadrature);
362 
368  template <class VectorType>
369  void
370  evaluate_field(const std::string & name,
371  const VectorType & solution,
372  const DataPostprocessor<dim> &data_postprocessor,
373  const Quadrature<dim> & quadrature);
374 
375 
388  template <class VectorType>
389  void
390  evaluate_field_at_requested_location(const std::string &name,
391  const VectorType & solution);
392 
393 
402  void
403  start_new_dataset(const double key);
404 
411  void
412  push_back_independent(const std::vector<double> &independent_values);
413 
414 
433  void
434  write_gnuplot(const std::string & base_name,
435  const std::vector<Point<dim>> &postprocessor_locations =
436  std::vector<Point<dim>>());
437 
438 
464 
472  void
473  get_support_locations(std::vector<std::vector<Point<dim>>> &locations);
474 
483  DEAL_II_DEPRECATED
484  void
485  get_points(std::vector<std::vector<Point<dim>>> &locations);
486 
496  void
497  get_postprocessor_locations(const Quadrature<dim> & quadrature,
498  std::vector<Point<dim>> &locations);
499 
511  void
512  close();
513 
514 
524  void
525  clear();
526 
532  void
533  status(std::ostream &out);
534 
535 
544  bool
545  deep_check(const bool strict);
546 
551  "A call has been made to push_back_independent() when "
552  "no independent values were requested.");
553 
559  "This error is thrown to indicate that the data sets appear to be out of "
560  "sync. The class requires that the number of dataset keys is the same as "
561  "the number of independent values sets and mesh linked value sets. The "
562  "number of each of these is allowed to differ by one to allow new values "
563  "to be added with out restricting the order the user choses to do so. "
564  "Special cases of no FHandler and no independent values should not "
565  "trigger this error.");
566 
567 
573  "A method which requires access to a @p DoFHandler to be meaningful has "
574  "been called when have_dof_handler is false (most likely due to default "
575  "constructor being called). Only independent variables may be logged with "
576  "no DoFHandler.");
577 
583  "The triangulation has been refined or coarsened in some way. This "
584  "suggests that the internal DoF indices stored by the current "
585  "object are no longer meaningful.");
586 
587 private:
592  std::vector<double> dataset_key;
593 
597  std::vector<std::vector<double>> independent_values;
598 
604  std::vector<std::string> indep_names;
605 
613  std::map<std::string, std::vector<std::vector<double>>> data_store;
614 
618  std::map<std::string, ComponentMask> component_mask;
619 
620 
625  std::map<std::string, std::vector<std::string>> component_names_map;
626 
630  std::vector<internal::PointValueHistoryImplementation::PointGeometryData<dim>>
632 
633 
637  bool closed;
638 
642  bool cleared;
643 
644 
650 
651 
658 
664 
668  boost::signals2::connection tria_listener;
669 
673  unsigned int n_indep;
674 
675 
683  void
685 };
686 
687 
688 DEAL_II_NAMESPACE_CLOSE
689 #endif /* dealii_point_value_history_h */
std::vector< std::vector< double > > independent_values
static ::ExceptionBase & ExcDataLostSync()
PointValueHistory(const unsigned int n_independent_variables=0)
std::vector< internal::PointValueHistoryImplementation::PointGeometryData< dim > > point_geometry_data
void get_support_locations(std::vector< std::vector< Point< dim >>> &locations)
static ::ExceptionBase & ExcDoFHandlerChanged()
std::vector< double > dataset_key
std::vector< std::string > indep_names
SmartPointer< const DoFHandler< dim >, PointValueHistory< dim > > dof_handler
static ::ExceptionBase & ExcNoIndependent()
void add_component_names(const std::string &vector_name, const std::vector< std::string > &component_names)
void add_field_name(const std::string &vector_name, const ComponentMask &component_mask=ComponentMask())
boost::signals2::connection tria_listener
void evaluate_field_at_requested_location(const std::string &name, const VectorType &solution)
PointGeometryData(const Point< dim > &new_requested_location, const std::vector< Point< dim >> &new_locations, const std::vector< types::global_dof_index > &new_sol_indices)
Only a constructor needed for this class (a struct really)
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
void get_postprocessor_locations(const Quadrature< dim > &quadrature, std::vector< Point< dim >> &locations)
std::map< std::string, std::vector< std::vector< double > > > data_store
void evaluate_field(const std::string &name, const VectorType &solution)
void start_new_dataset(const double key)
PointValueHistory & operator=(const PointValueHistory &point_value_history)
void add_independent_names(const std::vector< std::string > &independent_names)
Vector< double > mark_support_locations()
bool deep_check(const bool strict)
std::map< std::string, ComponentMask > component_mask
void write_gnuplot(const std::string &base_name, const std::vector< Point< dim >> &postprocessor_locations=std::vector< Point< dim >>())
void get_points(std::vector< std::vector< Point< dim >>> &locations)
void status(std::ostream &out)
void add_point(const Point< dim > &location)
std::map< std::string, std::vector< std::string > > component_names_map
void add_points(const std::vector< Point< dim >> &locations)
static ::ExceptionBase & ExcDoFHandlerRequired()
void push_back_independent(const std::vector< double > &independent_values)