5#ifndef DUNE_LOCALFUNCTIONS_NEDELEC_NEDELECSIMPLEX_NEDELECSIMPLEXPREBASIS_HH
6#define DUNE_LOCALFUNCTIONS_NEDELEC_NEDELECSIMPLEX_NEDELECSIMPLEXPREBASIS_HH
11#include <dune/geometry/type.hh>
17 template < GeometryType::Id geometryId,
class Field >
18 struct NedelecVecMatrix;
20 template <
unsigned int dim,
class Field>
24 typedef typename MBasisFactory::Object
MBasis;
29 typedef std::size_t
Key;
31 template <
unsigned int dd,
class FF>
37 template< GeometryType::Id geometryId >
51 MBasis *mbasis = MBasisFactory::template create<geometryId>(order+1);
52 std::remove_const_t<Object>* tmBasis =
new std::remove_const_t<Object>(*mbasis);
53 tmBasis->fill(vecMatrix);
59 template <GeometryType::Id geometryId,
class Field>
62 static constexpr GeometryType
geometry = geometryId;
88 DUNE_THROW(Dune::NotImplemented,
"High order nedelec elements are only supported by simplices in 2d and 3d");
91 FieldVector< MI, dim > x;
98 for(
unsigned int i = 0; i <
dim; ++i )
100 std::vector< MI > val( basis.
size() );
108 unsigned int notHomogen = 0;
110 notHomogen = basis.
sizes()[order-1];
113 unsigned int homogen = basis.
sizes()[order]-notHomogen;
145 int homogenTwoVariables = 0;
146 for(
int w = notHomogen; w<notHomogen + homogen; w++)
148 homogenTwoVariables++;
149 row_ = (notHomogen*
dim+homogen*(
dim+2) + homogenTwoVariables)*
dim;
158 for (
unsigned int i=0; i<notHomogen+homogen; ++i)
160 for (
unsigned int r=0; r<
dim; ++r)
162 for (
unsigned int rr=0; rr<
dim; ++rr)
166 for (
unsigned int j=0; j<
col_; ++j)
179 for (
unsigned int i=0; i<homogen; ++i)
182 MI xval = val[notHomogen+i];
185 for (
unsigned int r=0; r<
dim; ++r)
189 for (
unsigned int j=0; j<
col_; ++j)
197 for (
int w=homogen+notHomogen; w<val.size(); ++w)
199 if (val[w] == xval*x[0])
201 if (val[w] == xval*x[1])
208 int skipLastDim = xval.
z(0)>0;
214 for (
unsigned int r=0; r<
dim*(
dim-skipLastDim); ++r)
218 for (
unsigned int j=0; j<
col_; ++j)
230 for (
unsigned int r=0; r<
dim - skipLastDim; ++r)
232 int index = (r+
dim-1)%
dim;
233 for (
int w=homogen+notHomogen; w<val.size(); ++w)
235 if (val[w] == xval*x[index])
237 if (val[w] == xval*x[r])
249 for (
unsigned int i=0; i<
rows(); ++i) {
263 template <
class Vector>
264 void row(
const unsigned int row, Vector &vec )
const
266 const unsigned int N =
cols();
267 assert( vec.size() == N );
268 for (
unsigned int i=0; i<N; ++i)
Definition bdfmcube.hh:18
void field_cast(const F1 &f1, F2 &f2)
a helper class to cast from one field to another
Definition field.hh:159
Definition nedelecsimplexprebasis.hh:61
NedelecVecMatrix(std::size_t order)
Definition nedelecsimplexprebasis.hh:66
MultiIndex< dim, Field > MI
Definition nedelecsimplexprebasis.hh:64
unsigned int row_
Definition nedelecsimplexprebasis.hh:272
unsigned int cols() const
Definition nedelecsimplexprebasis.hh:255
~NedelecVecMatrix()
Definition nedelecsimplexprebasis.hh:247
MonomialBasis< geometryId, MI > MIBasis
Definition nedelecsimplexprebasis.hh:65
unsigned int col_
Definition nedelecsimplexprebasis.hh:272
static const unsigned int dim
Definition nedelecsimplexprebasis.hh:63
void row(const unsigned int row, Vector &vec) const
Definition nedelecsimplexprebasis.hh:264
static constexpr GeometryType geometry
Definition nedelecsimplexprebasis.hh:62
unsigned int rows() const
Definition nedelecsimplexprebasis.hh:259
Field ** mat_
Definition nedelecsimplexprebasis.hh:273
Definition nedelecsimplexprebasis.hh:22
static void release(Object *object)
Definition nedelecsimplexprebasis.hh:56
MBasisFactory::Object MBasis
Definition nedelecsimplexprebasis.hh:24
static Object * create(Key order)
Definition nedelecsimplexprebasis.hh:38
PolynomialBasisWithMatrix< EvalMBasis, SparseCoeffMatrix< Field, dim > > Basis
Definition nedelecsimplexprebasis.hh:26
const Basis Object
Definition nedelecsimplexprebasis.hh:28
StandardEvaluator< MBasis > EvalMBasis
Definition nedelecsimplexprebasis.hh:25
MonomialBasisProvider< dim, Field > MBasisFactory
Definition nedelecsimplexprebasis.hh:23
std::size_t Key
Definition nedelecsimplexprebasis.hh:29
Definition nedelecsimplexprebasis.hh:33
MonomialBasisProvider< dd, FF > Type
Definition nedelecsimplexprebasis.hh:34
Definition basisevaluator.hh:131
Definition monomialbasis.hh:440
unsigned int size() const
Definition monomialbasis.hh:476
void evaluate(const unsigned int deriv, const DomainVector &x, Field *const values) const
Definition monomialbasis.hh:498
const unsigned int * sizes(unsigned int order) const
Definition monomialbasis.hh:465
Definition monomialbasis.hh:780
Definition multiindex.hh:37
int z(int i) const
Definition multiindex.hh:91
Definition polynomialbasis.hh:348