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>
22namespace Dune {
namespace Impl
30 template<
class D,
class R,
unsigned int dim>
31 class CrouzeixRaviartLocalBasis
34 using Traits = LocalBasisTraits<D,dim,FieldVector<D,dim>,R,1,FieldVector<R,1>,FieldMatrix<R,1,dim> >;
38 static constexpr unsigned int size ()
45 std::vector<typename Traits::RangeType>& out)
const
49 std::fill(out.begin(), out.end()-1, 1.0);
52 for (
unsigned int i=0; i<dim; i++)
54 out[i] -= dim * x[dim-i-1];
55 out.back() += dim*x[i];
65 std::vector<typename Traits::JacobianType>& out)
const
69 for (
unsigned i=0; i<dim; i++)
70 for (
unsigned j=0; j<dim; j++)
71 out[i][0][j] = (i==(dim-1-j)) ? -(double)dim : 0;
73 std::fill(out.back()[0].begin(), out.back()[0].end(), dim);
82 void partial(
const std::array<unsigned int,dim>& order,
84 std::vector<typename Traits::RangeType>& out)
const
86 auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
90 if (totalOrder == 0) {
91 evaluateFunction(in, out);
97 auto direction = std::find(order.begin(), order.end(), 1)-order.begin();
99 for (
unsigned int i=0; i<dim; i++)
100 out[i] = (i==(dim-1-direction)) ? -(double)dim : 0.0;
105 std::fill(out.begin(), out.end(), 0);
109 static constexpr unsigned int order ()
119 template<
unsigned int dim>
120 class CrouzeixRaviartLocalCoefficients
124 CrouzeixRaviartLocalCoefficients ()
127 for (std::size_t i=0; i<size(); i++)
128 localKeys_[i] = LocalKey(i,dim-1,0);
132 static constexpr std::size_t size ()
138 const LocalKey& localKey (std::size_t i)
const
140 return localKeys_[i];
144 std::vector<LocalKey> localKeys_;
151 template<
class LocalBasis>
152 class CrouzeixRaviartLocalInterpolation
163 template<
typename F,
typename C>
164 void interpolate (
const F& ff, std::vector<C>& out)
const
166 constexpr auto dim = LocalBasis::Traits::dimDomain;
168 auto&& f = Impl::makeFunctionWithCallOperator<typename LocalBasis::Traits::DomainType>(ff);
170 out.resize(LocalBasis::size());
173 auto refSimplex = ReferenceElements<typename LocalBasis::Traits::DomainFieldType,dim>::simplex();
174 for (
int i=0; i<refSimplex.size(1); i++)
175 out[i] = f(refSimplex.template geometry<1>(i).center());
189 template<
class D,
class R,
int dim>
196 Impl::CrouzeixRaviartLocalCoefficients<dim>,
197 Impl::CrouzeixRaviartLocalInterpolation<Impl::CrouzeixRaviartLocalBasis<D,R,dim> > >;
210 return coefficients_;
217 return interpolation_;
221 static constexpr std::size_t
size()
228 static constexpr GeometryType
type()
230 return GeometryTypes::simplex(dim);
234 Impl::CrouzeixRaviartLocalBasis<D,R,dim> basis_;
235 Impl::CrouzeixRaviartLocalCoefficients<dim> coefficients_;
236 Impl::CrouzeixRaviartLocalInterpolation<Impl::CrouzeixRaviartLocalBasis<D,R,dim> > interpolation_;
Definition bdfmcube.hh:18
D DomainType
domain type
Definition common/localbasis.hh:42
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:191
static constexpr GeometryType type()
The reference element that the local finite element is defined on.
Definition crouzeixraviart.hh:228
const Traits::LocalCoefficientsType & localCoefficients() const
Returns the assignment of the degrees of freedom to the element subentities.
Definition crouzeixraviart.hh:208
const Traits::LocalInterpolationType & localInterpolation() const
Returns object that evaluates degrees of freedom.
Definition crouzeixraviart.hh:215
const Traits::LocalBasisType & localBasis() const
Returns the local basis, i.e., the set of shape functions.
Definition crouzeixraviart.hh:201
static constexpr std::size_t size()
The number of shape functions.
Definition crouzeixraviart.hh:221