Reference documentation for deal.II version Git bcb2fdb641 2020-09-17 14:58:24 -0400
\(\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 - 2020 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.");
70 
75  "For the operation you are attempting, you first need to "
76  "tell the DataOut or related object which DoFHandler "
77  "you would like to work on.");
78 
83  int,
84  int,
85  int,
86  << "The vector has size " << arg1
87  << " but the DoFHandler object says that there are " << arg2
88  << " degrees of freedom and there are " << arg3
89  << " active cells. The size of your vector needs to be"
90  << " either equal to the number of degrees of freedom (when"
91  << " the data is of type type_dof_data), or equal to the"
92  << " number of active cells (when the data is of type "
93  << " type_cell_data).");
99  std::string,
100  size_t,
101  << "Please use only the characters [a-zA-Z0-9_<>()] for" << std::endl
102  << "description strings since some graphics formats will only accept these."
103  << std::endl
104  << "The string you gave was <" << arg1
105  << ">, within which the invalid character is <" << arg1[arg2] << ">."
106  << std::endl);
112  "When attaching a triangulation or DoFHandler object, it is "
113  "not allowed if old data vectors are still referenced. If "
114  "you want to reuse an object of the current type, you first "
115  "need to call the 'clear_data_vector()' function.");
120  int,
121  int,
122  << "You have to give one name per component in your "
123  << "data vector. The number you gave was " << arg1
124  << ", but the number of components is " << arg2 << ".");
129  "While merging sets of patches, the two sets to be merged "
130  "need to refer to data that agrees on the names of the "
131  "various variables represented. In other words, you "
132  "cannot merge sets of patches that originate from "
133  "entirely unrelated simulations.");
138  "While merging sets of patches, the two sets to be merged "
139  "need to refer to data that agrees on the number of "
140  "subdivisions and other properties. In other words, you "
141  "cannot merge sets of patches that originate from "
142  "entirely unrelated simulations.");
143 
145  int,
146  std::string,
147  << "When declaring that a number of components in a data "
148  << "set to be output logically form a vector instead of "
149  << "simply a set of scalar fields, you need to specify "
150  << "this for all relevant components. Furthermore, "
151  << "vectors must always consist of exactly <dim> "
152  << "components. However, the vector component at "
153  << "position " << arg1 << " with name <" << arg2
154  << "> does not satisfy these conditions.");
155 
157  int,
158  std::string,
159  << "When declaring that a number of components in a data "
160  << "set to be output logically form a tensor instead of "
161  << "simply a set of scalar fields, you need to specify "
162  << "this for all relevant components. Furthermore, "
163  << "tensors must always consist of exactly <dim*dim> "
164  << "components. However, the tensor component at "
165  << "position " << arg1 << " with name <" << arg2
166  << "> does not satisfy these conditions.");
167 
168  } // namespace DataOutImplementation
169 } // namespace Exceptions
170 
171 
172 namespace internal
173 {
174  namespace DataOutImplementation
175  {
200  {
201  real_part,
203  };
204 
205 
221  template <typename DoFHandlerType>
223  {
224  public:
230  DataEntryBase(const DoFHandlerType * dofs,
231  const std::vector<std::string> &names,
232  const std::vector<
234  &data_component_interpretation);
235 
241  DataEntryBase(const DoFHandlerType *dofs,
243  *data_postprocessor);
244 
248  virtual ~DataEntryBase() = default;
249 
254  virtual double
255  get_cell_data_value(const unsigned int cell_number,
256  const ComponentExtractor extract_component) const = 0;
257 
262  virtual void
263  get_function_values(
264  const FEValuesBase<DoFHandlerType::dimension,
265  DoFHandlerType::space_dimension> &fe_patch_values,
266  const ComponentExtractor extract_component,
267  std::vector<double> &patch_values) const = 0;
268 
274  virtual void
275  get_function_values(
276  const FEValuesBase<DoFHandlerType::dimension,
277  DoFHandlerType::space_dimension> &fe_patch_values,
278  const ComponentExtractor extract_component,
279  std::vector<::Vector<double>> &patch_values_system) const = 0;
280 
285  virtual void
286  get_function_gradients(
287  const FEValuesBase<DoFHandlerType::dimension,
288  DoFHandlerType::space_dimension> &fe_patch_values,
289  const ComponentExtractor extract_component,
291  &patch_gradients) const = 0;
292 
298  virtual void
299  get_function_gradients(
300  const FEValuesBase<DoFHandlerType::dimension,
301  DoFHandlerType::space_dimension> &fe_patch_values,
302  const ComponentExtractor extract_component,
303  std::vector<std::vector<Tensor<1, DoFHandlerType::space_dimension>>>
304  &patch_gradients_system) const = 0;
305 
310  virtual void
311  get_function_hessians(
312  const FEValuesBase<DoFHandlerType::dimension,
313  DoFHandlerType::space_dimension> &fe_patch_values,
314  const ComponentExtractor extract_component,
315  std::vector<Tensor<2, DoFHandlerType::space_dimension>> &patch_hessians)
316  const = 0;
317 
323  virtual void
324  get_function_hessians(
325  const FEValuesBase<DoFHandlerType::dimension,
326  DoFHandlerType::space_dimension> &fe_patch_values,
327  const ComponentExtractor extract_component,
328  std::vector<std::vector<Tensor<2, DoFHandlerType::space_dimension>>>
329  &patch_hessians_system) const = 0;
330 
335  virtual bool
336  is_complex_valued() const = 0;
337 
341  virtual void
342  clear() = 0;
343 
348  virtual std::size_t
349  memory_consumption() const = 0;
350 
355 
359  const std::vector<std::string> names;
360 
366  const std::vector<
369 
374  SmartPointer<
375  const ::DataPostprocessor<DoFHandlerType::space_dimension>>
377 
384  unsigned int n_output_variables;
385  };
386 
387 
404  template <int dim, int spacedim>
406  {
408  const unsigned int n_datasets,
409  const unsigned int n_subdivisions,
410  const std::vector<unsigned int> &n_postprocessor_outputs,
411  const Mapping<dim, spacedim> & mapping,
412  const std::vector<
413  std::shared_ptr<::hp::FECollection<dim, spacedim>>>
414  & finite_elements,
415  const UpdateFlags update_flags,
416  const bool use_face_values);
417 
419  const unsigned int n_datasets,
420  const unsigned int n_subdivisions,
421  const std::vector<unsigned int> &n_postprocessor_outputs,
422  const ::hp::MappingCollection<dim, spacedim> &mapping,
423  const std::vector<
424  std::shared_ptr<::hp::FECollection<dim, spacedim>>>
425  & finite_elements,
426  const UpdateFlags update_flags,
427  const bool use_face_values);
428 
429  ParallelDataBase(const ParallelDataBase &data);
430 
431  template <typename DoFHandlerType>
432  void
433  reinit_all_fe_values(
434  std::vector<std::shared_ptr<DataEntryBase<DoFHandlerType>>> &dof_data,
435  const typename ::Triangulation<dim, spacedim>::cell_iterator
436  & cell,
437  const unsigned int face = numbers::invalid_unsigned_int);
438 
440  get_present_fe_values(const unsigned int dataset) const;
441 
442  void
443  resize_system_vectors(const unsigned int n_components);
444 
445  const unsigned int n_datasets;
446  const unsigned int n_subdivisions;
447 
450  std::vector<std::vector<::Vector<double>>> postprocessed_values;
451 
452  const ::hp::MappingCollection<dim, spacedim> mapping_collection;
453  const std::vector<
454  std::shared_ptr<::hp::FECollection<dim, spacedim>>>
457 
458  std::vector<std::shared_ptr<::hp::FEValues<dim, spacedim>>>
460  std::vector<std::shared_ptr<::hp::FEFaceValues<dim, spacedim>>>
462  };
463  } // namespace DataOutImplementation
464 } // namespace internal
465 
466 
467 // TODO: Most of the documentation of DataOut_DoFData applies to DataOut.
468 
609 template <typename DoFHandlerType,
610  int patch_dim,
611  int patch_space_dim = patch_dim>
612 class DataOut_DoFData : public DataOutInterface<patch_dim, patch_space_dim>
613 {
614 public:
619  using cell_iterator =
620  typename Triangulation<DoFHandlerType::dimension,
621  DoFHandlerType::space_dimension>::cell_iterator;
622 
623 public:
633  {
638 
643 
647  type_automatic
648  };
649 
653  DataOut_DoFData();
654 
658  virtual ~DataOut_DoFData() override;
659 
668  void
669  attach_dof_handler(const DoFHandlerType &);
670 
680  void
681  attach_triangulation(const Triangulation<DoFHandlerType::dimension,
682  DoFHandlerType::space_dimension> &);
683 
752  template <class VectorType>
753  void
754  add_data_vector(
755  const VectorType & data,
756  const std::vector<std::string> &names,
757  const DataVectorType type = type_automatic,
758  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
759  &data_component_interpretation = std::vector<
761 
778  template <class VectorType>
779  void
780  add_data_vector(
781  const VectorType & data,
782  const std::string & name,
783  const DataVectorType type = type_automatic,
784  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
785  &data_component_interpretation = std::vector<
787 
803  template <class VectorType>
804  void
805  add_data_vector(
806  const DoFHandlerType & dof_handler,
807  const VectorType & data,
808  const std::vector<std::string> &names,
809  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
810  &data_component_interpretation = std::vector<
812 
813 
818  template <class VectorType>
819  void
820  add_data_vector(
821  const DoFHandlerType &dof_handler,
822  const VectorType & data,
823  const std::string & name,
824  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
825  &data_component_interpretation = std::vector<
827 
854  template <class VectorType>
855  void
856  add_data_vector(const VectorType &data,
858  &data_postprocessor);
859 
866  template <class VectorType>
867  void
868  add_data_vector(const DoFHandlerType &dof_handler,
869  const VectorType & data,
871  &data_postprocessor);
872 
890  template <class VectorType>
891  void
892  add_mg_data_vector(
893  const DoFHandlerType & dof_handler,
894  const MGLevelObject<VectorType> &data,
895  const std::vector<std::string> & names,
896  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
897  &data_component_interpretation = std::vector<
899 
903  template <class VectorType>
904  void
905  add_mg_data_vector(const DoFHandlerType & dof_handler,
906  const MGLevelObject<VectorType> &data,
907  const std::string & name);
908 
915  void
916  clear_data_vectors();
917 
928  void
929  clear_input_data_references();
930 
954  template <typename DoFHandlerType2>
955  void
956  merge_patches(
959 
966  virtual void
967  clear();
968 
973  std::size_t
974  memory_consumption() const;
975 
976 protected:
981 
985  SmartPointer<const Triangulation<DoFHandlerType::dimension,
986  DoFHandlerType::space_dimension>>
988 
993 
997  std::vector<std::shared_ptr<
1000 
1004  std::vector<std::shared_ptr<
1007 
1013  std::vector<Patch> patches;
1014 
1019  virtual const std::vector<Patch> &
1020  get_patches() const override;
1021 
1026  virtual std::vector<std::string>
1027  get_dataset_names() const override;
1028 
1033  std::vector<
1034  std::shared_ptr<::hp::FECollection<DoFHandlerType::dimension,
1035  DoFHandlerType::space_dimension>>>
1036  get_fes() const;
1037 
1042  virtual std::vector<
1043  std::tuple<unsigned int,
1044  unsigned int,
1045  std::string,
1047  get_nonscalar_data_ranges() const override;
1048 
1049  // Make all template siblings friends. Needed for the merge_patches()
1050  // function.
1051  template <class, int, int>
1052  friend class DataOut_DoFData;
1053 
1056  template <int, class>
1057  friend class MGDataOut;
1058 
1059 private:
1063  template <class VectorType>
1064  void
1065  add_data_vector_internal(
1066  const DoFHandlerType * dof_handler,
1067  const VectorType & data,
1068  const std::vector<std::string> &names,
1069  const DataVectorType type,
1070  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
1071  & data_component_interpretation,
1072  const bool deduce_output_names);
1073 };
1074 
1075 
1076 
1077 // -------------------- template and inline functions ------------------------
1078 template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
1079 template <typename VectorType>
1080 void
1082  const VectorType & vec,
1083  const std::string & name,
1084  const DataVectorType type,
1085  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
1086  &data_component_interpretation)
1087 {
1088  Assert(triangulation != nullptr,
1090  std::vector<std::string> names(1, name);
1091  add_data_vector_internal(
1092  dofs, vec, names, type, data_component_interpretation, true);
1093 }
1094 
1095 
1096 
1097 template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
1098 template <typename VectorType>
1099 void
1101  const VectorType & vec,
1102  const std::vector<std::string> &names,
1103  const DataVectorType type,
1104  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
1105  &data_component_interpretation)
1106 {
1107  Assert(triangulation != nullptr,
1109  add_data_vector_internal(
1110  dofs, vec, names, type, data_component_interpretation, false);
1111 }
1112 
1113 
1114 
1115 template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
1116 template <typename VectorType>
1117 void
1119  const DoFHandlerType &dof_handler,
1120  const VectorType & data,
1121  const std::string & name,
1122  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
1123  &data_component_interpretation)
1124 {
1125  std::vector<std::string> names(1, name);
1126  add_data_vector_internal(&dof_handler,
1127  data,
1128  names,
1129  type_dof_data,
1130  data_component_interpretation,
1131  true);
1132 }
1133 
1134 
1135 
1136 template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
1137 template <typename VectorType>
1138 void
1140  const DoFHandlerType & dof_handler,
1141  const VectorType & data,
1142  const std::vector<std::string> &names,
1143  const std::vector<DataComponentInterpretation::DataComponentInterpretation>
1144  &data_component_interpretation)
1145 {
1146  add_data_vector_internal(&dof_handler,
1147  data,
1148  names,
1149  type_dof_data,
1150  data_component_interpretation,
1151  false);
1152 }
1153 
1154 
1155 
1156 template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
1157 template <typename VectorType>
1158 void
1160  const VectorType & vec,
1161  const DataPostprocessor<DoFHandlerType::space_dimension> &data_postprocessor)
1162 {
1163  Assert(dofs != nullptr,
1165  add_data_vector(*dofs, vec, data_postprocessor);
1166 }
1167 
1168 
1169 
1170 template <typename DoFHandlerType, int patch_dim, int patch_space_dim>
1171 template <typename DoFHandlerType2>
1172 void
1175  const Point<patch_space_dim> & shift)
1176 {
1177  const std::vector<Patch> &source_patches = source.get_patches();
1178  Assert((patches.size() != 0) && (source_patches.size() != 0),
1179  ExcMessage("When calling this function, both the current "
1180  "object and the one being merged need to have a "
1181  "nonzero number of patches associated with it. "
1182  "Either you called this function on objects that "
1183  "are empty, or you may have forgotten to call "
1184  "the 'build_patches()' function."));
1185  // Check equality of component names
1186  Assert(get_dataset_names() == source.get_dataset_names(),
1188 
1189  // Make sure patches are compatible. Ideally, we would check that all input
1190  // patches from both collections are all compatible, but we'll be content
1191  // with checking that just the first ones from both sources are.
1192  //
1193  // We check compatibility by testing that both sets of patches result
1194  // from the same number of subdivisions, and that they have the same
1195  // number of source vectors (they really should, since we already checked
1196  // that there are the same number of source components above, but you
1197  // never know). This implies that the data should have the same number of
1198  // columns. They should really have the same number of rows as well,
1199  // but depending on whether a patch has points included or not, the
1200  // number of rows may or may not include coordinates for the points,
1201  // and the comparison has to account for that because in each source
1202  // stream, the patches may include some that have points included.
1203  Assert(patches[0].n_subdivisions == source_patches[0].n_subdivisions,
1205  Assert(patches[0].data.n_cols() == source_patches[0].data.n_cols(),
1207  Assert((patches[0].data.n_rows() +
1208  (patches[0].points_are_available ? 0 : patch_space_dim)) ==
1209  (source_patches[0].data.n_rows() +
1210  (source_patches[0].points_are_available ? 0 : patch_space_dim)),
1212 
1213  // check equality of the vector data
1214  // specifications
1215  Assert(get_nonscalar_data_ranges().size() ==
1216  source.get_nonscalar_data_ranges().size(),
1217  ExcMessage("Both sources need to declare the same components "
1218  "as vectors."));
1219  for (unsigned int i = 0; i < get_nonscalar_data_ranges().size(); ++i)
1220  {
1221  Assert(std::get<0>(get_nonscalar_data_ranges()[i]) ==
1222  std::get<0>(source.get_nonscalar_data_ranges()[i]),
1223  ExcMessage("Both sources need to declare the same components "
1224  "as vectors."));
1225  Assert(std::get<1>(get_nonscalar_data_ranges()[i]) ==
1226  std::get<1>(source.get_nonscalar_data_ranges()[i]),
1227  ExcMessage("Both sources need to declare the same components "
1228  "as vectors."));
1229  Assert(std::get<2>(get_nonscalar_data_ranges()[i]) ==
1230  std::get<2>(source.get_nonscalar_data_ranges()[i]),
1231  ExcMessage("Both sources need to declare the same components "
1232  "as vectors."));
1233  }
1234 
1235  // merge patches. store old number
1236  // of elements, since we need to
1237  // adjust patch numbers, etc
1238  // afterwards
1239  const unsigned int old_n_patches = patches.size();
1240  patches.insert(patches.end(), source_patches.begin(), source_patches.end());
1241 
1242  // perform shift, if so desired
1243  if (shift != Point<patch_space_dim>())
1244  for (unsigned int i = old_n_patches; i < patches.size(); ++i)
1245  for (const unsigned int v : GeometryInfo<patch_dim>::vertex_indices())
1246  patches[i].vertices[v] += shift;
1247 
1248  // adjust patch numbers
1249  for (unsigned int i = old_n_patches; i < patches.size(); ++i)
1250  patches[i].patch_index += old_n_patches;
1251 
1252  // adjust patch neighbors
1253  for (unsigned int i = old_n_patches; i < patches.size(); ++i)
1254  for (const unsigned int n : GeometryInfo<patch_dim>::face_indices())
1255  if (patches[i].neighbors[n] != Patch::no_neighbor)
1256  patches[i].neighbors[n] += old_n_patches;
1257 }
1258 
1259 
1261 
1262 #endif
void merge_patches(const DataOut_DoFData< DoFHandlerType2, patch_dim, patch_space_dim > &source, const Point< patch_space_dim > &shift=Point< patch_space_dim >())
static const unsigned int invalid_unsigned_int
Definition: types.h:196
static ::ExceptionBase & ExcNoTriangulationSelected()
static ::ExceptionBase & ExcIncompatiblePatchLists()
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:538
const std::vector< std::shared_ptr<::hp::FECollection< dim, spacedim > > > finite_elements
virtual std::vector< std::string > get_dataset_names() const override
static ::ExceptionBase & ExcOldDataStillPresent()
static ::ExceptionBase & ExcInvalidVectorSize(int arg1, int arg2, int arg3)
static ::ExceptionBase & ExcInvalidNumberOfSubdivisions(int arg1)
Definition: point.h:110
typename Triangulation< DoFHandlerType::dimension, DoFHandlerType::space_dimension >::cell_iterator cell_iterator
static ::ExceptionBase & ExcIncompatibleDatasetNames()
std::vector< std::shared_ptr< internal::DataOutImplementation::DataEntryBase< DoFHandlerType > > > dof_data
static ::ExceptionBase & ExcMessage(std::string arg1)
static ::ExceptionBase & ExcNoDoFHandlerSelected()
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:515
#define Assert(cond, exc)
Definition: exceptions.h:1411
UpdateFlags
Abstract base class for mapping classes.
Definition: mapping.h:301
static ::ExceptionBase & ExcInvalidNumberOfNames(int arg1, int arg2)
static ::ExceptionBase & ExcInvalidVectorDeclaration(int arg1, std::string arg2)
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
static ::ExceptionBase & ExcInvalidCharacter(std::string arg1, size_t arg2)
std::vector< Patch > patches
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=std::vector< DataComponentInterpretation::DataComponentInterpretation >())
Point< 3 > vertices[4]
SmartPointer< const ::DataPostprocessor< DoFHandlerType::space_dimension > > postprocessor
DataPostprocessorInputs::Scalar< spacedim > patch_values_scalar
virtual const std::vector< Patch > & get_patches() const override
const ::hp::MappingCollection< dim, spacedim > mapping_collection
static ::ExceptionBase & ExcInvalidTensorDeclaration(int arg1, std::string arg2)
DataPostprocessorInputs::Vector< spacedim > patch_values_system
virtual std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > get_nonscalar_data_ranges() const override
SmartPointer< const DoFHandlerType > dof_handler
Definition: tensor.h:448
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition: exceptions.h:561
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
SmartPointer< const Triangulation< DoFHandlerType::dimension, DoFHandlerType::space_dimension > > triangulation
std::vector< std::shared_ptr<::hp::FEValues< dim, spacedim > > > x_fe_values
std::vector< std::shared_ptr< internal::DataOutImplementation::DataEntryBase< DoFHandlerType > > > cell_data
const std::vector< DataComponentInterpretation::DataComponentInterpretation > data_component_interpretation
SmartPointer< const DoFHandlerType > dofs
std::vector< std::vector<::Vector< double > > > postprocessed_values
std::vector< std::shared_ptr<::hp::FEFaceValues< dim, spacedim > > > x_fe_face_values
void shift(const Tensor< 1, spacedim > &shift_vector, Triangulation< dim, spacedim > &triangulation)
Definition: grid_tools.cc:815
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)