Loading [MathJax]/extensions/TeX/AMSsymbols.js
 deal.II version GIT relicensing-3012-g95e643124d 2025-04-03 05:00:00+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\}}\)
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
fe_coupling_values.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2023 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_fe_coupling_values_h
16#define dealii_fe_coupling_values_h
17
18#include <deal.II/base/config.h>
19
21
23#include <deal.II/base/point.h>
27
29
31
32// Forward declaration
33#ifndef DOXYGEN
34template <int dim, int spacedim>
35class FEValuesBase;
36#endif
37
39{
51 {
92 const unsigned int n_inner_quadrature_points =
94 const std::vector<unsigned int> &dof_renumbering = {},
95 const std::vector<unsigned int> &quadrature_renumbering = {});
96
100 RenumberingData(const RenumberingData &other) = delete;
101
105 const unsigned int n_inner_dofs;
106
110 const unsigned int n_dofs;
111
116 const unsigned int n_inner_quadrature_points;
117
121 const unsigned int n_quadrature_points;
122
126 const std::vector<unsigned int> dof_renumbering;
127
131 const std::vector<unsigned int> quadrature_renumbering;
132
145 };
146
168 template <typename ViewType>
170 {
171 public:
175 using value_type = typename ViewType::value_type;
176
180 using gradient_type = typename ViewType::gradient_type;
181
186 template <typename Number>
188 typename ViewType::template solution_value_type<Number>;
189
194 template <typename Number>
196 typename ViewType::template solution_gradient_type<Number>;
197
211 RenumberedView(const ViewType &view, const RenumberingData &data);
212
230 value(const unsigned int shape_function, const unsigned int q_point) const;
231
249 gradient(const unsigned int shape_function,
250 const unsigned int q_point) const;
251
266 template <typename Number>
267 void
269 std::vector<solution_value_type<Number>> &values) const;
270
291 template <class InputVector>
292 void
294 const InputVector &dof_values,
296 &values) const;
297
313 template <typename Number>
314 void
316 const ReadVector<Number> &fe_function,
317 std::vector<solution_gradient_type<Number>> &gradients) const;
318
339 template <class InputVector>
340 void
342 const InputVector &dof_values,
344 &gradients) const;
345
346 private:
352
353
370 template <typename Number>
371 std::string
372 get_unique_container_name(const std::string &prefix,
373 const unsigned int size,
374 const Number &exemplar_number) const;
375
380 template <typename InputVector>
381 const InputVector &
382 outer_to_inner_dofs(const InputVector &outer_vector) const;
383
388 template <typename ValueType>
389 std::vector<ValueType> &
390 outer_to_inner_values(std::vector<ValueType> &outer_values) const;
391
397 template <typename ValueType>
398 void
399 inner_to_outer_values(const std::vector<ValueType> &inner_values,
400 std::vector<ValueType> &outer_values) const;
401
405 const ViewType &view;
406 };
407} // namespace FEValuesViews
408
409
435{
444
455 unrolled,
456
467 matching,
468
477 reorder,
478
489};
490
535{
543
557};
558
730template <int dim1, int dim2 = dim1, int spacedim = dim1>
732{
733public:
739
776 const FEValuesBase<dim1, spacedim> &fe_values_1,
777 const FEValuesBase<dim2, spacedim> &fe_values_2,
781
816 void
818 const FEValuesBase<dim1, spacedim> &fe_values_1,
819 const FEValuesBase<dim2, spacedim> &fe_values_2,
823
824
829 template <typename Extractor>
831 get_first_extractor(const Extractor &extractor) const;
832
837 template <typename Extractor>
839 get_second_extractor(const Extractor &extractor) const;
840
846 double
847 JxW(const unsigned int quadrature_point) const;
848
856 std::pair<Point<spacedim>, Point<spacedim>>
857 quadrature_point(const unsigned int quadrature_point) const;
858
868
876
884
899
940 std::vector<types::global_dof_index>
942 const std::vector<types::global_dof_index> &dof_indices_1,
943 const std::vector<types::global_dof_index> &dof_indices_2,
944 const types::global_dof_index dofs_offset_1 = 0,
945 const types::global_dof_index dofs_offset_2 = 0) const;
946
952 std::pair<unsigned int, unsigned int>
953 coupling_dof_to_dof_indices(const unsigned int coupling_dof_index) const;
954
960 std::pair<unsigned int, unsigned int>
962 const unsigned int quadrature_point) const;
963
974 template <typename Extractor>
978 const FEValuesExtractors::FirstCoupling<Extractor> &extractor) const;
979
985 template <typename Extractor>
989 const FEValuesExtractors::SecondCoupling<Extractor> &extractor) const;
990
1000 unsigned int
1002
1007 unsigned int
1009
1014 unsigned int
1016
1020 unsigned int
1022
1023private:
1028
1033
1038
1043
1047 std::unique_ptr<const FEValuesViews::RenumberingData> first_renumbering_data;
1048
1052 std::unique_ptr<const FEValuesViews::RenumberingData> second_renumbering_data;
1053
1058
1064 unsigned int n_coupling_dofs_;
1065};
1066
1067
1068#ifndef DOXYGEN
1069
1070
1071/*------------------------ Inline functions: namespace FEValuesViews --------*/
1072
1073namespace FEValuesViews
1074{
1076 const unsigned int n_inner_dofs,
1077 const unsigned int n_inner_quadrature_points,
1078 const std::vector<unsigned int> &dof_renumbering,
1079 const std::vector<unsigned int> &quadrature_renumbering)
1080 : n_inner_dofs(n_inner_dofs)
1081 , n_dofs(dof_renumbering.empty() ? n_inner_dofs : dof_renumbering.size())
1082 , n_inner_quadrature_points(n_inner_quadrature_points)
1083 , n_quadrature_points(quadrature_renumbering.empty() ?
1084 n_inner_quadrature_points :
1085 quadrature_renumbering.size())
1086 , dof_renumbering(dof_renumbering)
1087 , quadrature_renumbering(quadrature_renumbering)
1088 {
1089 // Check that the renumbering vectors are valid.
1090 if constexpr (running_in_debug_mode())
1091 {
1092 // While for dofs we admit invalid values, this is not the case for
1093 // quadrature points.
1094 for (const auto i : dof_renumbering)
1095 Assert(i < n_inner_dofs || i == numbers::invalid_unsigned_int,
1096 ExcIndexRange(i, 0, n_inner_dofs));
1097
1098 for (const auto q : quadrature_renumbering)
1099 AssertIndexRange(q, n_inner_quadrature_points);
1100 }
1101 }
1102
1103
1104
1105 template <typename ViewType>
1106 RenumberedView<ViewType>::RenumberedView(const ViewType &view,
1107 const RenumberingData &data)
1108 : view(view)
1109 , data(data)
1110 {}
1111
1112
1113
1114 template <typename ViewType>
1115 template <typename Number>
1116 inline std::string
1117 RenumberedView<ViewType>::get_unique_container_name(
1118 const std::string &prefix,
1119 const unsigned int size,
1120 const Number &exemplar_number) const
1121 {
1122 return prefix + "_" + Utilities::int_to_string(size) + "_" +
1123 Utilities::type_to_string(exemplar_number);
1124 }
1125
1126
1127
1128 template <typename ViewType>
1129 typename RenumberedView<ViewType>::value_type
1130 RenumberedView<ViewType>::value(const unsigned int shape_function,
1131 const unsigned int q_point) const
1132 {
1133 AssertIndexRange(shape_function, data.n_dofs);
1134 AssertIndexRange(q_point, data.n_quadrature_points);
1135
1136 const auto inner_shape_function = data.dof_renumbering.empty() ?
1137 shape_function :
1138 data.dof_renumbering[shape_function];
1139 const auto inner_q_point = data.quadrature_renumbering.empty() ?
1140 q_point :
1141 data.quadrature_renumbering[q_point];
1142 if (inner_shape_function == numbers::invalid_unsigned_int)
1143 return value_type(0);
1144 else
1145 {
1146 AssertIndexRange(inner_shape_function, data.n_inner_dofs);
1147 AssertIndexRange(inner_q_point, data.n_inner_quadrature_points);
1148 return view.value(inner_shape_function, inner_q_point);
1149 }
1150 }
1151
1152
1153
1154 template <typename ViewType>
1155 typename RenumberedView<ViewType>::gradient_type
1156 RenumberedView<ViewType>::gradient(const unsigned int shape_function,
1157 const unsigned int q_point) const
1158 {
1159 AssertIndexRange(shape_function, data.n_dofs);
1160 AssertIndexRange(q_point, data.n_quadrature_points);
1161
1162 const auto inner_shape_function = data.dof_renumbering.empty() ?
1163 shape_function :
1164 data.dof_renumbering[shape_function];
1165 const auto inner_q_point = data.quadrature_renumbering.empty() ?
1166 q_point :
1167 data.quadrature_renumbering[q_point];
1168 if (inner_shape_function == numbers::invalid_unsigned_int)
1169 return gradient_type();
1170 else
1171 return view.gradient(inner_shape_function, inner_q_point);
1172 }
1173
1174
1175
1176 template <typename ViewType>
1177 template <typename ValueType>
1178 std::vector<ValueType> &
1179 RenumberedView<ViewType>::outer_to_inner_values(
1180 std::vector<ValueType> &outer_values) const
1181 {
1182 AssertDimension(outer_values.size(), data.n_quadrature_points);
1183 if (data.quadrature_renumbering.empty())
1184 {
1185 return outer_values;
1186 }
1187 else
1188 {
1189 const auto name =
1190 get_unique_container_name("RenumberedView::outer_to_inner_values",
1191 data.n_inner_quadrature_points,
1192 outer_values[0]);
1193 auto &inner_values =
1194 data.data_storage.get()
1195 .template get_or_add_object_with_name<std::vector<ValueType>>(
1196 name, data.n_inner_quadrature_points);
1197 return inner_values;
1198 }
1199 }
1200
1201
1202
1203 template <typename ViewType>
1204 template <typename VectorType>
1205 const VectorType &
1206 RenumberedView<ViewType>::outer_to_inner_dofs(
1207 const VectorType &outer_dofs) const
1208 {
1209 AssertDimension(outer_dofs.size(), data.n_dofs);
1210 if (data.dof_renumbering.empty())
1211 {
1212 return outer_dofs;
1213 }
1214 else
1215 {
1216 const auto name =
1217 get_unique_container_name("RenumberedView::outer_to_inner_dofs",
1218 data.n_inner_dofs,
1219 outer_dofs[0]);
1220
1221 auto &inner_dofs = data.data_storage.get()
1222 .template get_or_add_object_with_name<VectorType>(
1223 name, data.n_inner_dofs);
1224 for (unsigned int i = 0; i < data.n_dofs; ++i)
1225 {
1226 const auto inner_i = data.dof_renumbering[i];
1227 if (inner_i != numbers::invalid_unsigned_int)
1228 {
1229 AssertIndexRange(inner_i, data.n_inner_dofs);
1230 inner_dofs[inner_i] = outer_dofs[i];
1231 }
1232 }
1233 return inner_dofs;
1234 }
1235 }
1236
1237
1238
1239 template <typename ViewType>
1240 template <typename ValueType>
1241 void
1242 RenumberedView<ViewType>::inner_to_outer_values(
1243 const std::vector<ValueType> &inner_values,
1244 std::vector<ValueType> &outer_values) const
1245 {
1246 AssertDimension(outer_values.size(), data.n_quadrature_points);
1247 AssertDimension(inner_values.size(), data.n_inner_quadrature_points);
1248 if (data.quadrature_renumbering.empty())
1249 {
1250 Assert(&inner_values == &outer_values, ExcInternalError());
1251 return;
1252 }
1253 for (unsigned int i = 0; i < data.quadrature_renumbering.size(); ++i)
1254 {
1255 outer_values[i] = inner_values[data.quadrature_renumbering[i]];
1256 }
1257 }
1258
1259
1260
1261 template <typename ViewType>
1262 template <typename Number>
1263 void
1264 RenumberedView<ViewType>::get_function_values(
1265 const ReadVector<Number> &fe_function,
1266 std::vector<solution_value_type<Number>> &values) const
1267 {
1268 auto &inner_values = outer_to_inner_values(values);
1269 view.get_function_values(fe_function, inner_values);
1270 inner_to_outer_values(inner_values, values);
1271 }
1272
1273
1274
1275 template <typename ViewType>
1276 template <typename InputVector>
1277 void
1278 RenumberedView<ViewType>::get_function_values_from_local_dof_values(
1279 const InputVector &dof_values,
1280 std::vector<solution_value_type<typename InputVector::value_type>> &values)
1281 const
1282 {
1283 const auto &inner_dof_values = outer_to_inner_dofs(dof_values);
1284 auto &inner_values = outer_to_inner_values(values);
1285
1286 view.get_function_values_from_local_dof_values(inner_dof_values,
1287 inner_values);
1288 inner_to_outer_values(inner_values, values);
1289 }
1290
1291
1292 template <typename ViewType>
1293 template <typename Number>
1294 void
1295 RenumberedView<ViewType>::get_function_gradients(
1296 const ReadVector<Number> &fe_function,
1297 std::vector<solution_gradient_type<Number>> &gradients) const
1298 {
1299 auto &inner_gradients = outer_to_inner_values(gradients);
1300 view.get_function_gradients(fe_function, inner_gradients);
1301 inner_to_outer_values(inner_gradients, gradients);
1302 }
1303
1304
1305
1306 template <typename ViewType>
1307 template <typename InputVector>
1308 void
1309 RenumberedView<ViewType>::get_function_gradients_from_local_dof_values(
1310 const InputVector &dof_values,
1311 std::vector<solution_gradient_type<typename InputVector::value_type>>
1312 &gradients) const
1313 {
1314 const auto &inner_dof_values = outer_to_inner_dofs(dof_values);
1315 auto &inner_gradients = outer_to_inner_values(gradients);
1316
1317 view.get_function_gradients_from_local_dof_values(inner_dof_values,
1318 inner_gradients);
1319 inner_to_outer_values(inner_gradients, gradients);
1320 }
1321} // namespace FEValuesViews
1322
1323/*-------------- Inline functions FECouplingValues ---------------------*/
1324
1325template <int dim1, int dim2, int spacedim>
1326std::vector<types::global_dof_index>
1328 const std::vector<types::global_dof_index> &dof_indices_1,
1329 const std::vector<types::global_dof_index> &dof_indices_2,
1330 const types::global_dof_index dofs_offset_1,
1331 const types::global_dof_index dofs_offset_2) const
1332{
1334 n_coupling_dofs_ != numbers::invalid_unsigned_int,
1335 ExcMessage(
1336 "Dofs are independent. You cannot ask for coupling dof indices."));
1337 AssertDimension(dof_indices_1.size(), first_fe_values->dofs_per_cell);
1338 AssertDimension(dof_indices_2.size(), second_fe_values->dofs_per_cell);
1339
1340 std::vector<types::global_dof_index> coupling_dof_indices(
1341 dof_indices_1.size() + dof_indices_2.size());
1342 unsigned int idx = 0;
1343 for (const auto &i : dof_indices_1)
1344 coupling_dof_indices[idx++] = i + dofs_offset_1;
1345 for (const auto &i : dof_indices_2)
1346 coupling_dof_indices[idx++] = i + dofs_offset_2;
1347 return coupling_dof_indices;
1348}
1349
1350
1351
1352template <int dim1, int dim2, int spacedim>
1354 : first_fe_values(nullptr)
1355 , second_fe_values(nullptr)
1356 , quadrature_coupling_type(QuadratureCouplingType::unrolled)
1357 , n_quadrature_points_(0)
1358 , n_coupling_dofs_(numbers::invalid_unsigned_int)
1359{}
1360
1361
1362
1363template <int dim1, int dim2, int spacedim>
1365 const FEValuesBase<dim1, spacedim> &fe_values_1,
1366 const FEValuesBase<dim2, spacedim> &fe_values_2,
1367 const DoFCouplingType &dof_coupling_type,
1368 const QuadratureCouplingType &quadrature_coupling_type)
1369 : FECouplingValues() // delegate to other constructor
1370{
1371 reinit(fe_values_1, fe_values_2, dof_coupling_type, quadrature_coupling_type);
1372}
1373
1374
1375
1376template <int dim1, int dim2, int spacedim>
1377void
1379 const FEValuesBase<dim1, spacedim> &fe_values_1,
1380 const FEValuesBase<dim2, spacedim> &fe_values_2,
1381 const DoFCouplingType &dof_coupling_type,
1382 const QuadratureCouplingType &quadrature_coupling_type)
1383{
1384 first_fe_values = &fe_values_1;
1385 second_fe_values = &fe_values_2;
1386 this->dof_coupling_type = dof_coupling_type;
1387 this->quadrature_coupling_type = quadrature_coupling_type;
1388
1389 // Gather information about the inner objects
1390 unsigned int first_n_inner_dofs = fe_values_1.dofs_per_cell;
1391 unsigned int second_n_inner_dofs = fe_values_2.dofs_per_cell;
1392
1393 unsigned int first_n_inner_quadrature_points =
1394 fe_values_1.n_quadrature_points;
1395 unsigned int second_n_inner_quadrature_points =
1396 fe_values_2.n_quadrature_points;
1397
1398 // Initialize counters and renumbering vectors to zero
1399 std::vector<unsigned int> first_dofs_map;
1400 std::vector<unsigned int> second_dofs_map;
1401
1402 std::vector<unsigned int> first_quad_map;
1403 std::vector<unsigned int> second_quad_map;
1404
1405 // Local relative tolerance for comparing quadrature points
1406 const double tol = std::min(fe_values_1.get_cell()->diameter(),
1407 fe_values_2.get_cell()->diameter()) *
1408 1e-6;
1409
1410 // Now fill the renumbering vectors
1411 switch (dof_coupling_type)
1412 {
1414 {
1415 n_coupling_dofs_ = numbers::invalid_unsigned_int;
1416 first_dofs_map.clear();
1417 second_dofs_map.clear();
1418 break;
1419 }
1421 {
1422 n_coupling_dofs_ =
1423 fe_values_1.dofs_per_cell + fe_values_2.dofs_per_cell;
1424
1425 first_dofs_map.resize(n_coupling_dofs_);
1426 second_dofs_map.resize(n_coupling_dofs_);
1427
1428 unsigned int idx = 0;
1429 for (const unsigned int &i : fe_values_1.dof_indices())
1430 {
1431 first_dofs_map[idx] = i;
1432 second_dofs_map[idx++] = numbers::invalid_unsigned_int;
1433 }
1434 for (const unsigned int &i : fe_values_2.dof_indices())
1435 {
1436 first_dofs_map[idx] = numbers::invalid_unsigned_int;
1437 second_dofs_map[idx++] = i;
1438 }
1439 // Make sure we have the right number of dofs
1440 AssertDimension(idx, n_coupling_dofs_);
1441 break;
1442 }
1443 default:
1445 }
1446 // Compute the quadrature maps
1447 switch (quadrature_coupling_type)
1448 {
1450 {
1451 const auto &quadrature_points_1 = fe_values_1.get_quadrature_points();
1452 const auto &quadrature_points_2 = fe_values_2.get_quadrature_points();
1453
1454 n_quadrature_points_ =
1455 quadrature_points_1.size() * quadrature_points_2.size();
1456
1457 first_quad_map.resize(n_quadrature_points_);
1458 second_quad_map.resize(n_quadrature_points_);
1459
1460 unsigned int idx = 0;
1461 for (const unsigned int &i : fe_values_1.quadrature_point_indices())
1462 for (const unsigned int &j : fe_values_2.quadrature_point_indices())
1463 {
1464 first_quad_map[idx] = i;
1465 second_quad_map[idx] = j;
1466 ++idx;
1467 }
1468 break;
1469 }
1471 {
1472 Assert(fe_values_1.get_quadrature_points().size() ==
1473 fe_values_2.get_quadrature_points().size(),
1474 ExcMessage("The two FEValuesBase objects must have the same "
1475 "number of quadrature points"));
1476
1477 n_quadrature_points_ = fe_values_1.get_quadrature_points().size();
1478
1479 first_quad_map.clear();
1480 second_quad_map.clear();
1481
1482 break;
1483 }
1485 {
1486 const auto &quadrature_points_1 = fe_values_1.get_quadrature_points();
1487 const auto &quadrature_points_2 = fe_values_2.get_quadrature_points();
1488
1489 Assert(quadrature_points_1.size() == quadrature_points_2.size(),
1490 ExcMessage("The two FEValuesBase objects must have the same "
1491 "number of quadrature points"));
1492
1493 for (const unsigned int &i : fe_values_1.quadrature_point_indices())
1494 {
1495 Assert(quadrature_points_1[i].distance(quadrature_points_2[i]) <
1496 tol,
1497 ExcMessage(
1498 "The two FEValuesBase objects must have the same "
1499 "quadrature points"));
1500 }
1501
1502 n_quadrature_points_ = fe_values_1.get_quadrature_points().size();
1503
1504 first_quad_map.clear();
1505 second_quad_map.clear();
1506
1507 break;
1508 }
1510 {
1511 const auto &quadrature_points_1 = fe_values_1.get_quadrature_points();
1512 const auto &quadrature_points_2 = fe_values_2.get_quadrature_points();
1513
1514 Assert(quadrature_points_1.size() == quadrature_points_2.size(),
1515 ExcMessage("The two FEValuesBase objects must have the same "
1516 "number of quadrature points"));
1517
1518 n_quadrature_points_ = fe_values_1.get_quadrature_points().size();
1519
1520 // The first is the id. The second is renumbered.
1521 first_quad_map.clear();
1522 second_quad_map.resize(fe_values_1.get_quadrature_points().size());
1523
1524 // [TODO]: Avoid quadratic complexity here
1525 for (const unsigned int &i : fe_values_1.quadrature_point_indices())
1526 {
1528 for (const unsigned int &j :
1529 fe_values_2.quadrature_point_indices())
1530 if (quadrature_points_1[i].distance(quadrature_points_2[j]) <
1531 tol)
1532 {
1533 id = i;
1534 second_quad_map[i] = j;
1535 break;
1536 }
1538 ExcMessage(
1539 "The two FEValuesBase objects must have the same "
1540 "quadrature points, even if not in the same order."));
1541 }
1542 break;
1543 }
1545 {
1546 const auto &quadrature_points_1 = fe_values_1.get_quadrature_points();
1547 const auto &quadrature_points_2 = fe_values_2.get_quadrature_points();
1548
1549 // [TODO]: Avoid quadratic complexity here
1550 for (const unsigned int &i : fe_values_1.quadrature_point_indices())
1551 {
1552 for (const unsigned int &j :
1553 fe_values_2.quadrature_point_indices())
1554 if (quadrature_points_1[i].distance(quadrature_points_2[j]) <
1555 tol)
1556 {
1557 first_quad_map.emplace_back(i);
1558 second_quad_map.emplace_back(j);
1559 break;
1560 }
1561 }
1562 n_quadrature_points_ = first_quad_map.size();
1563
1564 break;
1565 }
1566 default:
1567 Assert(false, ExcInternalError());
1568 }
1569
1570 first_renumbering_data = std::make_unique<FEValuesViews::RenumberingData>(
1571 first_n_inner_dofs,
1572 first_n_inner_quadrature_points,
1573 first_dofs_map,
1574 first_quad_map);
1575
1576 second_renumbering_data = std::make_unique<FEValuesViews::RenumberingData>(
1577 second_n_inner_dofs,
1578 second_n_inner_quadrature_points,
1579 second_dofs_map,
1580 second_quad_map);
1581}
1582
1583
1584template <int dim1, int dim2, int spacedim>
1585inline double
1586FECouplingValues<dim1, dim2, spacedim>::JxW(const unsigned int q) const
1587{
1588 AssertIndexRange(q, n_quadrature_points_);
1589
1590 const auto first_q = first_renumbering_data->quadrature_renumbering.empty() ?
1591 q :
1592 first_renumbering_data->quadrature_renumbering[q];
1593
1594 if (quadrature_coupling_type == QuadratureCouplingType::tensor_product ||
1595 quadrature_coupling_type == QuadratureCouplingType::unrolled)
1596 {
1597 const auto second_q =
1598 second_renumbering_data->quadrature_renumbering.empty() ?
1599 q :
1600 second_renumbering_data->quadrature_renumbering[q];
1601 return first_fe_values->JxW(first_q) * second_fe_values->JxW(second_q);
1602 }
1603 else
1604 return first_fe_values->JxW(first_q);
1605}
1606
1607
1608
1609template <int dim1, int dim2, int spacedim>
1612{
1614 0U, n_quadrature_points_);
1615}
1616
1617
1618
1619template <int dim1, int dim2, int spacedim>
1622{
1623 AssertThrow(n_coupling_dofs_ != numbers::invalid_unsigned_int,
1624 ExcMessage(
1625 "Dofs are independent. You cannot ask for coupling dofs."));
1627 0U, n_coupling_dofs_);
1628}
1629
1630
1631
1632template <int dim1, int dim2, int spacedim>
1635{
1637 0U, n_first_dofs());
1638}
1639
1640
1641
1642template <int dim1, int dim2, int spacedim>
1645{
1647 0U, n_second_dofs());
1648}
1649
1650
1651
1652template <int dim1, int dim2, int spacedim>
1653inline unsigned int
1655{
1656 AssertThrow(n_coupling_dofs_ != numbers::invalid_unsigned_int,
1657 ExcMessage(
1658 "Dofs are independent. You cannot ask for coupling dofs."));
1659 return n_coupling_dofs_;
1660}
1661
1662
1663
1664template <int dim1, int dim2, int spacedim>
1665inline unsigned int
1667{
1668 return first_renumbering_data->n_dofs;
1669}
1670
1671
1672
1673template <int dim1, int dim2, int spacedim>
1674inline unsigned int
1676{
1677 return second_renumbering_data->n_dofs;
1678}
1679
1680
1681
1682template <int dim1, int dim2, int spacedim>
1683inline unsigned int
1685{
1686 return n_quadrature_points_;
1687}
1688
1689
1690
1691template <int dim1, int dim2, int spacedim>
1692std::pair<Point<spacedim>, Point<spacedim>>
1694 const unsigned int quadrature_point) const
1695{
1696 AssertIndexRange(quadrature_point, n_quadrature_points_);
1697 const auto first_q = first_fe_values->quadrature_point(
1698 first_renumbering_data->quadrature_renumbering.empty() ?
1699 quadrature_point :
1700 first_renumbering_data->quadrature_renumbering[quadrature_point]);
1701
1702 const auto second_q = second_fe_values->quadrature_point(
1703 second_renumbering_data->quadrature_renumbering.empty() ?
1704 quadrature_point :
1705 second_renumbering_data->quadrature_renumbering[quadrature_point]);
1706
1707 return {first_q, second_q};
1708}
1709
1710
1711
1712template <int dim1, int dim2, int spacedim>
1713std::pair<unsigned int, unsigned int>
1716 const unsigned int quadrature_point) const
1717{
1718 AssertIndexRange(quadrature_point, n_quadrature_points_);
1719 const auto first_id =
1720 first_renumbering_data->quadrature_renumbering.empty() ?
1721 quadrature_point :
1722 first_renumbering_data->quadrature_renumbering[quadrature_point];
1723
1724 const auto second_id =
1725 second_renumbering_data->quadrature_renumbering.empty() ?
1726 quadrature_point :
1727 second_renumbering_data->quadrature_renumbering[quadrature_point];
1728 return std::make_pair(first_id, second_id);
1729}
1730
1731
1732
1733template <int dim1, int dim2, int spacedim>
1734std::pair<unsigned int, unsigned int>
1736 const unsigned int coupling_dof_index) const
1737{
1738 AssertIndexRange(coupling_dof_index, n_coupling_dofs_);
1739 const auto first_id =
1740 first_renumbering_data->dof_renumbering.empty() ?
1741 coupling_dof_index :
1742 first_renumbering_data->dof_renumbering[coupling_dof_index];
1743
1744 const auto second_id =
1745 second_renumbering_data->dof_renumbering.empty() ?
1746 coupling_dof_index :
1747 second_renumbering_data->dof_renumbering[coupling_dof_index];
1748 return std::make_pair(first_id, second_id);
1749}
1750
1751
1752
1753template <int dim1, int dim2, int spacedim>
1754template <typename Extractor>
1757 const Extractor &extractor) const
1758{
1760}
1761
1762
1763
1764template <int dim1, int dim2, int spacedim>
1765template <typename Extractor>
1768 const Extractor &extractor) const
1769{
1771}
1772
1773
1774
1775template <int dim1, int dim2, int spacedim>
1776template <typename Extractor>
1780 const FEValuesExtractors::FirstCoupling<Extractor> &extractor) const
1781{
1784 (*first_fe_values)[extractor.extractor], *first_renumbering_data);
1785}
1786
1787
1788
1789template <int dim1, int dim2, int spacedim>
1790template <typename Extractor>
1794 const FEValuesExtractors::SecondCoupling<Extractor> &extractor) const
1795{
1798 (*second_fe_values)[extractor.extractor], *second_renumbering_data);
1799}
1800
1801#endif // DOXYGEN
1802
1804
1805#endif
unsigned int n_quadrature_points() const
std::pair< Point< spacedim >, Point< spacedim > > quadrature_point(const unsigned int quadrature_point) const
const FEValuesViews::RenumberedView< FEValuesViews::View< dim2, spacedim, Extractor > > operator[](const FEValuesExtractors::SecondCoupling< Extractor > &extractor) const
std::vector< types::global_dof_index > get_coupling_dof_indices(const std::vector< types::global_dof_index > &dof_indices_1, const std::vector< types::global_dof_index > &dof_indices_2, const types::global_dof_index dofs_offset_1=0, const types::global_dof_index dofs_offset_2=0) const
FEValuesExtractors::FirstCoupling< Extractor > get_first_extractor(const Extractor &extractor) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > first_dof_indices() const
FECouplingValues(const FEValuesBase< dim1, spacedim > &fe_values_1, const FEValuesBase< dim2, spacedim > &fe_values_2, const DoFCouplingType &dof_coupling_type=DoFCouplingType::independent, const QuadratureCouplingType &quadrature_coupling_type=QuadratureCouplingType::tensor_product)
unsigned int n_quadrature_points_
std_cxx20::ranges::iota_view< unsigned int, unsigned int > quadrature_point_indices() const
unsigned int n_coupling_dofs_
unsigned int n_coupling_dofs() const
QuadratureCouplingType quadrature_coupling_type
const FEValuesViews::RenumberedView< FEValuesViews::View< dim1, spacedim, Extractor > > operator[](const FEValuesExtractors::FirstCoupling< Extractor > &extractor) const
std::pair< unsigned int, unsigned int > coupling_quadrature_to_quadrature_indices(const unsigned int quadrature_point) const
ObserverPointer< const FEValuesBase< dim1, spacedim > > first_fe_values
double JxW(const unsigned int quadrature_point) const
ObserverPointer< const FEValuesBase< dim2, spacedim > > second_fe_values
std::unique_ptr< const FEValuesViews::RenumberingData > first_renumbering_data
std_cxx20::ranges::iota_view< unsigned int, unsigned int > coupling_dof_indices() const
unsigned int n_second_dofs() const
FEValuesExtractors::SecondCoupling< Extractor > get_second_extractor(const Extractor &extractor) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > second_dof_indices() const
DoFCouplingType dof_coupling_type
unsigned int n_first_dofs() const
void reinit(const FEValuesBase< dim1, spacedim > &fe_values_1, const FEValuesBase< dim2, spacedim > &fe_values_2, const DoFCouplingType &dof_coupling_type=DoFCouplingType::independent, const QuadratureCouplingType &quadrature_coupling_type=QuadratureCouplingType::tensor_product)
std::pair< unsigned int, unsigned int > coupling_dof_to_dof_indices(const unsigned int coupling_dof_index) const
std::unique_ptr< const FEValuesViews::RenumberingData > second_renumbering_data
Triangulation< dim, spacedim >::cell_iterator get_cell() const
const std::vector< Point< spacedim > > & get_quadrature_points() const
const unsigned int dofs_per_cell
const unsigned int n_quadrature_points
void get_function_values_from_local_dof_values(const InputVector &dof_values, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
void get_function_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< solution_gradient_type< typename InputVector::value_type > > &gradients) const
const InputVector & outer_to_inner_dofs(const InputVector &outer_vector) const
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
typename ViewType::value_type value_type
void get_function_gradients(const ReadVector< Number > &fe_function, std::vector< solution_gradient_type< Number > > &gradients) const
void inner_to_outer_values(const std::vector< ValueType > &inner_values, std::vector< ValueType > &outer_values) const
std::vector< ValueType > & outer_to_inner_values(std::vector< ValueType > &outer_values) const
typename ViewType::gradient_type gradient_type
void get_function_values(const ReadVector< Number > &fe_function, std::vector< solution_value_type< Number > > &values) const
RenumberedView(const ViewType &view, const RenumberingData &data)
typename ViewType::template solution_value_type< Number > solution_value_type
std::string get_unique_container_name(const std::string &prefix, const unsigned int size, const Number &exemplar_number) const
typename ViewType::template solution_gradient_type< Number > solution_gradient_type
value_type value(const unsigned int shape_function, const unsigned int q_point) const
Definition point.h:113
A class that provides a separate storage location on each thread that accesses the object.
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:35
constexpr bool running_in_debug_mode()
Definition config.h:73
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
QuadratureCouplingType
DoFCouplingType
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
std::vector< index_type > data
Definition mpi.cc:746
std::size_t size
Definition mpi.cc:745
typename ::internal::FEValuesViews::ViewType< dim, spacedim, Extractor >::type View
Tpetra::Vector< Number, LO, GO, NodeType< MemorySpace > > VectorType
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
std::string type_to_string(const T &t)
Definition utilities.h:1024
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition utilities.cc:466
constexpr unsigned int invalid_unsigned_int
Definition types.h:232
boost::integer_range< IncrementableType > iota_view
Definition iota_view.h:45
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
const std::vector< unsigned int > dof_renumbering
RenumberingData(const unsigned int n_inner_dofs=numbers::invalid_unsigned_int, const unsigned int n_inner_quadrature_points=numbers::invalid_unsigned_int, const std::vector< unsigned int > &dof_renumbering={}, const std::vector< unsigned int > &quadrature_renumbering={})
const unsigned int n_inner_quadrature_points
const std::vector< unsigned int > quadrature_renumbering
RenumberingData(const RenumberingData &other)=delete
Threads::ThreadLocalStorage< GeneralDataStorage > data_storage