Reference documentation for deal.II version Git 4abc4a1666 2020-07-04 19:58:34 +0200
\(\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  // CurrentCellEvaluation
530 
560  template <typename VectorType, typename Number = double>
561  void
562  extract_local_dof_values(const std::string &global_vector_name,
563  const VectorType & input_vector,
564  const Number dummy = Number(0));
565 
574  template <typename Number = double>
575  const std::vector<Number> &
576  get_local_dof_values(const std::string &global_vector_name,
577  Number dummy = Number(0)) const;
578 
600  template <typename Extractor, typename Number = double>
601  const std::vector<
603  value_type> &
604  get_values(const std::string &global_vector_name,
605  const Extractor & variable,
606  const Number dummy = Number(0));
607 
629  template <typename Extractor, typename Number = double>
630  const std::vector<
632  gradient_type> &
633  get_gradients(const std::string &global_vector_name,
634  const Extractor & variable,
635  const Number dummy = Number(0));
636 
659  template <typename Extractor, typename Number = double>
660  const std::vector<
662  symmetric_gradient_type> &
663  get_symmetric_gradients(const std::string &global_vector_name,
664  const Extractor & variable,
665  const Number dummy = Number(0));
666 
688  template <typename Extractor, typename Number = double>
689  const std::vector<
691  divergence_type> &
692  get_divergences(const std::string &global_vector_name,
693  const Extractor & variable,
694  const Number dummy = Number(0));
695 
717  template <typename Extractor, typename Number = double>
718  const std::vector<typename internal::
719  OutputType<dim, spacedim, Number, Extractor>::curl_type>
720  &
721  get_curls(const std::string &global_vector_name,
722  const Extractor & variable,
723  const Number dummy = Number(0));
724 
746  template <typename Extractor, typename Number = double>
747  const std::vector<
749  hessian_type> &
750  get_hessians(const std::string &global_vector_name,
751  const Extractor & variable,
752  const Number dummy = Number(0));
753 
775  template <typename Extractor, typename Number = double>
776  const std::vector<
778  third_derivative_type> &
779  get_third_derivatives(const std::string &global_vector_name,
780  const Extractor & variable,
781  const Number dummy = Number(0));
782  // CurrentCellEvaluation
784 
788  const Mapping<dim, spacedim> &
789  get_mapping() const;
790 
791  private:
796  template <typename Extractor, typename Number = double>
797  std::string
798  get_unique_name(const std::string &global_vector_name,
799  const Extractor & variable,
800  const std::string &object_type,
801  const unsigned int size,
802  const Number & exemplar_number) const;
803 
807  template <typename Number = double>
808  std::string
809  get_unique_dofs_name(const std::string &global_vector_name,
810  const unsigned int size,
811  const Number & exemplar_number) const;
812 
818 
824 
830 
836 
841 
846 
852 
858 
862  std::unique_ptr<FEValues<dim, spacedim>> fe_values;
863 
867  std::unique_ptr<FEFaceValues<dim, spacedim>> fe_face_values;
868 
872  std::unique_ptr<FESubfaceValues<dim, spacedim>> fe_subface_values;
873 
877  std::unique_ptr<FEValues<dim, spacedim>> neighbor_fe_values;
878 
882  std::unique_ptr<FEFaceValues<dim, spacedim>> neighbor_fe_face_values;
883 
887  std::unique_ptr<FESubfaceValues<dim, spacedim>> neighbor_fe_subface_values;
888 
892  std::unique_ptr<FEInterfaceValues<dim, spacedim>> interface_fe_values;
893 
897  std::vector<types::global_dof_index> local_dof_indices;
898 
902  std::vector<types::global_dof_index> neighbor_dof_indices;
903 
908 
913 
919 
925  };
926 
927 #ifndef DOXYGEN
928  template <int dim, int spacedim>
929  template <typename Extractor, typename Number>
930  std::string
932  const std::string &global_vector_name,
933  const Extractor & variable,
934  const std::string &object_type,
935  const unsigned int size,
936  const Number & exemplar_number) const
937  {
938  return global_vector_name + "_" + variable.get_name() + "_" + object_type +
939  "_" + Utilities::int_to_string(size) + "_" +
940  Utilities::type_to_string(exemplar_number);
941  }
942 
943 
944 
945  template <int dim, int spacedim>
946  template <typename Number>
947  std::string
949  const std::string &global_vector_name,
950  const unsigned int size,
951  const Number & exemplar_number) const
952  {
953  return global_vector_name + "_independent_local_dofs_" +
954  Utilities::int_to_string(size) + "_" +
955  Utilities::type_to_string(exemplar_number);
956  }
957 
958 
959 
960  template <int dim, int spacedim>
961  template <typename VectorType, typename Number>
962  void
964  const std::string &global_vector_name,
965  const VectorType & input_vector,
966  const Number dummy)
967  {
968  const unsigned int n_dofs = local_dof_indices.size();
969 
970  const std::string name =
971  get_unique_dofs_name(global_vector_name, n_dofs, dummy);
972 
973  auto &independent_local_dofs =
974  internal_data_storage
975  .template get_or_add_object_with_name<std::vector<Number>>(name,
976  n_dofs);
977 
978  AssertDimension(independent_local_dofs.size(), n_dofs);
979 
981  for (unsigned int i = 0; i < n_dofs; ++i)
983  input_vector(local_dof_indices[i]),
984  i,
985  n_dofs,
986  independent_local_dofs[i]);
987  else
988  for (unsigned int i = 0; i < n_dofs; ++i)
989  independent_local_dofs[i] = input_vector(local_dof_indices[i]);
990  }
991 
992 
993 
994  template <int dim, int spacedim>
995  template <typename Number>
996  const std::vector<Number> &
998  const std::string &global_vector_name,
999  Number dummy) const
1000  {
1001  const unsigned int n_dofs = get_current_fe_values().get_fe().dofs_per_cell;
1002 
1003  const std::string dofs_name =
1004  get_unique_dofs_name(global_vector_name, n_dofs, dummy);
1005 
1006  Assert(
1007  internal_data_storage.stores_object_with_name(dofs_name),
1008  ExcMessage(
1009  "You did not call yet extract_local_dof_values with the right types!"));
1010 
1011  return internal_data_storage
1012  .template get_object_with_name<std::vector<Number>>(dofs_name);
1013  }
1014 
1015 
1016 
1017  template <int dim, int spacedim>
1018  template <typename Extractor, typename Number>
1019  const std::vector<
1021  &
1023  const std::string &global_vector_name,
1024  const Extractor & variable,
1025  const Number dummy)
1026  {
1027  const std::vector<Number> &independent_local_dofs =
1028  get_local_dof_values(global_vector_name, dummy);
1029 
1030  const FEValuesBase<dim, spacedim> &fev = get_current_fe_values();
1031 
1032  const unsigned int n_q_points = fev.n_quadrature_points;
1033 
1034  const std::string name = get_unique_name(
1035  global_vector_name, variable, "_values_q", n_q_points, dummy);
1036 
1037  // Now build the return type
1038  using RetType =
1039  std::vector<typename internal::
1040  OutputType<dim, spacedim, Number, Extractor>::value_type>;
1041 
1042  RetType &ret =
1043  internal_data_storage.template get_or_add_object_with_name<RetType>(
1044  name, n_q_points);
1045 
1046  AssertDimension(ret.size(), n_q_points);
1047 
1048  fev[variable].get_function_values_from_local_dof_values(
1049  independent_local_dofs, ret);
1050  return ret;
1051  }
1052 
1053 
1054 
1055  template <int dim, int spacedim>
1056  template <typename Extractor, typename Number>
1057  const std::vector<
1059  gradient_type> &
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, "_gradients_q", n_q_points, dummy);
1074 
1075  // Now build the return type
1076  using RetType = std::vector<
1078  gradient_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_gradients_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  hessian_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, "_hessians_q", n_q_points, dummy);
1112 
1113  // Now build the return type
1114  using RetType =
1115  std::vector<typename internal::
1116  OutputType<dim, spacedim, Number, Extractor>::hessian_type>;
1117 
1118  RetType &ret =
1119  internal_data_storage.template get_or_add_object_with_name<RetType>(
1120  name, n_q_points);
1121 
1122 
1123  AssertDimension(ret.size(), n_q_points);
1124 
1125  fev[variable].get_function_hessians_from_local_dof_values(
1126  independent_local_dofs, ret);
1127  return ret;
1128  }
1129 
1130 
1131 
1132  template <int dim, int spacedim>
1133  template <typename Extractor, typename Number>
1134  const std::vector<
1136  third_derivative_type> &
1138  const std::string &global_vector_name,
1139  const Extractor & variable,
1140  const Number dummy)
1141  {
1142  const std::vector<Number> &independent_local_dofs =
1143  get_local_dof_values(global_vector_name, dummy);
1144 
1145  const FEValuesBase<dim, spacedim> &fev = get_current_fe_values();
1146 
1147  const unsigned int n_q_points = fev.n_quadrature_points;
1148 
1149  const std::string name = get_unique_name(
1150  global_vector_name, variable, "_third_derivatives_q", n_q_points, dummy);
1151 
1152  // Now build the return type
1153  using RetType = std::vector<
1155  third_derivative_type>;
1156 
1157  RetType &ret =
1158  internal_data_storage.template get_or_add_object_with_name<RetType>(
1159  name, n_q_points);
1160 
1161 
1162  AssertDimension(ret.size(), n_q_points);
1163 
1164  fev[variable].get_function_third_derivatives_from_local_dof_values(
1165  independent_local_dofs, ret);
1166  return ret;
1167  }
1168 
1169 
1170 
1171  template <int dim, int spacedim>
1172  template <typename Extractor, typename Number>
1173  const std::vector<
1175  symmetric_gradient_type> &
1177  const std::string &global_vector_name,
1178  const Extractor & variable,
1179  const Number dummy)
1180  {
1181  const std::vector<Number> &independent_local_dofs =
1182  get_local_dof_values(global_vector_name, dummy);
1183 
1184  const FEValuesBase<dim, spacedim> &fev = get_current_fe_values();
1185 
1186  const unsigned int n_q_points = fev.n_quadrature_points;
1187 
1188  const std::string name = get_unique_name(
1189  global_vector_name, variable, "_symmetric_gradient_q", n_q_points, dummy);
1190 
1191 
1192  // Now build the return type
1193  using RetType = std::vector<
1195  symmetric_gradient_type>;
1196 
1197  RetType &ret =
1198  internal_data_storage.template get_or_add_object_with_name<RetType>(
1199  name, n_q_points);
1200 
1201 
1202  AssertDimension(ret.size(), n_q_points);
1203 
1204  fev[variable].get_function_symmetric_gradients_from_local_dof_values(
1205  independent_local_dofs, ret);
1206  return ret;
1207  }
1208 
1209 
1210  template <int dim, int spacedim>
1211  template <typename Extractor, typename Number>
1212  const std::vector<
1214  divergence_type> &
1216  const std::string &global_vector_name,
1217  const Extractor & variable,
1218  const Number dummy)
1219  {
1220  const std::vector<Number> &independent_local_dofs =
1221  get_local_dof_values(global_vector_name, dummy);
1222 
1223  const FEValuesBase<dim, spacedim> &fev = get_current_fe_values();
1224 
1225  const unsigned int n_q_points = fev.n_quadrature_points;
1226 
1227  const std::string name = get_unique_name(
1228  global_vector_name, variable, "_divergence_q", n_q_points, dummy);
1229 
1230  // Now build the return type
1231  using RetType = std::vector<
1233  divergence_type>;
1234 
1235  RetType &ret =
1236  internal_data_storage.template get_or_add_object_with_name<RetType>(
1237  name, n_q_points);
1238 
1239 
1240  AssertDimension(ret.size(), n_q_points);
1241 
1242  fev[variable].get_function_divergences_from_local_dof_values(
1243  independent_local_dofs, ret);
1244  return ret;
1245  }
1246 
1247 
1248 
1249  template <int dim, int spacedim>
1250  template <typename Extractor, typename Number>
1251  const std::vector<
1253  &
1254  ScratchData<dim, spacedim>::get_curls(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, "_curl_q", n_q_points, dummy);
1267 
1268  // Now build the return type
1269  using RetType =
1270  std::vector<typename internal::
1271  OutputType<dim, spacedim, Number, Extractor>::curl_type>;
1272 
1273  RetType &ret =
1274  internal_data_storage.template get_or_add_object_with_name<RetType>(
1275  name, n_q_points);
1276 
1277  AssertDimension(ret.size(), n_q_points);
1278 
1279  fev[variable].get_function_curls_from_local_dof_values(
1280  independent_local_dofs, ret);
1281  return ret;
1282  }
1283 
1284 #endif
1285 
1286 } // namespace MeshWorker
1287 
1289 
1290 #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:912
SmartPointer< const FEValuesBase< dim, spacedim > > current_fe_values
Definition: scratch_data.h:918
std::string type_to_string(const T &t)
Definition: utilities.h:1066
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1560
std::unique_ptr< FEFaceValues< dim, spacedim > > fe_face_values
Definition: scratch_data.h:867
typename ::internal::FEValuesViews::ViewType< dim, spacedim, Extractor >::type View
Definition: fe_values.h:1970
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:829
SmartPointer< const Mapping< dim, spacedim > > mapping
Definition: scratch_data.h:817
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:877
UpdateFlags neighbor_face_update_flags
Definition: scratch_data.h:857
SmartPointer< const FEValuesBase< dim, spacedim > > current_neighbor_fe_values
Definition: scratch_data.h:924
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:845
std::vector< types::global_dof_index > neighbor_dof_indices
Definition: scratch_data.h:902
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))
std::unique_ptr< FESubfaceValues< dim, spacedim > > fe_subface_values
Definition: scratch_data.h:872
std::unique_ptr< FESubfaceValues< dim, spacedim > > neighbor_fe_subface_values
Definition: scratch_data.h:887
static ::ExceptionBase & ExcMessage(std::string arg1)
No update.
#define Assert(cond, exc)
Definition: exceptions.h:1403
UpdateFlags
Abstract base class for mapping classes.
Definition: mapping.h:301
std::unique_ptr< FEFaceValues< dim, spacedim > > neighbor_fe_face_values
Definition: scratch_data.h:882
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
UpdateFlags face_update_flags
Definition: scratch_data.h:851
GeneralDataStorage user_data_storage
Definition: scratch_data.h:907
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:474
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:840
const unsigned int n_quadrature_points
Definition: fe_values.h:2090
Quadrature< dim - 1 > face_quadrature
Definition: scratch_data.h:835
std::unique_ptr< FEValues< dim, spacedim > > fe_values
Definition: scratch_data.h:862
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:892
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
typename ActiveSelector::active_cell_iterator active_cell_iterator
Definition: dof_handler.h:343
std::vector< types::global_dof_index > local_dof_indices
Definition: scratch_data.h:897
SmartPointer< const FiniteElement< dim, spacedim > > fe
Definition: scratch_data.h:823
const std::vector< Number > & get_local_dof_values(const std::string &global_vector_name, Number dummy=Number(0)) const