3#ifndef DUNE_LOCALFUNCTIONS_NEDELEC_NEDELECSIMPLEX_NEDELECSIMPLEXPREBASIS_HH
4#define DUNE_LOCALFUNCTIONS_NEDELEC_NEDELECSIMPLEX_NEDELECSIMPLEXPREBASIS_HH
9#include <dune/geometry/type.hh>
15 template < GeometryType::Id geometryId,
class Field >
16 struct NedelecVecMatrix;
18 template <
unsigned int dim,
class Field>
22 typedef typename MBasisFactory::Object
MBasis;
27 typedef std::size_t
Key;
29 template <
unsigned int dd,
class FF>
35 template< GeometryType::Id geometryId >
49 MBasis *mbasis = MBasisFactory::template create<geometryId>(order+1);
50 std::remove_const_t<Object>* tmBasis =
new std::remove_const_t<Object>(*mbasis);
51 tmBasis->fill(vecMatrix);
57 template <GeometryType::Id geometryId,
class Field>
60 static constexpr GeometryType
geometry = geometryId;
86 DUNE_THROW(Dune::NotImplemented,
"High order nedelec elements are only supported by simplices in 2d and 3d");
89 FieldVector< MI, dim > x;
96 for(
unsigned int i = 0; i <
dim; ++i )
98 std::vector< MI > val( basis.
size() );
106 unsigned int notHomogen = 0;
108 notHomogen = basis.
sizes()[order-1];
111 unsigned int homogen = basis.
sizes()[order]-notHomogen;
143 int homogenTwoVariables = 0;
144 for(
int w = notHomogen; w<notHomogen + homogen; w++)
146 homogenTwoVariables++;
147 row_ = (notHomogen*
dim+homogen*(
dim+2) + homogenTwoVariables)*
dim;
156 for (
unsigned int i=0; i<notHomogen+homogen; ++i)
158 for (
unsigned int r=0; r<
dim; ++r)
160 for (
unsigned int rr=0; rr<
dim; ++rr)
164 for (
unsigned int j=0; j<
col_; ++j)
177 for (
unsigned int i=0; i<homogen; ++i)
180 MI xval = val[notHomogen+i];
183 for (
unsigned int r=0; r<
dim; ++r)
187 for (
unsigned int j=0; j<
col_; ++j)
195 for (
int w=homogen+notHomogen; w<val.size(); ++w)
197 if (val[w] == xval*x[0])
199 if (val[w] == xval*x[1])
206 int skipLastDim = xval.
z(0)>0;
212 for (
unsigned int r=0; r<
dim*(
dim-skipLastDim); ++r)
216 for (
unsigned int j=0; j<
col_; ++j)
228 for (
unsigned int r=0; r<
dim - skipLastDim; ++r)
230 int index = (r+
dim-1)%
dim;
231 for (
int w=homogen+notHomogen; w<val.size(); ++w)
233 if (val[w] == xval*x[index])
235 if (val[w] == xval*x[r])
247 for (
unsigned int i=0; i<
rows(); ++i) {
261 template <
class Vector>
262 void row(
const unsigned int row, Vector &vec )
const
264 const unsigned int N =
cols();
265 assert( vec.size() == N );
266 for (
unsigned int i=0; i<N; ++i)
Definition: bdfmcube.hh:16
void field_cast(const F1 &f1, F2 &f2)
a helper class to cast from one field to another
Definition: field.hh:157
Definition: nedelecsimplexprebasis.hh:59
NedelecVecMatrix(std::size_t order)
Definition: nedelecsimplexprebasis.hh:64
MultiIndex< dim, Field > MI
Definition: nedelecsimplexprebasis.hh:62
unsigned int row_
Definition: nedelecsimplexprebasis.hh:270
unsigned int cols() const
Definition: nedelecsimplexprebasis.hh:253
~NedelecVecMatrix()
Definition: nedelecsimplexprebasis.hh:245
MonomialBasis< geometryId, MI > MIBasis
Definition: nedelecsimplexprebasis.hh:63
unsigned int col_
Definition: nedelecsimplexprebasis.hh:270
static const unsigned int dim
Definition: nedelecsimplexprebasis.hh:61
void row(const unsigned int row, Vector &vec) const
Definition: nedelecsimplexprebasis.hh:262
static constexpr GeometryType geometry
Definition: nedelecsimplexprebasis.hh:60
unsigned int rows() const
Definition: nedelecsimplexprebasis.hh:257
Field ** mat_
Definition: nedelecsimplexprebasis.hh:271
Definition: nedelecsimplexprebasis.hh:20
static void release(Object *object)
Definition: nedelecsimplexprebasis.hh:54
MBasisFactory::Object MBasis
Definition: nedelecsimplexprebasis.hh:22
static Object * create(Key order)
Definition: nedelecsimplexprebasis.hh:36
PolynomialBasisWithMatrix< EvalMBasis, SparseCoeffMatrix< Field, dim > > Basis
Definition: nedelecsimplexprebasis.hh:24
const Basis Object
Definition: nedelecsimplexprebasis.hh:26
StandardEvaluator< MBasis > EvalMBasis
Definition: nedelecsimplexprebasis.hh:23
MonomialBasisProvider< dim, Field > MBasisFactory
Definition: nedelecsimplexprebasis.hh:21
std::size_t Key
Definition: nedelecsimplexprebasis.hh:27
Definition: nedelecsimplexprebasis.hh:31
MonomialBasisProvider< dd, FF > Type
Definition: nedelecsimplexprebasis.hh:32
Definition: basisevaluator.hh:129
Definition: monomialbasis.hh:438
unsigned int size() const
Definition: monomialbasis.hh:474
void evaluate(const unsigned int deriv, const DomainVector &x, Field *const values) const
Definition: monomialbasis.hh:496
const unsigned int * sizes(unsigned int order) const
Definition: monomialbasis.hh:463
Definition: monomialbasis.hh:789
Definition: multiindex.hh:35
int z(int i) const
Definition: multiindex.hh:89
Definition: polynomialbasis.hh:335