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;
23 template <
unsigned int dim,
class Field>
27 typedef typename MBasisFactory::Object
MBasis;
32 typedef std::size_t
Key;
34 template <
unsigned int dd,
class FF>
40 template< GeometryType::Id geometryId >
54 MBasis *mbasis = MBasisFactory::template create<geometryId>(order+1);
55 std::remove_const_t<Object>* tmBasis =
new std::remove_const_t<Object>(*mbasis);
56 tmBasis->fill(vecMatrix);
62 template <GeometryType::Id geometryId,
class Field>
65 static constexpr GeometryType
geometry = geometryId;
91 DUNE_THROW(Dune::NotImplemented,
"High order nedelec elements are only supported by simplices in 2d and 3d");
94 FieldVector< MI, dim > x;
101 for(
unsigned int i = 0; i <
dim; ++i )
103 std::vector< MI > val( basis.
size() );
111 unsigned int notHomogen = 0;
113 notHomogen = basis.
sizes()[order-1];
116 unsigned int homogen = basis.
sizes()[order]-notHomogen;
148 int homogenTwoVariables = 0;
149 for(
int w = notHomogen; w<notHomogen + homogen; w++)
151 homogenTwoVariables++;
152 row_ = (notHomogen*
dim+homogen*(
dim+2) + homogenTwoVariables)*
dim;
161 for (
unsigned int i=0; i<notHomogen+homogen; ++i)
163 for (
unsigned int r=0; r<
dim; ++r)
165 for (
unsigned int rr=0; rr<
dim; ++rr)
169 for (
unsigned int j=0; j<
col_; ++j)
182 for (
unsigned int i=0; i<homogen; ++i)
185 MI xval = val[notHomogen+i];
188 for (
unsigned int r=0; r<
dim; ++r)
192 for (
unsigned int j=0; j<
col_; ++j)
200 for (
int w=homogen+notHomogen; w<val.size(); ++w)
202 if (val[w] == xval*x[0])
204 if (val[w] == xval*x[1])
211 int skipLastDim = xval.
z(0)>0;
217 for (
unsigned int r=0; r<
dim*(
dim-skipLastDim); ++r)
221 for (
unsigned int j=0; j<
col_; ++j)
233 for (
unsigned int r=0; r<
dim - skipLastDim; ++r)
235 int index = (r+
dim-1)%
dim;
236 for (
int w=homogen+notHomogen; w<val.size(); ++w)
238 if (val[w] == xval*x[index])
240 if (val[w] == xval*x[r])
252 for (
unsigned int i=0; i<
rows(); ++i) {
266 template <
class Vector>
267 void row(
const unsigned int row, Vector &vec )
const
269 const unsigned int N =
cols();
270 assert( vec.size() == N );
271 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:64
NedelecVecMatrix(std::size_t order)
Definition nedelecsimplexprebasis.hh:69
MultiIndex< dim, Field > MI
Definition nedelecsimplexprebasis.hh:67
unsigned int row_
Definition nedelecsimplexprebasis.hh:275
unsigned int cols() const
Definition nedelecsimplexprebasis.hh:258
~NedelecVecMatrix()
Definition nedelecsimplexprebasis.hh:250
MonomialBasis< geometryId, MI > MIBasis
Definition nedelecsimplexprebasis.hh:68
unsigned int col_
Definition nedelecsimplexprebasis.hh:275
static const unsigned int dim
Definition nedelecsimplexprebasis.hh:66
void row(const unsigned int row, Vector &vec) const
Definition nedelecsimplexprebasis.hh:267
static constexpr GeometryType geometry
Definition nedelecsimplexprebasis.hh:65
unsigned int rows() const
Definition nedelecsimplexprebasis.hh:262
Field ** mat_
Definition nedelecsimplexprebasis.hh:276
Definition nedelecsimplexprebasis.hh:25
static void release(Object *object)
Definition nedelecsimplexprebasis.hh:59
MBasisFactory::Object MBasis
Definition nedelecsimplexprebasis.hh:27
static Object * create(Key order)
Definition nedelecsimplexprebasis.hh:41
PolynomialBasisWithMatrix< EvalMBasis, SparseCoeffMatrix< Field, dim > > Basis
Definition nedelecsimplexprebasis.hh:29
const Basis Object
Definition nedelecsimplexprebasis.hh:31
StandardEvaluator< MBasis > EvalMBasis
Definition nedelecsimplexprebasis.hh:28
MonomialBasisProvider< dim, Field > MBasisFactory
Definition nedelecsimplexprebasis.hh:26
std::size_t Key
Definition nedelecsimplexprebasis.hh:32
Definition nedelecsimplexprebasis.hh:36
MonomialBasisProvider< dd, FF > Type
Definition nedelecsimplexprebasis.hh:37
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:38
int z(int i) const
Definition multiindex.hh:92
Definition polynomialbasis.hh:348