Reference documentation for deal.II version Git 78a8940608 2021-04-17 09:24:19 -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\}}\)
scratch_data.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2019 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_meshworker_scratch_data_h
17 #define dealii_meshworker_scratch_data_h
18 
19 #include <deal.II/base/config.h>
20 
22 
24 
26 #include <deal.II/fe/fe_values.h>
27 #include <deal.II/fe/mapping_q1.h>
28 
29 #include <boost/any.hpp>
30 
31 #include <algorithm>
32 #include <map>
33 #include <vector>
34 
36 
37 namespace MeshWorker
38 {
39  namespace internal
40  {
45  template <int dim, int spacedim, typename NumberType, typename Extractor>
48  } // namespace internal
49 
230  template <int dim, int spacedim = dim>
232  {
233  public:
253  ScratchData(
254  const Mapping<dim, spacedim> & mapping,
256  const Quadrature<dim> & quadrature,
257  const UpdateFlags & update_flags,
258  const Quadrature<dim - 1> &face_quadrature = Quadrature<dim - 1>(),
259  const UpdateFlags & face_update_flags = update_default);
260 
277  ScratchData(
278  const Mapping<dim, spacedim> & mapping,
280  const Quadrature<dim> & quadrature,
281  const UpdateFlags & update_flags,
282  const UpdateFlags & neighbor_update_flags,
283  const Quadrature<dim - 1> &face_quadrature = Quadrature<dim - 1>(),
284  const UpdateFlags & face_update_flags = update_default,
285  const UpdateFlags & neighbor_face_update_flags = update_default);
286 
299  ScratchData(
301  const Quadrature<dim> & quadrature,
302  const UpdateFlags & update_flags,
303  const Quadrature<dim - 1> &face_quadrature = Quadrature<dim - 1>(),
304  const UpdateFlags & face_update_flags = update_default);
305 
320  ScratchData(
322  const Quadrature<dim> & quadrature,
323  const UpdateFlags & update_flags,
324  const UpdateFlags & neighbor_update_flags,
325  const Quadrature<dim - 1> &face_quadrature = Quadrature<dim - 1>(),
326  const UpdateFlags & face_update_flags = update_default,
327  const UpdateFlags & neighbor_face_update_flags = update_default);
328 
332  ScratchData(const ScratchData<dim, spacedim> &scratch);
333  // CurrentCellMethods
338 
347  reinit(
349 
359  const unsigned int face_no);
360 
373  const unsigned int face_no,
374  const unsigned int subface_no);
375 
392  const unsigned int face_no,
393  const unsigned int sub_face_no,
395  & cell_neighbor,
396  const unsigned int face_no_neighbor,
397  const unsigned int sub_face_no_neighbor);
398 
410  get_current_fe_values() const;
411 
415  const std::vector<Point<spacedim>> &
416  get_quadrature_points() const;
417 
421  const std::vector<double> &
422  get_JxW_values() const;
423 
427  const std::vector<Tensor<1, spacedim>> &
428  get_normal_vectors() const;
429 
434  const std::vector<types::global_dof_index> &
435  get_local_dof_indices() const;
436  // CurrentCellMethods
438  // NeighborCellMethods
443 
452  reinit_neighbor(
454 
463  reinit_neighbor(
465  const unsigned int face_no);
466 
478  reinit_neighbor(
480  const unsigned int face_no,
481  const unsigned int subface_no);
482 
495  get_current_neighbor_fe_values() const;
496 
500  const std::vector<double> &
501  get_neighbor_JxW_values() const;
502 
506  const std::vector<Tensor<1, spacedim>> &
507  get_neighbor_normal_vectors();
508 
513  const std::vector<types::global_dof_index> &
514  get_neighbor_dof_indices() const;
515  // NeighborCellMethods
517 
524  get_general_data_storage();
525 
531  const GeneralDataStorage &
532  get_general_data_storage() const;
533  // CurrentCellEvaluation
538 
568  template <typename VectorType, typename Number = double>
569  void
570  extract_local_dof_values(const std::string &global_vector_name,
571  const VectorType & input_vector,
572  const Number dummy = Number(0));
573 
582  template <typename Number = double>
583  const std::vector<Number> &
584  get_local_dof_values(const std::string &global_vector_name,
585  Number dummy = Number(0)) const;
586 
608  template <typename Extractor, typename Number = double>
609  const std::vector<
611  value_type> &
612  get_values(const std::string &global_vector_name,
613  const Extractor & variable,
614  const Number dummy = Number(0));
615 
637  template <typename Extractor, typename Number = double>
638  const std::vector<
640  gradient_type> &
641  get_gradients(const std::string &global_vector_name,
642  const Extractor & variable,
643  const Number dummy = Number(0));
644 
667  template <typename Extractor, typename Number = double>
668  const std::vector<
670  symmetric_gradient_type> &
671  get_symmetric_gradients(const std::string &global_vector_name,
672  const Extractor & variable,
673  const Number dummy = Number(0));
674 
696  template <typename Extractor, typename Number = double>
697  const std::vector<
699  divergence_type> &
700  get_divergences(const std::string &global_vector_name,
701  const Extractor & variable,
702  const Number dummy = Number(0));
703 
725  template <typename Extractor, typename Number = double>
726  const std::vector<typename internal::
727  OutputType<dim, spacedim, Number, Extractor>::curl_type>
728  &
729  get_curls(const std::string &global_vector_name,
730  const Extractor & variable,
731  const Number dummy = Number(0));
732 
754  template <typename Extractor, typename Number = double>
755  const std::vector<
757  hessian_type> &
758  get_hessians(const std::string &global_vector_name,
759  const Extractor & variable,
760  const Number dummy = Number(0));
761 
783  template <typename Extractor, typename Number = double>
784  const std::vector<
786  laplacian_type> &
787  get_laplacians(const std::string &global_vector_name,
788  const Extractor & variable,
789  const Number dummy = Number(0));
790 
812  template <typename Extractor, typename Number = double>
813  const std::vector<
815  third_derivative_type> &
816  get_third_derivatives(const std::string &global_vector_name,
817  const Extractor & variable,
818  const Number dummy = Number(0));
819  // CurrentCellEvaluation
821 
825  const Mapping<dim, spacedim> &
826  get_mapping() const;
827 
828  private:
833  template <typename Extractor, typename Number = double>
834  std::string
835  get_unique_name(const std::string &global_vector_name,
836  const Extractor & variable,
837  const std::string &object_type,
838  const unsigned int size,
839  const Number & exemplar_number) const;
840 
844  template <typename Number = double>
845  std::string
846  get_unique_dofs_name(const std::string &global_vector_name,
847  const unsigned int size,
848  const Number & exemplar_number) const;
849 
855 
861 
867 
873 
878 
883 
889 
895 
899  std::unique_ptr<FEValues<dim, spacedim>> fe_values;
900 
904  std::unique_ptr<FEFaceValues<dim, spacedim>> fe_face_values;
905 
909  std::unique_ptr<FESubfaceValues<dim, spacedim>> fe_subface_values;
910 
914  std::unique_ptr<FEValues<dim, spacedim>> neighbor_fe_values;
915 
919  std::unique_ptr<FEFaceValues<dim, spacedim>> neighbor_fe_face_values;
920 
924  std::unique_ptr<FESubfaceValues<dim, spacedim>> neighbor_fe_subface_values;
925 
929  std::unique_ptr<FEInterfaceValues<dim, spacedim>> interface_fe_values;
930 
934  std::vector<types::global_dof_index> local_dof_indices;
935 
939  std::vector<types::global_dof_index> neighbor_dof_indices;
940 
945 
950 
956 
962  };
963 
964 #ifndef DOXYGEN
965  template <int dim, int spacedim>
966  template <typename Extractor, typename Number>
967  std::string
969  const std::string &global_vector_name,
970  const Extractor & variable,
971  const std::string &object_type,
972  const unsigned int size,
973  const Number & exemplar_number) const
974  {
975  return global_vector_name + "_" + variable.get_name() + "_" + object_type +
976  "_" + Utilities::int_to_string(size) + "_" +
977  Utilities::type_to_string(exemplar_number);
978  }
979 
980 
981 
982  template <int dim, int spacedim>
983  template <typename Number>
984  std::string
986  const std::string &global_vector_name,
987  const unsigned int size,
988  const Number & exemplar_number) const
989  {
990  return global_vector_name + "_independent_local_dofs_" +
991  Utilities::int_to_string(size) + "_" +
992  Utilities::type_to_string(exemplar_number);
993  }
994 
995 
996 
997  template <int dim, int spacedim>
998  template <typename VectorType, typename Number>
999  void
1001  const std::string &global_vector_name,
1002  const VectorType & input_vector,
1003  const Number dummy)
1004  {
1005  const unsigned int n_dofs = local_dof_indices.size();
1006 
1007  const std::string name =
1008  get_unique_dofs_name(global_vector_name, n_dofs, dummy);
1009 
1010  auto &independent_local_dofs =
1011  internal_data_storage
1012  .template get_or_add_object_with_name<std::vector<Number>>(name,
1013  n_dofs);
1014 
1015  AssertDimension(independent_local_dofs.size(), n_dofs);
1016 
1018  for (unsigned int i = 0; i < n_dofs; ++i)
1020  input_vector(local_dof_indices[i]),
1021  i,
1022  n_dofs,
1023  independent_local_dofs[i]);
1024  else
1025  for (unsigned int i = 0; i < n_dofs; ++i)
1026  independent_local_dofs[i] = input_vector(local_dof_indices[i]);
1027  }
1028 
1029 
1030 
1031  template <int dim, int spacedim>
1032  template <typename Number>
1033  const std::vector<Number> &
1035  const std::string &global_vector_name,
1036  Number dummy) const
1037  {
1038  const unsigned int n_dofs =
1039  get_current_fe_values().get_fe().n_dofs_per_cell();
1040 
1041  const std::string dofs_name =
1042  get_unique_dofs_name(global_vector_name, n_dofs, dummy);
1043 
1044  Assert(
1045  internal_data_storage.stores_object_with_name(dofs_name),
1046  ExcMessage(
1047  "You did not call yet extract_local_dof_values with the right types!"));
1048 
1049  return internal_data_storage
1050  .template get_object_with_name<std::vector<Number>>(dofs_name);
1051  }
1052 
1053 
1054 
1055  template <int dim, int spacedim>
1056  template <typename Extractor, typename Number>
1057  const std::vector<
1059  &
1061  const std::string &global_vector_name,
1062  const Extractor & variable,
1063  const Number dummy)
1064  {
1065  const std::vector<Number> &independent_local_dofs =
1066  get_local_dof_values(global_vector_name, dummy);
1067 
1068  const FEValuesBase<dim, spacedim> &fev = get_current_fe_values();
1069 
1070  const unsigned int n_q_points = fev.n_quadrature_points;
1071 
1072  const std::string name = get_unique_name(
1073  global_vector_name, variable, "_values_q", n_q_points, dummy);
1074 
1075  // Now build the return type
1076  using RetType =
1077  std::vector<typename internal::
1078  OutputType<dim, spacedim, Number, Extractor>::value_type>;
1079 
1080  RetType &ret =
1081  internal_data_storage.template get_or_add_object_with_name<RetType>(
1082  name, n_q_points);
1083 
1084  AssertDimension(ret.size(), n_q_points);
1085 
1086  fev[variable].get_function_values_from_local_dof_values(
1087  independent_local_dofs, ret);
1088  return ret;
1089  }
1090 
1091 
1092 
1093  template <int dim, int spacedim>
1094  template <typename Extractor, typename Number>
1095  const std::vector<
1097  gradient_type> &
1099  const std::string &global_vector_name,
1100  const Extractor & variable,
1101  const Number dummy)
1102  {
1103  const std::vector<Number> &independent_local_dofs =
1104  get_local_dof_values(global_vector_name, dummy);
1105 
1106  const FEValuesBase<dim, spacedim> &fev = get_current_fe_values();
1107 
1108  const unsigned int n_q_points = fev.n_quadrature_points;
1109 
1110  const std::string name = get_unique_name(
1111  global_vector_name, variable, "_gradients_q", n_q_points, dummy);
1112 
1113  // Now build the return type
1114  using RetType = std::vector<
1116  gradient_type>;
1117 
1118  RetType &ret =
1119  internal_data_storage.template get_or_add_object_with_name<RetType>(
1120  name, n_q_points);
1121 
1122  AssertDimension(ret.size(), n_q_points);
1123 
1124  fev[variable].get_function_gradients_from_local_dof_values(
1125  independent_local_dofs, ret);
1126  return ret;
1127  }
1128 
1129 
1130 
1131  template <int dim, int spacedim>
1132  template <typename Extractor, typename Number>
1133  const std::vector<
1135  hessian_type> &
1137  const std::string &global_vector_name,
1138  const Extractor & variable,
1139  const Number dummy)
1140  {
1141  const std::vector<Number> &independent_local_dofs =
1142  get_local_dof_values(global_vector_name, dummy);
1143 
1144  const FEValuesBase<dim, spacedim> &fev = get_current_fe_values();
1145 
1146  const unsigned int n_q_points = fev.n_quadrature_points;
1147 
1148  const std::string name = get_unique_name(
1149  global_vector_name, variable, "_hessians_q", n_q_points, dummy);
1150 
1151  // Now build the return type
1152  using RetType =
1153  std::vector<typename internal::
1154  OutputType<dim, spacedim, Number, Extractor>::hessian_type>;
1155 
1156  RetType &ret =
1157  internal_data_storage.template get_or_add_object_with_name<RetType>(
1158  name, n_q_points);
1159 
1160 
1161  AssertDimension(ret.size(), n_q_points);
1162 
1163  fev[variable].get_function_hessians_from_local_dof_values(
1164  independent_local_dofs, ret);
1165  return ret;
1166  }
1167 
1168 
1169 
1170  template <int dim, int spacedim>
1171  template <typename Extractor, typename Number>
1172  const std::vector<
1174  laplacian_type> &
1176  const std::string &global_vector_name,
1177  const Extractor & variable,
1178  const Number dummy)
1179  {
1180  const std::vector<Number> &independent_local_dofs =
1181  get_local_dof_values(global_vector_name, dummy);
1182 
1183  const FEValuesBase<dim, spacedim> &fev = get_current_fe_values();
1184 
1185  const unsigned int n_q_points = fev.n_quadrature_points;
1186 
1187  const std::string name = get_unique_name(
1188  global_vector_name, variable, "_laplacians_q", n_q_points, dummy);
1189 
1190  // Now build the return type
1191  using RetType = std::vector<
1193  laplacian_type>;
1194 
1195  RetType &ret =
1196  internal_data_storage.template get_or_add_object_with_name<RetType>(
1197  name, n_q_points);
1198 
1199 
1200  AssertDimension(ret.size(), n_q_points);
1201 
1202  fev[variable].get_function_laplacians_from_local_dof_values(
1203  independent_local_dofs, ret);
1204  return ret;
1205  }
1206 
1207 
1208 
1209  template <int dim, int spacedim>
1210  template <typename Extractor, typename Number>
1211  const std::vector<
1213  third_derivative_type> &
1215  const std::string &global_vector_name,
1216  const Extractor & variable,
1217  const Number dummy)
1218  {
1219  const std::vector<Number> &independent_local_dofs =
1220  get_local_dof_values(global_vector_name, dummy);
1221 
1222  const FEValuesBase<dim, spacedim> &fev = get_current_fe_values();
1223 
1224  const unsigned int n_q_points = fev.n_quadrature_points;
1225 
1226  const std::string name = get_unique_name(
1227  global_vector_name, variable, "_third_derivatives_q", n_q_points, dummy);
1228 
1229  // Now build the return type
1230  using RetType = std::vector<
1232  third_derivative_type>;
1233 
1234  RetType &ret =
1235  internal_data_storage.template get_or_add_object_with_name<RetType>(
1236  name, n_q_points);
1237 
1238 
1239  AssertDimension(ret.size(), n_q_points);
1240 
1241  fev[variable].get_function_third_derivatives_from_local_dof_values(
1242  independent_local_dofs, ret);
1243  return ret;
1244  }
1245 
1246 
1247 
1248  template <int dim, int spacedim>
1249  template <typename Extractor, typename Number>
1250  const std::vector<
1252  symmetric_gradient_type> &
1254  const std::string &global_vector_name,
1255  const Extractor & variable,
1256  const Number dummy)
1257  {
1258  const std::vector<Number> &independent_local_dofs =
1259  get_local_dof_values(global_vector_name, dummy);
1260 
1261  const FEValuesBase<dim, spacedim> &fev = get_current_fe_values();
1262 
1263  const unsigned int n_q_points = fev.n_quadrature_points;
1264 
1265  const std::string name = get_unique_name(
1266  global_vector_name, variable, "_symmetric_gradient_q", n_q_points, dummy);
1267 
1268 
1269  // Now build the return type
1270  using RetType = std::vector<
1272  symmetric_gradient_type>;
1273 
1274  RetType &ret =
1275  internal_data_storage.template get_or_add_object_with_name<RetType>(
1276  name, n_q_points);
1277 
1278 
1279  AssertDimension(ret.size(), n_q_points);
1280 
1281  fev[variable].get_function_symmetric_gradients_from_local_dof_values(
1282  independent_local_dofs, ret);
1283  return ret;
1284  }
1285 
1286 
1287  template <int dim, int spacedim>
1288  template <typename Extractor, typename Number>
1289  const std::vector<
1291  divergence_type> &
1293  const std::string &global_vector_name,
1294  const Extractor & variable,
1295  const Number dummy)
1296  {
1297  const std::vector<Number> &independent_local_dofs =
1298  get_local_dof_values(global_vector_name, dummy);
1299 
1300  const FEValuesBase<dim, spacedim> &fev = get_current_fe_values();
1301 
1302  const unsigned int n_q_points = fev.n_quadrature_points;
1303 
1304  const std::string name = get_unique_name(
1305  global_vector_name, variable, "_divergence_q", n_q_points, dummy);
1306 
1307  // Now build the return type
1308  using RetType = std::vector<
1310  divergence_type>;
1311 
1312  RetType &ret =
1313  internal_data_storage.template get_or_add_object_with_name<RetType>(
1314  name, n_q_points);
1315 
1316 
1317  AssertDimension(ret.size(), n_q_points);
1318 
1319  fev[variable].get_function_divergences_from_local_dof_values(
1320  independent_local_dofs, ret);
1321  return ret;
1322  }
1323 
1324 
1325 
1326  template <int dim, int spacedim>
1327  template <typename Extractor, typename Number>
1328  const std::vector<
1330  &
1331  ScratchData<dim, spacedim>::get_curls(const std::string &global_vector_name,
1332  const Extractor & variable,
1333  const Number dummy)
1334  {
1335  const std::vector<Number> &independent_local_dofs =
1336  get_local_dof_values(global_vector_name, dummy);
1337 
1338  const FEValuesBase<dim, spacedim> &fev = get_current_fe_values();
1339 
1340  const unsigned int n_q_points = fev.n_quadrature_points;
1341 
1342  const std::string name = get_unique_name(
1343  global_vector_name, variable, "_curl_q", n_q_points, dummy);
1344 
1345  // Now build the return type
1346  using RetType =
1347  std::vector<typename internal::
1348  OutputType<dim, spacedim, Number, Extractor>::curl_type>;
1349 
1350  RetType &ret =
1351  internal_data_storage.template get_or_add_object_with_name<RetType>(
1352  name, n_q_points);
1353 
1354  AssertDimension(ret.size(), n_q_points);
1355 
1356  fev[variable].get_function_curls_from_local_dof_values(
1357  independent_local_dofs, ret);
1358  return ret;
1359  }
1360 
1361 #endif
1362 
1363 } // namespace MeshWorker
1364 
1366 
1367 #endif
typename FEValuesViews::View< dim, spacedim, Extractor >::template OutputType< NumberType > OutputType
Definition: scratch_data.h:47
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
GeneralDataStorage internal_data_storage
Definition: scratch_data.h:949
SmartPointer< const FEValuesBase< dim, spacedim > > current_fe_values
Definition: scratch_data.h:955
std::string type_to_string(const T &t)
Definition: utilities.h:1102
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1623
std::unique_ptr< FEFaceValues< dim, spacedim > > fe_face_values
Definition: scratch_data.h:904
typename ::internal::FEValuesViews::ViewType< dim, spacedim, Extractor >::type View
Definition: fe_values.h:2062
const std::vector< typename internal::OutputType< dim, spacedim, Number, Extractor >::gradient_type > & get_gradients(const std::string &global_vector_name, const Extractor &variable, const Number dummy=Number(0))
Quadrature< dim > cell_quadrature
Definition: scratch_data.h:866
SmartPointer< const Mapping< dim, spacedim > > mapping
Definition: scratch_data.h:854
std::string get_unique_dofs_name(const std::string &global_vector_name, const unsigned int size, const Number &exemplar_number) const
std::unique_ptr< FEValues< dim, spacedim > > neighbor_fe_values
Definition: scratch_data.h:914
UpdateFlags neighbor_face_update_flags
Definition: scratch_data.h:894
SmartPointer< const FEValuesBase< dim, spacedim > > current_neighbor_fe_values
Definition: scratch_data.h:961
const std::vector< typename internal::OutputType< dim, spacedim, Number, Extractor >::value_type > & get_values(const std::string &global_vector_name, const Extractor &variable, const Number dummy=Number(0))
std::vector< double > get_quadrature_points(const unsigned int n)
const std::vector< typename internal::OutputType< dim, spacedim, Number, Extractor >::symmetric_gradient_type > & get_symmetric_gradients(const std::string &global_vector_name, const Extractor &variable, const Number dummy=Number(0))
UpdateFlags neighbor_cell_update_flags
Definition: scratch_data.h:882
std::vector< types::global_dof_index > neighbor_dof_indices
Definition: scratch_data.h:939
const std::vector< typename internal::OutputType< dim, spacedim, Number, Extractor >::hessian_type > & get_hessians(const std::string &global_vector_name, const Extractor &variable, const Number dummy=Number(0))
const std::vector< typename internal::OutputType< dim, spacedim, Number, Extractor >::laplacian_type > & get_laplacians(const std::string &global_vector_name, const Extractor &variable, const Number dummy=Number(0))
std::unique_ptr< FESubfaceValues< dim, spacedim > > fe_subface_values
Definition: scratch_data.h:909
std::unique_ptr< FESubfaceValues< dim, spacedim > > neighbor_fe_subface_values
Definition: scratch_data.h:924
static ::ExceptionBase & ExcMessage(std::string arg1)
No update.
#define Assert(cond, exc)
Definition: exceptions.h:1466
UpdateFlags
Abstract base class for mapping classes.
Definition: mapping.h:303
std::unique_ptr< FEFaceValues< dim, spacedim > > neighbor_fe_face_values
Definition: scratch_data.h:919
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:396
UpdateFlags face_update_flags
Definition: scratch_data.h:888
GeneralDataStorage user_data_storage
Definition: scratch_data.h:944
std::string get_unique_name(const std::string &global_vector_name, const Extractor &variable, const std::string &object_type, const unsigned int size, const Number &exemplar_number) const
const std::vector< typename internal::OutputType< dim, spacedim, Number, Extractor >::third_derivative_type > & get_third_derivatives(const std::string &global_vector_name, const Extractor &variable, const Number dummy=Number(0))
const std::vector< typename internal::OutputType< dim, spacedim, Number, Extractor >::curl_type > & get_curls(const std::string &global_vector_name, const Extractor &variable, const Number dummy=Number(0))
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:473
void extract_local_dof_values(const std::string &global_vector_name, const VectorType &input_vector, const Number dummy=Number(0))
UpdateFlags cell_update_flags
Definition: scratch_data.h:877
const unsigned int n_quadrature_points
Definition: fe_values.h:2186
Quadrature< dim - 1 > face_quadrature
Definition: scratch_data.h:872
std::unique_ptr< FEValues< dim, spacedim > > fe_values
Definition: scratch_data.h:899
const std::vector< typename internal::OutputType< dim, spacedim, Number, Extractor >::divergence_type > & get_divergences(const std::string &global_vector_name, const Extractor &variable, const Number dummy=Number(0))
std::unique_ptr< FEInterfaceValues< dim, spacedim > > interface_fe_values
Definition: scratch_data.h:929
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:395
typename ActiveSelector::active_cell_iterator active_cell_iterator
Definition: dof_handler.h:438
std::vector< types::global_dof_index > local_dof_indices
Definition: scratch_data.h:934
SmartPointer< const FiniteElement< dim, spacedim > > fe
Definition: scratch_data.h:860
const std::vector< Number > & get_local_dof_values(const std::string &global_vector_name, Number dummy=Number(0)) const