DUNE-FEM (unstable)

localinterpolation.hh
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3// SPDX-FileCopyrightInfo: Copyright © DUNE Project contributors, see file LICENSE.md in module root
4// SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
5#ifndef DUNE_LOCALFUNCTIONS_ENRICHED_SIMPLEXP1BUBBLE_LOCALINTERPOLATION_HH
6#define DUNE_LOCALFUNCTIONS_ENRICHED_SIMPLEXP1BUBBLE_LOCALINTERPOLATION_HH
7
8#include <type_traits>
9#include <vector>
10
11namespace Dune
12{
32 template<class LB>
34 {
35 static const int dim = LB::dimension;
36 static const int numVertices = dim+1;
37
38 using DomainType = typename LB::Traits::DomainType;
39 using RangeType = typename LB::Traits::RangeType;
40
41 public:
49 template<class F, class C,
50 class R = std::invoke_result_t<F, DomainType>,
51 std::enable_if_t<std::is_convertible_v<R, C>, int> = 0>
52 static constexpr void interpolate (const F& f, std::vector<C>& out)
53 {
54 out.resize(numVertices+1);
55
56 // vertices
57 DomainType x(0);
58 out[0] = f(x);
59
60 for (int i = 0; i < dim; ++i) {
61 x = 0;
62 x[i] = 1;
63 out[i+1] = f(x);
64 }
65
66 // element bubble
67 x = 1.0/(dim+1);
68 R y = f(x);
69
70 // evaluate the other shape functions in x and subtract this value
71 std::vector<RangeType> sfValues;
72 LB::evaluateFunction(x, sfValues);
73
74 out[numVertices] = y;
75 for (int i = 0; i < numVertices; ++i)
76 out[numVertices] -= out[i]*sfValues[i];
77 }
78 };
79
80} // end namespace Dune
81
82#endif // DUNE_LOCALFUNCTIONS_ENRICHED_SIMPLEXP1BUBBLE_LOCALINTERPOLATION_HH
Interpolation into the SimplexP1BubbleLocalBasis.
Definition: localinterpolation.hh:34
static constexpr void interpolate(const F &f, std::vector< C > &out)
Local interpolation of the function f.
Definition: localinterpolation.hh:52
Dune namespace.
Definition: alignedallocator.hh:13
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 12, 23:30, 2024)