15#ifndef dealii_fe_coupling_values_h
16#define dealii_fe_coupling_values_h
34template <
int dim,
int spacedim>
168 template <
typename ViewType>
186 template <
typename Number>
194 template <
typename Number>
230 value(
const unsigned int shape_function,
const unsigned int q_point)
const;
250 const unsigned int q_point)
const;
266 template <
typename Number>
291 template <
class InputVector>
294 const InputVector &dof_values,
313 template <
typename Number>
339 template <
class InputVector>
342 const InputVector &dof_values,
370 template <
typename Number>
373 const unsigned int size,
374 const Number &exemplar_number)
const;
380 template <
typename InputVector>
388 template <
typename ValueType>
389 std::vector<ValueType> &
397 template <
typename ValueType>
400 std::vector<ValueType> &outer_values)
const;
730template <
int dim1,
int dim2 = dim1,
int spacedim = dim1>
829 template <
typename Extractor>
837 template <
typename Extractor>
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,
952 std::pair<unsigned int, unsigned int>
960 std::pair<unsigned int, unsigned int>
974 template <
typename Extractor>
985 template <
typename Extractor>
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)
1094 for (
const auto i : dof_renumbering)
1098 for (
const auto q : quadrature_renumbering)
1105 template <
typename ViewType>
1107 const RenumberingData &
data)
1114 template <
typename ViewType>
1115 template <
typename Number>
1117 RenumberedView<ViewType>::get_unique_container_name(
1118 const std::string &prefix,
1119 const unsigned int size,
1120 const Number &exemplar_number)
const
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
1136 const auto inner_shape_function =
data.dof_renumbering.empty() ?
1138 data.dof_renumbering[shape_function];
1139 const auto inner_q_point =
data.quadrature_renumbering.empty() ?
1141 data.quadrature_renumbering[q_point];
1143 return value_type(0);
1148 return view.value(inner_shape_function, inner_q_point);
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
1162 const auto inner_shape_function =
data.dof_renumbering.empty() ?
1164 data.dof_renumbering[shape_function];
1165 const auto inner_q_point =
data.quadrature_renumbering.empty() ?
1167 data.quadrature_renumbering[q_point];
1169 return gradient_type();
1171 return view.gradient(inner_shape_function, inner_q_point);
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
1183 if (
data.quadrature_renumbering.empty())
1185 return outer_values;
1190 get_unique_container_name(
"RenumberedView::outer_to_inner_values",
1191 data.n_inner_quadrature_points,
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;
1203 template <
typename ViewType>
1204 template <
typename VectorType>
1206 RenumberedView<ViewType>::outer_to_inner_dofs(
1207 const VectorType &outer_dofs)
const
1210 if (
data.dof_renumbering.empty())
1217 get_unique_container_name(
"RenumberedView::outer_to_inner_dofs",
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)
1226 const auto inner_i =
data.dof_renumbering[i];
1230 inner_dofs[inner_i] = outer_dofs[i];
1239 template <
typename ViewType>
1240 template <
typename ValueType>
1242 RenumberedView<ViewType>::inner_to_outer_values(
1243 const std::vector<ValueType> &inner_values,
1244 std::vector<ValueType> &outer_values)
const
1248 if (
data.quadrature_renumbering.empty())
1253 for (
unsigned int i = 0; i <
data.quadrature_renumbering.size(); ++i)
1255 outer_values[i] = inner_values[
data.quadrature_renumbering[i]];
1261 template <
typename ViewType>
1262 template <
typename Number>
1264 RenumberedView<ViewType>::get_function_values(
1266 std::vector<solution_value_type<Number>> &values)
const
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);
1275 template <
typename ViewType>
1276 template <
typename InputVector>
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)
1283 const auto &inner_dof_values = outer_to_inner_dofs(dof_values);
1284 auto &inner_values = outer_to_inner_values(values);
1286 view.get_function_values_from_local_dof_values(inner_dof_values,
1288 inner_to_outer_values(inner_values, values);
1292 template <
typename ViewType>
1293 template <
typename Number>
1295 RenumberedView<ViewType>::get_function_gradients(
1297 std::vector<solution_gradient_type<Number>> &gradients)
const
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);
1306 template <
typename ViewType>
1307 template <
typename InputVector>
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>>
1314 const auto &inner_dof_values = outer_to_inner_dofs(dof_values);
1315 auto &inner_gradients = outer_to_inner_values(gradients);
1317 view.get_function_gradients_from_local_dof_values(inner_dof_values,
1319 inner_to_outer_values(inner_gradients, gradients);
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,
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);
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;
1352template <
int dim1,
int dim2,
int spacedim>
1354 : first_fe_values(nullptr)
1355 , second_fe_values(nullptr)
1357 , n_quadrature_points_(0)
1363template <
int dim1,
int dim2,
int spacedim>
1376template <
int dim1,
int dim2,
int spacedim>
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;
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;
1393 unsigned int first_n_inner_quadrature_points =
1395 unsigned int second_n_inner_quadrature_points =
1399 std::vector<unsigned int> first_dofs_map;
1400 std::vector<unsigned int> second_dofs_map;
1402 std::vector<unsigned int> first_quad_map;
1403 std::vector<unsigned int> second_quad_map;
1407 fe_values_2.
get_cell()->diameter()) *
1411 switch (dof_coupling_type)
1416 first_dofs_map.clear();
1417 second_dofs_map.clear();
1425 first_dofs_map.resize(n_coupling_dofs_);
1426 second_dofs_map.resize(n_coupling_dofs_);
1428 unsigned int idx = 0;
1429 for (
const unsigned int &i : fe_values_1.dof_indices())
1431 first_dofs_map[idx] = i;
1434 for (
const unsigned int &i : fe_values_2.dof_indices())
1437 second_dofs_map[idx++] = i;
1447 switch (quadrature_coupling_type)
1454 n_quadrature_points_ =
1455 quadrature_points_1.size() * quadrature_points_2.size();
1457 first_quad_map.resize(n_quadrature_points_);
1458 second_quad_map.resize(n_quadrature_points_);
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())
1464 first_quad_map[idx] = i;
1465 second_quad_map[idx] = j;
1474 ExcMessage(
"The two FEValuesBase objects must have the same "
1475 "number of quadrature points"));
1479 first_quad_map.clear();
1480 second_quad_map.clear();
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"));
1493 for (
const unsigned int &i : fe_values_1.quadrature_point_indices())
1495 Assert(quadrature_points_1[i].distance(quadrature_points_2[i]) <
1498 "The two FEValuesBase objects must have the same "
1499 "quadrature points"));
1504 first_quad_map.clear();
1505 second_quad_map.clear();
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"));
1521 first_quad_map.clear();
1525 for (
const unsigned int &i : fe_values_1.quadrature_point_indices())
1528 for (
const unsigned int &j :
1529 fe_values_2.quadrature_point_indices())
1530 if (quadrature_points_1[i].distance(quadrature_points_2[j]) <
1534 second_quad_map[i] = j;
1539 "The two FEValuesBase objects must have the same "
1540 "quadrature points, even if not in the same order."));
1550 for (
const unsigned int &i : fe_values_1.quadrature_point_indices())
1552 for (
const unsigned int &j :
1553 fe_values_2.quadrature_point_indices())
1554 if (quadrature_points_1[i].distance(quadrature_points_2[j]) <
1557 first_quad_map.emplace_back(i);
1558 second_quad_map.emplace_back(j);
1562 n_quadrature_points_ = first_quad_map.size();
1570 first_renumbering_data = std::make_unique<FEValuesViews::RenumberingData>(
1572 first_n_inner_quadrature_points,
1576 second_renumbering_data = std::make_unique<FEValuesViews::RenumberingData>(
1577 second_n_inner_dofs,
1578 second_n_inner_quadrature_points,
1584template <
int dim1,
int dim2,
int spacedim>
1590 const auto first_q = first_renumbering_data->quadrature_renumbering.empty() ?
1592 first_renumbering_data->quadrature_renumbering[q];
1597 const auto second_q =
1598 second_renumbering_data->quadrature_renumbering.empty() ?
1600 second_renumbering_data->quadrature_renumbering[q];
1601 return first_fe_values->JxW(first_q) * second_fe_values->JxW(second_q);
1604 return first_fe_values->JxW(first_q);
1609template <
int dim1,
int dim2,
int spacedim>
1614 0U, n_quadrature_points_);
1619template <
int dim1,
int dim2,
int spacedim>
1625 "Dofs are independent. You cannot ask for coupling dofs."));
1627 0U, n_coupling_dofs_);
1632template <
int dim1,
int dim2,
int spacedim>
1637 0U, n_first_dofs());
1642template <
int dim1,
int dim2,
int spacedim>
1647 0U, n_second_dofs());
1652template <
int dim1,
int dim2,
int spacedim>
1658 "Dofs are independent. You cannot ask for coupling dofs."));
1659 return n_coupling_dofs_;
1664template <
int dim1,
int dim2,
int spacedim>
1668 return first_renumbering_data->n_dofs;
1673template <
int dim1,
int dim2,
int spacedim>
1677 return second_renumbering_data->n_dofs;
1682template <
int dim1,
int dim2,
int spacedim>
1686 return n_quadrature_points_;
1691template <
int dim1,
int dim2,
int spacedim>
1694 const unsigned int quadrature_point)
const
1697 const auto first_q = first_fe_values->quadrature_point(
1698 first_renumbering_data->quadrature_renumbering.empty() ?
1700 first_renumbering_data->quadrature_renumbering[quadrature_point]);
1702 const auto second_q = second_fe_values->quadrature_point(
1703 second_renumbering_data->quadrature_renumbering.empty() ?
1705 second_renumbering_data->quadrature_renumbering[quadrature_point]);
1707 return {first_q, second_q};
1712template <
int dim1,
int dim2,
int spacedim>
1713std::pair<unsigned int, unsigned int>
1716 const unsigned int quadrature_point)
const
1719 const auto first_id =
1720 first_renumbering_data->quadrature_renumbering.empty() ?
1722 first_renumbering_data->quadrature_renumbering[quadrature_point];
1724 const auto second_id =
1725 second_renumbering_data->quadrature_renumbering.empty() ?
1727 second_renumbering_data->quadrature_renumbering[quadrature_point];
1728 return std::make_pair(first_id, second_id);
1733template <
int dim1,
int dim2,
int spacedim>
1734std::pair<unsigned int, unsigned int>
1736 const unsigned int coupling_dof_index)
const
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];
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);
1753template <
int dim1,
int dim2,
int spacedim>
1754template <
typename Extractor>
1757 const Extractor &extractor)
const
1764template <
int dim1,
int dim2,
int spacedim>
1765template <
typename Extractor>
1768 const Extractor &extractor)
const
1775template <
int dim1,
int dim2,
int spacedim>
1776template <
typename Extractor>
1784 (*first_fe_values)[extractor.
extractor], *first_renumbering_data);
1789template <
int dim1,
int dim2,
int spacedim>
1790template <
typename Extractor>
1798 (*second_fe_values)[extractor.
extractor], *second_renumbering_data);
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
const RenumberingData & data
A class that provides a separate storage location on each thread that accesses the object.
#define DEAL_II_NAMESPACE_OPEN
constexpr bool running_in_debug_mode()
#define DEAL_II_NAMESPACE_CLOSE
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
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)
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
constexpr unsigned int invalid_unsigned_int
boost::integer_range< IncrementableType > iota_view
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
const unsigned int n_inner_dofs
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
const unsigned int n_dofs
const unsigned int n_quadrature_points