Reference documentation for deal.II version 9.5.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// Copyright (C) 2010 - 2021 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
18
20namespace internal
21{
22 template <int dim>
25 const unsigned int component_in_dim_plus_1,
26 const double coordinate_value)
27 {
28 AssertIndexRange(component_in_dim_plus_1, dim + 1);
29
30 Point<dim + 1> output;
31 output(component_in_dim_plus_1) = coordinate_value;
32 for (int d = 0; d < dim; ++d)
33 {
34 const unsigned int component_to_write_to =
35 ::internal::coordinate_to_one_dim_higher<dim>(
36 component_in_dim_plus_1, d);
37 output(component_to_write_to) = point(d);
38 }
39
40 return output;
41 }
42} // namespace internal
43
44
45
46namespace Functions
47{
48 template <int dim>
50 const Function<dim + 1> &function,
51 const unsigned int direction,
52 const double coordinate_value)
53 : function(&function)
54 , restricted_direction(direction)
55 , coordinate_value(coordinate_value)
56 {
58 }
59
60
61
62 template <int dim>
63 double
65 const unsigned int component) const
66 {
67 const Point<dim + 1> full_point =
69 restricted_direction,
70 coordinate_value);
71
72 return function->value(full_point, component);
73 }
74
75
76
77 template <int dim>
80 const unsigned int component) const
81 {
82 const Point<dim + 1> full_point =
84 restricted_direction,
85 coordinate_value);
86
87 const Tensor<1, dim + 1> full_gradient =
88 function->gradient(full_point, component);
89
90 // CoordinateRestriction is constant in restricted direction. Through away
91 // the derivatives with respect to this direction and copy the other
92 // values.
93 Tensor<1, dim> grad;
94 for (unsigned int d = 0; d < dim; ++d)
95 {
96 const unsigned int index_to_write_from =
97 internal::coordinate_to_one_dim_higher<dim>(restricted_direction, d);
98 grad[d] = full_gradient[index_to_write_from];
99 }
100 return grad;
101 }
102
103
104
105 template <int dim>
108 const unsigned int component) const
109 {
110 const Point<dim + 1> full_point =
112 restricted_direction,
113 coordinate_value);
114
115 const Tensor<2, dim + 1> full_hessian =
116 function->hessian(full_point, component);
117
118 // CoordinateRestriction is constant in restricted direction. Through away
119 // the derivatives with respect to this direction and copy the other
120 // values.
122 for (unsigned int i = 0; i < dim; ++i)
123 {
124 const unsigned int i_to_write_from =
125 internal::coordinate_to_one_dim_higher<dim>(restricted_direction, i);
126 for (unsigned int j = 0; j < dim; ++j)
127 {
128 const unsigned int j_to_write_from =
129 internal::coordinate_to_one_dim_higher<dim>(restricted_direction,
130 j);
131 hess[i][j] = full_hessian[i_to_write_from][j_to_write_from];
132 }
134 return hess;
135 }
136
137
138
139 template <int dim>
141 const unsigned int open_direction,
142 const Point<dim> & point)
143 : function(&function)
144 , open_direction(open_direction)
145 , point(point)
146 {
148 }
149
150
151
152 template <int dim>
153 double
155 const unsigned int component) const
156 {
157 const Point<dim + 1> full_point =
158 internal::create_higher_dim_point(point, open_direction, point_1D(0));
159 return function->value(full_point, component);
160 }
161
162
163
164 template <int dim>
167 const unsigned int component) const
168 {
169 const Point<dim + 1> full_point =
170 internal::create_higher_dim_point(point, open_direction, point_1D(0));
171 const Tensor<1, dim + 1> full_gradient =
172 function->gradient(full_point, component);
173
174 // The PointRestrictions is constant in all but the open direction. Throw
175 // away the derivatives in all but this direction.
176 Tensor<1, 1> grad;
177 grad[0] = full_gradient[open_direction];
178 return grad;
179 }
180
181
182
183 template <int dim>
186 const unsigned int component) const
187 {
188 const Point<dim + 1> full_point =
189 internal::create_higher_dim_point(point, open_direction, point_1D(0));
190 const Tensor<2, dim + 1> full_hessian =
191 function->hessian(full_point, component);
192
193 // The PointRestrictions is constant in all but the open direction. Throw
194 // away the derivatives in all but this direction.
196 hess[0][0] = full_hessian[open_direction][open_direction];
197 return hess;
198 }
199
200
201
202} // namespace Functions
203#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:112
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#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)