Reference documentation for deal.II version GIT b065060f03 2023-05-30 23:50:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
data_out_dof_data.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1999 - 2022 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 #ifndef dealii_data_out_dof_data_h
17 #define dealii_data_out_dof_data_h
18 
19 
20 
21 #include <deal.II/base/config.h>
22 
26 
28 
29 #include <deal.II/fe/mapping.h>
30 
31 #include <deal.II/grid/tria.h>
32 
34 #include <deal.II/hp/fe_values.h>
37 
40 
41 #include <memory>
42 
44 
45 namespace Exceptions
46 {
51  namespace DataOutImplementation
52  {
57  int,
58  << "The number of subdivisions per patch, " << arg1
59  << ", is not valid. It needs to be greater or equal to "
60  "one, or zero if you want it to be determined "
61  "automatically.");
62 
67  "For the operation you are attempting, you first need to "
68  "tell the DataOut or related object which DoFHandler or "
69  "triangulation you would like to work on. This is "
70  "generally done using the 'attach_dof_handler()' or "
71  "'attach_triangulation()' member functions.");
72 
77  "For the operation you are attempting, you first need to "
78  "tell the DataOut or related object which DoFHandler "
79  "you would like to work on. This is "
80  "generally done using the 'attach_dof_handler()' "
81  "member function.");
82 
87  int,
88  int,
89  int,
90  << "The vector has size " << arg1
91  << " but the DoFHandler object says that there are " << arg2
92  << " degrees of freedom and there are " << arg3
93  << " active cells. The size of your vector needs to be"
94  << " either equal to the number of degrees of freedom (when"
95  << " the data is of type type_dof_data), or equal to the"
96  << " number of active cells (when the data is of type "
97  << " type_cell_data).");
103  std::string,
104  size_t,
105  << "Please use only the characters [a-zA-Z0-9_<>()] for" << std::endl
106  << "description strings since some graphics formats will only accept these."
107  << std::endl
108  << "The string you gave was <" << arg1
109  << ">, within which the invalid character is <" << arg1[arg2] << ">."
110  << std::endl);
116  "When attaching a triangulation or DoFHandler object, it is "
117  "not allowed if old data vectors are still referenced. If "
118  "you want to reuse an object of the current type, you first "
119  "need to call the 'clear_data_vector()' function.");
124  int,
125  int,
126  << "You have to give one name per component in your "
127  << "data vector. The number you gave was " << arg1
128  << ", but the number of components is " << arg2 << '.');
133  "While merging sets of patches, the two sets to be merged "
134  "need to refer to data that agrees on the names of the "
135  "various variables represented. In other words, you "
136  "cannot merge sets of patches that originate from "
137  "entirely unrelated simulations.");
142  "While merging sets of patches, the two sets to be merged "
143  "need to refer to data that agrees on the number of "
144  "subdivisions and other properties. In other words, you "
145  "cannot merge sets of patches that originate from "
146  "entirely unrelated simulations.");
147 
149  int,
150  std::string,
151  << "When declaring that a number of components in a data "
152  << "set to be output logically form a vector instead of "
153  << "simply a set of scalar fields, you need to specify "
154  << "this for all relevant components. Furthermore, "
155  << "vectors must always consist of exactly <dim> "
156  << "components. However, the vector component at "
157  << "position " << arg1 << " with name <" << arg2
158  << "> does not satisfy these conditions.");
159 
161  int,
162  std::string,
163  << "When declaring that a number of components in a data "
164  << "set to be output logically form a tensor instead of "
165  << "simply a set of scalar fields, you need to specify "
166  << "this for all relevant components. Furthermore, "
167  << "tensors must always consist of exactly <dim*dim> "
168  << "components. However, the tensor component at "
169  << "position " << arg1 << " with name <" << arg2
170  << "> does not satisfy these conditions.");
171 
172  } // namespace DataOutImplementation
173 } // namespace Exceptions
174 
175 
176 namespace internal
177 {
178  namespace DataOutImplementation
179  {
204  {
205  real_part,
207  };
208 
209 
225  template <int dim, int spacedim>
227  {
228  public:
235  const std::vector<std::string> & names,
236  const std::vector<
239 
246  const DataPostprocessor<spacedim> *data_postprocessor);
247 
251  virtual ~DataEntryBase() = default;
252 
257  virtual double
258  get_cell_data_value(const unsigned int cell_number,
259  const ComponentExtractor extract_component) const = 0;
260 
265  virtual void
267  const ComponentExtractor extract_component,
268  std::vector<double> &patch_values) const = 0;
269 
275  virtual void
277  const FEValuesBase<dim, spacedim> & fe_patch_values,
278  const ComponentExtractor extract_component,
279  std::vector<::Vector<double>> &patch_values_system) const = 0;
280 
285  virtual void
287  const FEValuesBase<dim, spacedim> &fe_patch_values,
288  const ComponentExtractor extract_component,
289  std::vector<Tensor<1, spacedim>> & patch_gradients) const = 0;
290 
296  virtual void
298  const ComponentExtractor extract_component,
299  std::vector<std::vector<Tensor<1, spacedim>>>
300  &patch_gradients_system) const = 0;
301 
306  virtual void
308  const FEValuesBase<dim, spacedim> &fe_patch_values,
309  const ComponentExtractor extract_component,
310  std::vector<Tensor<2, spacedim>> & patch_hessians) const = 0;
311 
317  virtual void
319  const ComponentExtractor extract_component,
320  std::vector<std::vector<Tensor<2, spacedim>>>
321  &patch_hessians_system) const = 0;
322 
327  virtual bool
328  is_complex_valued() const = 0;
329 
333  virtual void
334  clear() = 0;
335 
340  virtual std::size_t
341  memory_consumption() const = 0;
342 
347 
351  const std::vector<std::string> names;
352 
358  const std::vector<
361 
367 
374  unsigned int n_output_variables;
375  };
376 
377 
394  template <int dim, int spacedim>
396  {
398  const unsigned int n_datasets,
399  const unsigned int n_subdivisions,
400  const std::vector<unsigned int> &n_postprocessor_outputs,
401  const Mapping<dim, spacedim> & mapping,
402  const std::vector<
403  std::shared_ptr<::hp::FECollection<dim, spacedim>>>
404  & finite_elements,
406  const bool use_face_values);
407 
409  const unsigned int n_datasets,
410  const unsigned int n_subdivisions,
411  const std::vector<unsigned int> &n_postprocessor_outputs,
412  const ::hp::MappingCollection<dim, spacedim> &mapping,
413  const std::vector<
414  std::shared_ptr<::hp::FECollection<dim, spacedim>>>
415  & finite_elements,
417  const bool use_face_values);
418 
420 
421  void
423  std::vector<std::shared_ptr<DataEntryBase<dim, spacedim>>> &dof_data,
424  const typename ::Triangulation<dim, spacedim>::cell_iterator
425  & cell,
426  const unsigned int face = numbers::invalid_unsigned_int);
427 
429  get_present_fe_values(const unsigned int dataset) const;
430 
431  void
432  resize_system_vectors(const unsigned int n_components);
433 
434  const unsigned int n_datasets;
435  const unsigned int n_subdivisions;
436 
439  std::vector<std::vector<::Vector<double>>> postprocessed_values;
440 
441  const ::hp::MappingCollection<dim, spacedim> mapping_collection;
442  const std::vector<
443  std::shared_ptr<::hp::FECollection<dim, spacedim>>>
446 
447  std::vector<std::shared_ptr<::hp::FEValues<dim, spacedim>>>
449  std::vector<std::shared_ptr<::hp::FEFaceValues<dim, spacedim>>>
451  };
452  } // namespace DataOutImplementation
453 } // namespace internal
454 
455 
456 // TODO: Most of the documentation of DataOut_DoFData applies to DataOut.
457 
593 template <int dim,
594  int patch_dim,
595  int spacedim = dim,
596  int patch_spacedim = patch_dim>
597 class DataOut_DoFData : public DataOutInterface<patch_dim, patch_spacedim>
598 {
599 public:
605 
606 public:
616  {
621 
626 
631  };
632 
637 
641  virtual ~DataOut_DoFData() override;
642 
651  void
653 
663  void
665 
729  template <class VectorType>
730  void
732  const VectorType & data,
733  const std::vector<std::string> &names,
734  const DataVectorType type = type_automatic,
735  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
736  &data_component_interpretation = {});
737 
754  template <class VectorType>
755  void
757  const VectorType & data,
758  const std::string & name,
759  const DataVectorType type = type_automatic,
760  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
761  &data_component_interpretation = {});
762 
778  template <class VectorType>
779  void
781  const DoFHandler<dim, spacedim> &dof_handler,
782  const VectorType & data,
783  const std::vector<std::string> & names,
784  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
785  &data_component_interpretation = {});
786 
787 
792  template <class VectorType>
793  void
795  const DoFHandler<dim, spacedim> &dof_handler,
796  const VectorType & data,
797  const std::string & name,
798  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
799  &data_component_interpretation = {});
800 
827  template <class VectorType>
828  void
829  add_data_vector(const VectorType & data,
830  const DataPostprocessor<spacedim> &data_postprocessor);
831 
838  template <class VectorType>
839  void
841  const VectorType & data,
842  const DataPostprocessor<spacedim> &data_postprocessor);
843 
861  template <class VectorType>
862  void
864  const DoFHandler<dim, spacedim> &dof_handler,
865  const MGLevelObject<VectorType> &data,
866  const std::vector<std::string> & names,
867  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
868  &data_component_interpretation = std::vector<
870 
874  template <class VectorType>
875  void
877  const MGLevelObject<VectorType> &data,
878  const std::string & name);
879 
886  void
888 
899  void
901 
925  template <int dim2, int spacedim2>
926  void
930 
937  virtual void
938  clear();
939 
944  std::size_t
946 
947 protected:
952 
957 
962 
966  std::vector<std::shared_ptr<
969 
973  std::vector<std::shared_ptr<
976 
982  std::vector<Patch> patches;
983 
988 public:
989  virtual const std::vector<Patch> &
990  get_patches() const override;
991 
992 protected:
997  virtual std::vector<std::string>
998  get_dataset_names() const override;
999 
1004  virtual std::vector<
1005  std::tuple<unsigned int,
1006  unsigned int,
1007  std::string,
1009  get_nonscalar_data_ranges() const override;
1010 
1011 protected:
1016  std::vector<std::shared_ptr<::hp::FECollection<dim, spacedim>>>
1017  get_fes() const;
1018 
1019  // Make all template siblings friends. Needed for the merge_patches()
1020  // function.
1021  template <int, int, int, int>
1022  friend class DataOut_DoFData;
1023 
1026  template <int, class>
1027  friend class MGDataOut;
1028 
1029 private:
1033  template <class VectorType>
1034  void
1036  const DoFHandler<dim, spacedim> *dof_handler,
1037  const VectorType & data,
1038  const std::vector<std::string> & names,
1039  const DataVectorType type,
1040  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
1041  & data_component_interpretation,
1042  const bool deduce_output_names);
1043 };
1044 
1045 
1046 
1047 // -------------------- template and inline functions ------------------------
1048 template <int dim, int patch_dim, int spacedim, int patch_spacedim>
1049 template <typename VectorType>
1050 void
1052  const VectorType & vec,
1053  const std::string & name,
1054  const DataVectorType type,
1055  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
1056  &data_component_interpretation)
1057 {
1058  Assert(triangulation != nullptr,
1060  std::vector<std::string> names(1, name);
1061  add_data_vector_internal(
1062  dofs, vec, names, type, data_component_interpretation, true);
1063 }
1064 
1065 
1066 
1067 template <int dim, int patch_dim, int spacedim, int patch_spacedim>
1068 template <typename VectorType>
1069 void
1071  const VectorType & vec,
1072  const std::vector<std::string> &names,
1073  const DataVectorType type,
1074  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
1075  &data_component_interpretation)
1076 {
1077  Assert(triangulation != nullptr,
1079  add_data_vector_internal(
1080  dofs, vec, names, type, data_component_interpretation, false);
1081 }
1082 
1083 
1084 
1085 template <int dim, int patch_dim, int spacedim, int patch_spacedim>
1086 template <typename VectorType>
1087 void
1089  const DoFHandler<dim, spacedim> &dof_handler,
1090  const VectorType & data,
1091  const std::string & name,
1092  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
1093  &data_component_interpretation)
1094 {
1095  std::vector<std::string> names(1, name);
1096  add_data_vector_internal(&dof_handler,
1097  data,
1098  names,
1099  type_dof_data,
1100  data_component_interpretation,
1101  true);
1102 }
1103 
1104 
1105 
1106 template <int dim, int patch_dim, int spacedim, int patch_spacedim>
1107 template <typename VectorType>
1108 void
1110  const DoFHandler<dim, spacedim> &dof_handler,
1111  const VectorType & data,
1112  const std::vector<std::string> & names,
1113  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
1114  &data_component_interpretation)
1115 {
1116  add_data_vector_internal(&dof_handler,
1117  data,
1118  names,
1119  type_dof_data,
1120  data_component_interpretation,
1121  false);
1122 }
1123 
1124 
1125 
1126 template <int dim, int patch_dim, int spacedim, int patch_spacedim>
1127 template <typename VectorType>
1128 void
1130  const VectorType & vec,
1131  const DataPostprocessor<spacedim> &data_postprocessor)
1132 {
1133  Assert(dofs != nullptr,
1135  add_data_vector(*dofs, vec, data_postprocessor);
1136 }
1137 
1138 
1139 
1140 template <int dim, int patch_dim, int spacedim, int patch_spacedim>
1141 template <int dim2, int spacedim2>
1142 void
1145  const Point<patch_spacedim> & shift)
1146 {
1147  const std::vector<Patch> &source_patches = source.get_patches();
1148  Assert((patches.size() != 0) && (source_patches.size() != 0),
1149  ExcMessage("When calling this function, both the current "
1150  "object and the one being merged need to have a "
1151  "nonzero number of patches associated with it. "
1152  "Either you called this function on objects that "
1153  "are empty, or you may have forgotten to call "
1154  "the 'build_patches()' function."));
1155  // Check equality of component names
1156  Assert(get_dataset_names() == source.get_dataset_names(),
1158 
1159  // Make sure patches are compatible. Ideally, we would check that all input
1160  // patches from both collections are all compatible, but we'll be content
1161  // with checking that just the first ones from both sources are.
1162  //
1163  // We check compatibility by testing that both sets of patches result
1164  // from the same number of subdivisions, and that they have the same
1165  // number of source vectors (they really should, since we already checked
1166  // that there are the same number of source components above, but you
1167  // never know). This implies that the data should have the same number of
1168  // columns. They should really have the same number of rows as well,
1169  // but depending on whether a patch has points included or not, the
1170  // number of rows may or may not include coordinates for the points,
1171  // and the comparison has to account for that because in each source
1172  // stream, the patches may include some that have points included.
1173  Assert(patches[0].n_subdivisions == source_patches[0].n_subdivisions,
1175  Assert(patches[0].data.n_cols() == source_patches[0].data.n_cols(),
1177  Assert((patches[0].data.n_rows() +
1178  (patches[0].points_are_available ? 0 : patch_spacedim)) ==
1179  (source_patches[0].data.n_rows() +
1180  (source_patches[0].points_are_available ? 0 : patch_spacedim)),
1182 
1183  // check equality of the vector data
1184  // specifications
1185  Assert(get_nonscalar_data_ranges().size() ==
1186  source.get_nonscalar_data_ranges().size(),
1187  ExcMessage("Both sources need to declare the same components "
1188  "as vectors."));
1189  for (unsigned int i = 0; i < get_nonscalar_data_ranges().size(); ++i)
1190  {
1191  Assert(std::get<0>(get_nonscalar_data_ranges()[i]) ==
1192  std::get<0>(source.get_nonscalar_data_ranges()[i]),
1193  ExcMessage("Both sources need to declare the same components "
1194  "as vectors."));
1195  Assert(std::get<1>(get_nonscalar_data_ranges()[i]) ==
1196  std::get<1>(source.get_nonscalar_data_ranges()[i]),
1197  ExcMessage("Both sources need to declare the same components "
1198  "as vectors."));
1199  Assert(std::get<2>(get_nonscalar_data_ranges()[i]) ==
1200  std::get<2>(source.get_nonscalar_data_ranges()[i]),
1201  ExcMessage("Both sources need to declare the same components "
1202  "as vectors."));
1203  }
1204 
1205  // merge patches. store old number
1206  // of elements, since we need to
1207  // adjust patch numbers, etc
1208  // afterwards
1209  const unsigned int old_n_patches = patches.size();
1210  patches.insert(patches.end(), source_patches.begin(), source_patches.end());
1211 
1212  // perform shift, if so desired
1213  if (shift != Point<patch_spacedim>())
1214  for (unsigned int i = old_n_patches; i < patches.size(); ++i)
1215  for (const unsigned int v : GeometryInfo<patch_dim>::vertex_indices())
1216  patches[i].vertices[v] += shift;
1217 
1218  // adjust patch numbers
1219  for (unsigned int i = old_n_patches; i < patches.size(); ++i)
1220  patches[i].patch_index += old_n_patches;
1221 
1222  // adjust patch neighbors
1223  for (unsigned int i = old_n_patches; i < patches.size(); ++i)
1224  for (const unsigned int n : GeometryInfo<patch_dim>::face_indices())
1225  if (patches[i].neighbors[n] != Patch::no_neighbor)
1226  patches[i].neighbors[n] += old_n_patches;
1227 }
1228 
1229 
1231 
1232 #endif
void add_data_vector(const VectorType &data, const std::string &name, const DataVectorType type=type_automatic, const std::vector< DataComponentInterpretation::DataComponentInterpretation > &data_component_interpretation={})
std::vector< std::shared_ptr<::hp::FECollection< dim, spacedim > > > get_fes() const
void clear_data_vectors()
void add_data_vector(const DoFHandler< dim, spacedim > &dof_handler, const VectorType &data, const std::string &name, const std::vector< DataComponentInterpretation::DataComponentInterpretation > &data_component_interpretation={})
std::vector< std::shared_ptr< internal::DataOutImplementation::DataEntryBase< dim, spacedim > > > cell_data
std::vector< Patch > patches
void add_data_vector(const DoFHandler< dim, spacedim > &dof_handler, const VectorType &data, const DataPostprocessor< spacedim > &data_postprocessor)
void merge_patches(const DataOut_DoFData< dim2, patch_dim, spacedim2, patch_spacedim > &source, const Point< patch_spacedim > &shift=Point< patch_spacedim >())
void add_data_vector_internal(const DoFHandler< dim, spacedim > *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)
virtual std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > get_nonscalar_data_ranges() const override
void add_data_vector(const VectorType &data, const DataPostprocessor< spacedim > &data_postprocessor)
std::vector< std::shared_ptr< internal::DataOutImplementation::DataEntryBase< dim, spacedim > > > dof_data
SmartPointer< const DoFHandler< dim, spacedim > > dofs
virtual void clear()
virtual std::vector< std::string > get_dataset_names() const override
void attach_dof_handler(const DoFHandler< dim, spacedim > &)
void add_data_vector(const DoFHandler< dim, spacedim > &dof_handler, const VectorType &data, const std::vector< std::string > &names, const std::vector< DataComponentInterpretation::DataComponentInterpretation > &data_component_interpretation={})
std::size_t memory_consumption() const
virtual const std::vector< Patch > & get_patches() const override
void add_mg_data_vector(const DoFHandler< dim, spacedim > &dof_handler, const MGLevelObject< VectorType > &data, const std::vector< std::string > &names, const std::vector< DataComponentInterpretation::DataComponentInterpretation > &data_component_interpretation=std::vector< DataComponentInterpretation::DataComponentInterpretation >())
void clear_input_data_references()
typename Triangulation< dim, spacedim >::cell_iterator cell_iterator
void add_mg_data_vector(const DoFHandler< dim, spacedim > &dof_handler, const MGLevelObject< VectorType > &data, const std::string &name)
void attach_triangulation(const Triangulation< dim, spacedim > &)
SmartPointer< const Triangulation< dim, spacedim > > triangulation
virtual ~DataOut_DoFData() override
void add_data_vector(const VectorType &data, const std::vector< std::string > &names, const DataVectorType type=type_automatic, const std::vector< DataComponentInterpretation::DataComponentInterpretation > &data_component_interpretation={})
friend class MGDataOut
Definition: point.h:112
Definition: vector.h:109
const std::vector< DataComponentInterpretation::DataComponentInterpretation > data_component_interpretation
virtual void get_function_values(const FEValuesBase< dim, spacedim > &fe_patch_values, const ComponentExtractor extract_component, std::vector< double > &patch_values) const =0
virtual void get_function_hessians(const FEValuesBase< dim, spacedim > &fe_patch_values, const ComponentExtractor extract_component, std::vector< std::vector< Tensor< 2, spacedim >>> &patch_hessians_system) const =0
virtual void get_function_hessians(const FEValuesBase< dim, spacedim > &fe_patch_values, const ComponentExtractor extract_component, std::vector< Tensor< 2, spacedim >> &patch_hessians) const =0
SmartPointer< const DoFHandler< dim, spacedim > > dof_handler
virtual void get_function_values(const FEValuesBase< dim, spacedim > &fe_patch_values, const ComponentExtractor extract_component, std::vector<::Vector< double >> &patch_values_system) const =0
virtual double get_cell_data_value(const unsigned int cell_number, const ComponentExtractor extract_component) const =0
virtual std::size_t memory_consumption() const =0
SmartPointer< const ::DataPostprocessor< spacedim > > postprocessor
virtual void get_function_gradients(const FEValuesBase< dim, spacedim > &fe_patch_values, const ComponentExtractor extract_component, std::vector< Tensor< 1, spacedim >> &patch_gradients) const =0
DataEntryBase(const DoFHandler< dim, spacedim > *dofs, const std::vector< std::string > &names, const std::vector< DataComponentInterpretation::DataComponentInterpretation > &data_component_interpretation)
DataEntryBase(const DoFHandler< dim, spacedim > *dofs, const DataPostprocessor< spacedim > *data_postprocessor)
virtual void get_function_gradients(const FEValuesBase< dim, spacedim > &fe_patch_values, const ComponentExtractor extract_component, std::vector< std::vector< Tensor< 1, spacedim >>> &patch_gradients_system) const =0
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:474
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:475
static ::ExceptionBase & ExcIncompatiblePatchLists()
static ::ExceptionBase & ExcInvalidVectorSize(int arg1, int arg2, int arg3)
static ::ExceptionBase & ExcOldDataStillPresent()
static ::ExceptionBase & ExcNoDoFHandlerSelected()
static ::ExceptionBase & ExcInvalidTensorDeclaration(int arg1, std::string arg2)
static ::ExceptionBase & ExcInvalidNumberOfSubdivisions(int arg1)
static ::ExceptionBase & ExcInvalidCharacter(std::string arg1, size_t arg2)
static ::ExceptionBase & ExcNoTriangulationSelected()
#define Assert(cond, exc)
Definition: exceptions.h:1614
static ::ExceptionBase & ExcIncompatibleDatasetNames()
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:533
static ::ExceptionBase & ExcInvalidNumberOfNames(int arg1, int arg2)
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:488
static ::ExceptionBase & ExcInvalidVectorDeclaration(int arg1, std::string arg2)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition: exceptions.h:556
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:510
static ::ExceptionBase & ExcMessage(std::string arg1)
UpdateFlags
void shift(const Tensor< 1, spacedim > &shift_vector, Triangulation< dim, spacedim > &triangulation)
Definition: grid_tools.cc:2131
static const unsigned int invalid_unsigned_int
Definition: types.h:213
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
std::vector< std::shared_ptr<::hp::FEFaceValues< dim, spacedim > > > x_fe_face_values
DataPostprocessorInputs::Scalar< spacedim > patch_values_scalar
const std::vector< std::shared_ptr<::hp::FECollection< dim, spacedim > > > finite_elements
std::vector< std::vector<::Vector< double > > > postprocessed_values
const ::hp::MappingCollection< dim, spacedim > mapping_collection
std::vector< std::shared_ptr<::hp::FEValues< dim, spacedim > > > x_fe_values
DataPostprocessorInputs::Vector< spacedim > patch_values_system
ParallelDataBase(const ParallelDataBase &data)
void reinit_all_fe_values(std::vector< std::shared_ptr< DataEntryBase< dim, spacedim >>> &dof_data, const typename ::Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face=numbers::invalid_unsigned_int)
const FEValuesBase< dim, spacedim > & get_present_fe_values(const unsigned int dataset) const
ParallelDataBase(const unsigned int n_datasets, const unsigned int n_subdivisions, const std::vector< unsigned int > &n_postprocessor_outputs, const Mapping< dim, spacedim > &mapping, const std::vector< std::shared_ptr<::hp::FECollection< dim, spacedim >>> &finite_elements, const UpdateFlags update_flags, const bool use_face_values)
ParallelDataBase(const unsigned int n_datasets, const unsigned int n_subdivisions, const std::vector< unsigned int > &n_postprocessor_outputs, const ::hp::MappingCollection< dim, spacedim > &mapping, const std::vector< std::shared_ptr<::hp::FECollection< dim, spacedim >>> &finite_elements, const UpdateFlags update_flags, const bool use_face_values)
void resize_system_vectors(const unsigned int n_components)