Reference documentation for deal.II version GIT 4efb66ecd0 2023-02-07 13:45:02+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\}}\)
function_tools.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 
17 
18 #include <cmath>
19 
20 
22 namespace FunctionTools
23 {
24  template <int dim>
25  void
27  const Function<dim> & function,
28  const BoundingBox<dim> & box,
29  std::pair<double, double> & value_bounds,
30  std::array<std::pair<double, double>, dim> &gradient_bounds,
31  const unsigned int component)
32  {
33  const Point<dim> center = box.center();
34  const double value = function.value(center, component);
35  const Tensor<1, dim> gradient = function.gradient(center, component);
36  const SymmetricTensor<2, dim> hessian = function.hessian(center, component);
37 
38  // Deviation from function value at the center, based on the
39  // Taylor-expansion: |f'| * dx + 1/2 * |f''| * dx^2, (in 1D). dx is half
40  // the side-length of the box.
41  double taylor_bound_f = 0;
42 
43  for (unsigned int i = 0; i < dim; ++i)
44  {
45  const double dx_i = .5 * box.side_length(i);
46 
47  taylor_bound_f += std::abs(gradient[i]) * dx_i;
48 
49  // Deviation from value of df/dx_i at the center,
50  // |f''| * dx, (in 1D).
51  double taylor_bound_dfdxi = 0;
52 
53  for (unsigned int j = 0; j < dim; ++j)
54  {
55  const double dx_j = .5 * box.side_length(j);
56 
57  taylor_bound_dfdxi += std::abs(hessian[i][j]) * dx_j;
58  taylor_bound_f += .5 * std::abs(hessian[i][j]) * dx_i * dx_j;
59  }
60 
61  gradient_bounds[i].first = gradient[i] - taylor_bound_dfdxi;
62  gradient_bounds[i].second = gradient[i] + taylor_bound_dfdxi;
63  }
64 
65  value_bounds.first = value - taylor_bound_f;
66  value_bounds.second = value + taylor_bound_f;
67  }
68 
69 } // namespace FunctionTools
70 
71 #include "function_tools.inst"
72 
Point< spacedim, Number > center() const
Number side_length(const unsigned int direction) const
Definition: point.h:111
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:461
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:462
Point< 3 > center
void taylor_estimate_function_bounds(const Function< dim > &function, const BoundingBox< dim > &box, std::pair< double, double > &value_bounds, std::array< std::pair< double, double >, dim > &gradient_bounds, const unsigned int component=0)