Caching could be done as you suggest - we don't yet have restriction/prolongation for functionals available. I'll hopefully be able to merge the branch from dune-localfunctions then this branch can also be merged. Perhaps you could then test this with acfem (if you haven't done that yet). Concerning caching: I will not get around to that for a bit - so perhaps we should just put it into an issue for now, I don't think it should be a show stopper for this MR?

I'd like to merge this by end of week

**Andreas Dedner**
(12a9f35b)
*
at
10 Apr 17:04
*

make it possible to use a general vector class to store interpolation

*
... and
7 more commits
*

**Andreas Dedner**
(edfc42b3)
*
at
09 Apr 22:57
*

Correct the return type of the interface methods `size`

and `order`

to return `size_t`

instead of `unsigned int`

Also fix two unsigned int warnings

Correct the return type of the interface methods `size`

and `order`

to return `size_t`

instead of `unsigned int`

**Andreas Dedner**
(73be92d8)
*
at
09 Apr 22:52
*

fixed some unsinged/signed comparison warnings

*
... and
1 more commit
*

I had a look at the issue with the static variable. The problem seems to be the following code (dune/localfunctions/raviartthomas/raviartthomassimplex/raviartthomassimplexinterpolation.hh:311):

```
typedef Dune::QuadratureRule<Field, dimension> Quadrature;
typedef Dune::QuadratureRules<Field, dimension> QuadratureRules;
const Quadrature &elemQuad = QuadratureRules::rule
( geoType, 2*order_+1 );
```

Looking at `elemQuad`

in gdb leads to

```
(gdb) print elemQuad
$26 = (const Dune::RaviartThomasL2Interpolation<2u, double>::Quadrature &) @0x75e678: {<std::vector<Dune::QuadraturePoint<double, 2>, std::allocator<Dune::QuadraturePoint<double, 2> > >> = std::vector of length 4, capacity 4 = {{
local = {<Dune::DenseVector<Dune::FieldVector<double, 2> >> = {<No data fields>}, _data = {
_M_elems = {0, 0}}}, weight_ = 0}, {
local = {<Dune::DenseVector<Dune::FieldVector<double, 2> >> = {<No data fields>}, _data = {
_M_elems = {0, 0}}}, weight_ = 0}, {
local = {<Dune::DenseVector<Dune::FieldVector<double, 2> >> = {<No data fields>}, _data = {
_M_elems = {0, 0}}}, weight_ = 0}, {
local = {<Dune::DenseVector<Dune::FieldVector<double, 2> >> = {<No data fields>}, _data = {
_M_elems = {0, 0}}}, weight_ = 0}},
_vptr.QuadratureRule = 0x73b380 <vtable for Dune::QuadratureRule<double, 2>+16>,
geometry_type = {dim_ = 2 '\002', none_ = false, topologyId_ = 0}, delivered_order = 0}
```

i.e. the set up of the quadratures in `dune-geometry`

has not been carried out (also static). The following version of the code works because the `dune-geometry`

quadratures are initialized before the `instance`

method is called. I don't know if there is anything I can do about this....

```
#include <config.h>
#include <dune/geometry/referenceelements.hh>
#include <dune/localfunctions/raviartthomas.hh>
constexpr int dim = 2;
constexpr int k = 1;
//some class
struct RTstatic {
// rt as static member: do not work
static auto instance()
{
static const Dune::RaviartThomasSimplexLocalFiniteElement<dim,double,double> rt(Dune::GeometryTypes::simplex(dim),k);
return rt;
}
static const Dune::RaviartThomasSimplexLocalFiniteElement<dim,double,double> rt;
};
// initialize the static member (this causes an error)
// const Dune::RaviartThomasSimplexLocalFiniteElement<dim,double,double> RTstatic::rt(Dune::GeometryTypes::simplex(dim),k);
int main(int argc, char const *argv[])
{
// rt as normal variable: works
Dune::RaviartThomasSimplexLocalFiniteElement<dim,double,double> rt(Dune::GeometryTypes::simplex(dim),k);
RTstatic::instance().size();
return 0;
}
```

Sorry, I lost track of this issue: I (hopefully) fixed the return type issue (see !121 (closed))

`size`

and `order`

to return `size_t`

instead of `unsigned int`

**Andreas Dedner**
(edfc42b3)
*
at
09 Apr 11:18
*

fixed some unsinged/signed comparison warnings

*
... and
8 more commits
*

I'm not quite sure how you would be planning to do this. Take a 1d case [0,2h] split into [0,h],[h,2h] and as interpolation on the father the point evaluation in the midpoint, e.g., in 0.5. The values on the children could just be two constant 1 and 2 then

```
\sum_C I_F(u_C\circ (G_{CF})^{-1}(0.5)) =
u_{C_1}((G^{C_1F})^{-1}(0.5)) +
u_{C_2}((G^{C_2F})^{-1}(0.5)) =
u_{C_1}(1) + u_{C_2}(0) = 1 + 2 = 3
```

If the quadrature is at two points 0.25,0.75 with weights 0.5 then I get

```
0.5(\sum_C I_F(u_C\circ (G_{CF})^{-1}(0.25)) +
0.5(\sum_C I_F(u_C\circ (G_{CF})^{-1}(0.75)) =
0.5(u_{C_1}((G^{C_1F})^{-1}(0.25)) + 0) +
0.5(0+u_{C_2}((G^{C_2F})^{-1}(0.75))) =
0.5\cdot 1 + 0.5\cdot 2 = 1.5
```

The second result is what I would expect but I don't see how to fix the first case. Or are we talking at cross purpose?

**Andreas Dedner**
(e93d43ec)
*
at
08 Apr 22:46
*

avoid usage of dune-fempy global/local grid function instead use du...

**Andreas Dedner**
(8a67e439)
*
at
08 Apr 22:39
*

added name to grid function class decorator

**Andreas Dedner**
(33aae228)
*
at
08 Apr 22:37
*

add specialization for DGLFE space - used same base implementation ...

I just realized that there is a `DiscontinuousLocalFiniteElementSpace`

for that the `DefaultLocalRestrictProlong`

has to also be specialized as well.
I'm just not quite sure the summation works - at least not guaranteed. Just take for example a quadrature that uses the three edge midpoints (fine for quadratics). The values at those points would be added up multiple times if we did the restriction in `restrictLocal`

.

That would work as well and I had considered that. The reason was discontinuous functions the approach now at least gives an averaging where the result is independent of the order in which the children are traversed. Also one doesn't really have the option of a custom quadrature at least as long as the provided interpolation is used. But my final reason for using the above implementation was that the additional cost isn't really that bad (could be testing on 8 children for y inside I guess...)

@martin.alkaemper those were typos, sorry. I've corrected them. Does it make sense now?

**Andreas Dedner**
(df6e904f)
*
at
07 Apr 23:18
*

added a standard constructor (gv,interface,direction) to the LocalF...