Reference documentation for deal.II version 9.6.0
\(\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\}}\)
Loading...
Searching...
No Matches
function_restriction.cc
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2021 - 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
17
19namespace internal
20{
21 template <int dim>
24 const unsigned int component_in_dim_plus_1,
25 const double coordinate_value)
26 {
27 AssertIndexRange(component_in_dim_plus_1, dim + 1);
28
29 Point<dim + 1> output;
30 output[component_in_dim_plus_1] = coordinate_value;
31 for (int d = 0; d < dim; ++d)
32 {
33 const unsigned int component_to_write_to =
34 ::internal::coordinate_to_one_dim_higher<dim>(
35 component_in_dim_plus_1, d);
36 output[component_to_write_to] = point[d];
37 }
38
39 return output;
40 }
41} // namespace internal
42
43
44
45namespace Functions
46{
47 template <int dim>
49 const Function<dim + 1> &function,
50 const unsigned int direction,
51 const double coordinate_value)
52 : function(&function)
53 , restricted_direction(direction)
54 , coordinate_value(coordinate_value)
55 {
57 }
58
59
60
61 template <int dim>
62 double
64 const unsigned int component) const
65 {
66 const Point<dim + 1> full_point =
68 restricted_direction,
69 coordinate_value);
70
71 return function->value(full_point, component);
72 }
73
74
75
76 template <int dim>
79 const unsigned int component) const
80 {
81 const Point<dim + 1> full_point =
83 restricted_direction,
84 coordinate_value);
85
86 const Tensor<1, dim + 1> full_gradient =
87 function->gradient(full_point, component);
88
89 // CoordinateRestriction is constant in restricted direction. Through away
90 // the derivatives with respect to this direction and copy the other
91 // values.
92 Tensor<1, dim> grad;
93 for (unsigned int d = 0; d < dim; ++d)
94 {
95 const unsigned int index_to_write_from =
96 internal::coordinate_to_one_dim_higher<dim>(restricted_direction, d);
97 grad[d] = full_gradient[index_to_write_from];
98 }
99 return grad;
100 }
101
102
103
104 template <int dim>
107 const unsigned int component) const
108 {
109 const Point<dim + 1> full_point =
111 restricted_direction,
112 coordinate_value);
113
114 const Tensor<2, dim + 1> full_hessian =
115 function->hessian(full_point, component);
116
117 // CoordinateRestriction is constant in restricted direction. Through away
118 // the derivatives with respect to this direction and copy the other
119 // values.
121 for (unsigned int i = 0; i < dim; ++i)
122 {
123 const unsigned int i_to_write_from =
124 internal::coordinate_to_one_dim_higher<dim>(restricted_direction, i);
125 for (unsigned int j = 0; j < dim; ++j)
126 {
127 const unsigned int j_to_write_from =
128 internal::coordinate_to_one_dim_higher<dim>(restricted_direction,
129 j);
130 hess[i][j] = full_hessian[i_to_write_from][j_to_write_from];
131 }
133 return hess;
134 }
135
136
137
138 template <int dim>
140 const unsigned int open_direction,
141 const Point<dim> &point)
142 : function(&function)
143 , open_direction(open_direction)
144 , point(point)
145 {
147 }
148
149
150
151 template <int dim>
152 double
154 const unsigned int component) const
155 {
156 const Point<dim + 1> full_point =
157 internal::create_higher_dim_point(point, open_direction, point_1D[0]);
158 return function->value(full_point, component);
159 }
160
161
162
163 template <int dim>
166 const unsigned int component) const
167 {
168 const Point<dim + 1> full_point =
169 internal::create_higher_dim_point(point, open_direction, point_1D[0]);
170 const Tensor<1, dim + 1> full_gradient =
171 function->gradient(full_point, component);
172
173 // The PointRestrictions is constant in all but the open direction. Throw
174 // away the derivatives in all but this direction.
175 Tensor<1, 1> grad;
176 grad[0] = full_gradient[open_direction];
177 return grad;
178 }
179
180
181
182 template <int dim>
185 const unsigned int component) const
186 {
187 const Point<dim + 1> full_point =
188 internal::create_higher_dim_point(point, open_direction, point_1D[0]);
189 const Tensor<2, dim + 1> full_hessian =
190 function->hessian(full_point, component);
191
192 // The PointRestrictions is constant in all but the open direction. Throw
193 // away the derivatives in all but this direction.
195 hess[0][0] = full_hessian[open_direction][open_direction];
196 return hess;
197 }
198
199
200
201} // namespace Functions
202#include "function_restriction.inst"
double value(const Point< dim > &point, const unsigned int component) const override
Tensor< 1, dim > gradient(const Point< dim > &point, const unsigned int component) const override
SymmetricTensor< 2, dim > hessian(const Point< dim > &point, const unsigned int component) const override
CoordinateRestriction(const Function< dim+1 > &function, const unsigned int direction, const double coordinate_value)
PointRestriction(const Function< dim+1 > &function, const unsigned int open_direction, const Point< dim > &point)
double value(const Point< 1 > &point, const unsigned int component) const override
Tensor< 1, 1 > gradient(const Point< 1 > &point, const unsigned int component) const override
SymmetricTensor< 2, 1 > hessian(const Point< 1 > &point, const unsigned int component) const override
Definition point.h:111
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:503
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:504
#define AssertIndexRange(index, range)
Point< dim+1 > create_higher_dim_point(const Point< dim > &point, const unsigned int component_in_dim_plus_1, const double coordinate_value)