5#ifndef DUNE_LOCALFUNCTIONS_CROUZEIXRAVIART_HH
6#define DUNE_LOCALFUNCTIONS_CROUZEIXRAVIART_HH
11#include <dune/common/fmatrix.hh>
12#include <dune/common/fvector.hh>
14#include <dune/geometry/type.hh>
15#include <dune/geometry/referenceelements.hh>
21namespace Dune {
namespace Impl
29 template<
class D,
class R,
unsigned int dim>
30 class CrouzeixRaviartLocalBasis
33 using Traits = LocalBasisTraits<D,dim,FieldVector<D,dim>,R,1,FieldVector<R,1>,FieldMatrix<R,1,dim> >;
37 static constexpr unsigned int size ()
44 std::vector<typename Traits::RangeType>& out)
const
48 std::fill(out.begin(), out.end()-1, 1.0);
51 for (
unsigned int i=0; i<dim; i++)
53 out[i] -= dim * x[dim-i-1];
54 out.back() += dim*x[i];
64 std::vector<typename Traits::JacobianType>& out)
const
68 for (
unsigned i=0; i<dim; i++)
69 for (
unsigned j=0; j<dim; j++)
70 out[i][0][j] = (i==(dim-1-j)) ? -(double)dim : 0;
72 std::fill(out.back()[0].begin(), out.back()[0].end(), dim);
81 void partial(
const std::array<unsigned int,dim>& order,
83 std::vector<typename Traits::RangeType>& out)
const
85 auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
89 if (totalOrder == 0) {
90 evaluateFunction(in, out);
96 auto direction = std::find(order.begin(), order.end(), 1)-order.begin();
98 for (
unsigned int i=0; i<dim; i++)
99 out[i] = (i==(dim-1-direction)) ? -(double)dim : 0.0;
104 std::fill(out.begin(), out.end(), 0);
108 static constexpr unsigned int order ()
118 template<
unsigned int dim>
119 class CrouzeixRaviartLocalCoefficients
123 CrouzeixRaviartLocalCoefficients ()
126 for (std::size_t i=0; i<size(); i++)
127 localKeys_[i] = LocalKey(i,1,0);
131 static constexpr std::size_t size ()
137 const LocalKey& localKey (std::size_t i)
const
139 return localKeys_[i];
143 std::vector<LocalKey> localKeys_;
150 template<
class LocalBasis>
151 class CrouzeixRaviartLocalInterpolation
162 template<
typename F,
typename C>
163 void interpolate (
const F& f, std::vector<C>& out)
const
165 constexpr auto dim = LocalBasis::Traits::dimDomain;
167 out.resize(LocalBasis::size());
170 auto refSimplex = ReferenceElements<typename LocalBasis::Traits::DomainFieldType,dim>::simplex();
171 for (
int i=0; i<refSimplex.size(1); i++)
172 out[i] = f(refSimplex.template geometry<1>(i).center());
186 template<
class D,
class R,
int dim>
193 Impl::CrouzeixRaviartLocalCoefficients<dim>,
194 Impl::CrouzeixRaviartLocalInterpolation<Impl::CrouzeixRaviartLocalBasis<D,R,dim> > >;
207 return coefficients_;
214 return interpolation_;
218 static constexpr std::size_t
size()
225 static constexpr GeometryType
type()
227 return GeometryTypes::simplex(dim);
231 Impl::CrouzeixRaviartLocalBasis<D,R,dim> basis_;
232 Impl::CrouzeixRaviartLocalCoefficients<dim> coefficients_;
233 Impl::CrouzeixRaviartLocalInterpolation<Impl::CrouzeixRaviartLocalBasis<D,R,dim> > interpolation_;
Definition bdfmcube.hh:18
D DomainType
domain type
Definition common/localbasis.hh:43
traits helper struct
Definition localfiniteelementtraits.hh:13
LB LocalBasisType
Definition localfiniteelementtraits.hh:16
LC LocalCoefficientsType
Definition localfiniteelementtraits.hh:20
LI LocalInterpolationType
Definition localfiniteelementtraits.hh:24
Crouzeix-Raviart finite element.
Definition crouzeixraviart.hh:188
static constexpr GeometryType type()
The reference element that the local finite element is defined on.
Definition crouzeixraviart.hh:225
const Traits::LocalCoefficientsType & localCoefficients() const
Returns the assignment of the degrees of freedom to the element subentities.
Definition crouzeixraviart.hh:205
const Traits::LocalInterpolationType & localInterpolation() const
Returns object that evaluates degrees of freedom.
Definition crouzeixraviart.hh:212
const Traits::LocalBasisType & localBasis() const
Returns the local basis, i.e., the set of shape functions.
Definition crouzeixraviart.hh:198
static constexpr std::size_t size()
The number of shape functions.
Definition crouzeixraviart.hh:218