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_tools.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2010 - 2023 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
22namespace 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
virtual SymmetricTensor< 2, dim, RangeNumberType > hessian(const Point< dim > &p, const unsigned int component=0) const
virtual Tensor< 1, dim, RangeNumberType > gradient(const Point< dim > &p, const unsigned int component=0) const
virtual RangeNumberType value(const Point< dim > &p, const unsigned int component=0) const
Definition point.h:112
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
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)
::VectorizedArray< Number, width > abs(const ::VectorizedArray< Number, width > &)