3 #ifndef DUNE_LOCALFUNCTIONS_CROUZEIXRAVIART_HH
4 #define DUNE_LOCALFUNCTIONS_CROUZEIXRAVIART_HH
9 #include <dune/common/fmatrix.hh>
10 #include <dune/common/fvector.hh>
12 #include <dune/geometry/type.hh>
13 #include <dune/geometry/referenceelements.hh>
20 namespace Dune {
namespace Impl
28 template<
class D,
class R,
unsigned int dim>
29 class CrouzeixRaviartLocalBasis
32 using Traits = LocalBasisTraits<D,dim,FieldVector<D,dim>,R,1,FieldVector<R,1>,FieldMatrix<R,1,dim> >;
36 static constexpr
unsigned int size ()
43 std::vector<typename Traits::RangeType>& out)
const
47 std::fill(out.begin(), out.end()-1, 1.0);
50 for (
unsigned int i=0; i<dim; i++)
52 out[i] -= dim * x[dim-i-1];
53 out.back() += dim*x[i];
63 std::vector<typename Traits::JacobianType>& out)
const
67 for (
unsigned i=0; i<dim; i++)
68 for (
unsigned j=0; j<dim; j++)
69 out[i][0][j] = (i==(dim-1-j)) ? -(double)dim : 0;
71 std::fill(out.back()[0].begin(), out.back()[0].end(), dim);
80 void partial(
const std::array<unsigned int,dim>& order,
82 std::vector<typename Traits::RangeType>& out)
const
84 auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
88 if (totalOrder == 0) {
89 evaluateFunction(in, out);
95 auto direction = std::find(order.begin(), order.end(), 1)-order.begin();
97 for (
unsigned int i=0; i<dim; i++)
98 out[i] = (i==(dim-1-direction)) ? -(double)dim : 0.0;
103 std::fill(out.begin(), out.end(), 0);
107 static constexpr
unsigned int order ()
117 template<
unsigned int dim>
118 class CrouzeixRaviartLocalCoefficients
122 CrouzeixRaviartLocalCoefficients ()
125 for (std::size_t i=0; i<size(); i++)
126 localKeys_[i] = LocalKey(i,dim-1,0);
130 static constexpr std::size_t size ()
136 const LocalKey& localKey (std::size_t i)
const
138 return localKeys_[i];
142 std::vector<LocalKey> localKeys_;
149 template<
class LocalBasis>
150 class CrouzeixRaviartLocalInterpolation
161 template<
typename F,
typename C>
162 void interpolate (
const F& ff, std::vector<C>& out)
const
164 constexpr
auto dim = LocalBasis::Traits::dimDomain;
166 auto&& f = Impl::makeFunctionWithCallOperator<typename LocalBasis::Traits::DomainType>(ff);
168 out.resize(LocalBasis::size());
171 auto refSimplex = ReferenceElements<typename LocalBasis::Traits::DomainFieldType,dim>::simplex();
172 for (
int i=0; i<refSimplex.size(1); i++)
173 out[i] = f(refSimplex.template geometry<1>(i).center());
187 template<
class D,
class R,
int dim>
194 Impl::CrouzeixRaviartLocalCoefficients<dim>,
195 Impl::CrouzeixRaviartLocalInterpolation<Impl::CrouzeixRaviartLocalBasis<D,R,dim> > >;
208 return coefficients_;
215 return interpolation_;
219 static constexpr std::size_t
size()
226 static constexpr GeometryType
type()
228 return GeometryTypes::simplex(dim);
232 Impl::CrouzeixRaviartLocalBasis<D,R,dim> basis_;
233 Impl::CrouzeixRaviartLocalCoefficients<dim> coefficients_;
234 Impl::CrouzeixRaviartLocalInterpolation<Impl::CrouzeixRaviartLocalBasis<D,R,dim> > interpolation_;
239 #endif // DUNE_LOCALFUNCTIONS_CROUZEIXRAVIART_HH