dune-localfunctions 2.10
Loading...
Searching...
No Matches
raviartthomas0prismlocalbasis.hh
Go to the documentation of this file.
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_RAVIARTTHOMAS0_PRISM_LOCALBASIS_HH
6#define DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS0_PRISM_LOCALBASIS_HH
7
8#include <numeric>
9#include <vector>
10
11#include <dune/common/fmatrix.hh>
12
13#include "../../common/localbasis.hh"
14
15namespace Dune
16{
26 template<class D, class R>
28 {
29
30 public:
31 typedef LocalBasisTraits<D,3,Dune::FieldVector<D,3>,R,3,Dune::FieldVector<R,3>,
32 Dune::FieldMatrix<R,3,3> > Traits;
33
39 RT0PrismLocalBasis (std::bitset<5> s = 0)
40 {
41 for (size_t i=0; i<size(); i++)
42 sign[i] = s[i] ? -1.0 : 1.0;
43 }
44
46 unsigned int size () const
47 {
48 return 5;
49 }
50
57 inline void evaluateFunction (const typename Traits::DomainType& in,
58 std::vector<typename Traits::RangeType>& out) const
59 {
60 out.resize(5);
61
62 out[0] = { in[0], -1.0 + in[1], 0.0};
63
64 out[1] = { -1.0 + in[0], in[1], 0.0};
65
66 out[2] = { in[0], in[1], 0.0};
67
68 out[3] = { 0.0, 0.0, -2.0 + 2.0*in[2]};
69
70 out[4] = { 0.0, 0.0, 2.0*in[2]};
71
72 for (std::size_t i=0; i<out.size(); i++)
73 out[i] *= sign[i];
74
75 }
76
83 inline void evaluateJacobian (const typename Traits::DomainType& in,
84 std::vector<typename Traits::JacobianType>& out) const
85 {
86 out.resize(5);
87
88 for(int i=0; i<size(); i++)
89 for(int j=0; j<3; j++)
90 out[i][j] = {0.0, 0.0, 0.0};
91
92 out[0][0][0] = sign[0];
93 out[0][1][1] = sign[0];
94
95 out[1][0][0] = sign[1];
96 out[1][1][1] = sign[1];
97
98 out[2][0][0] = sign[2];
99 out[2][1][1] = sign[2];
100
101 out[3][2][2] = sign[3]*(2.0);
102
103 out[4][2][2] = sign[4]*(2.0);
104 }
105
107 void partial (const std::array<unsigned int, 3>& order,
108 const typename Traits::DomainType& in, // position
109 std::vector<typename Traits::RangeType>& out) const // return value
110 {
111 auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
112 if (totalOrder == 0) {
113 evaluateFunction(in, out);
114 } else {
115 DUNE_THROW(NotImplemented, "Desired derivative order is not implemented");
116 }
117 }
118
120 unsigned int order () const
121 {
122 return 1;
123 }
124
125 private:
126 std::array<R,5> sign;
127 };
128}
129#endif // DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS0_PRISM_LOCALBASIS_HH
Definition bdfmcube.hh:18
Type traits for LocalBasisVirtualInterface.
Definition common/localbasis.hh:35
D DomainType
domain type
Definition common/localbasis.hh:43
First order Raviart-Thomas shape functions on the reference prism.
Definition raviartthomas0prismlocalbasis.hh:28
RT0PrismLocalBasis(std::bitset< 5 > s=0)
Make set number s, where 0 <= s < 32.
Definition raviartthomas0prismlocalbasis.hh:39
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition raviartthomas0prismlocalbasis.hh:57
unsigned int size() const
number of shape functions
Definition raviartthomas0prismlocalbasis.hh:46
unsigned int order() const
Polynomial order of the shape functions.
Definition raviartthomas0prismlocalbasis.hh:120
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition raviartthomas0prismlocalbasis.hh:83
void partial(const std::array< unsigned int, 3 > &order, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate partial derivatives of all shape functions.
Definition raviartthomas0prismlocalbasis.hh:107
LocalBasisTraits< D, 3, Dune::FieldVector< D, 3 >, R, 3, Dune::FieldVector< R, 3 >, Dune::FieldMatrix< R, 3, 3 > > Traits
Definition raviartthomas0prismlocalbasis.hh:32