[dune-functions] Binding to non-leaf entities in dune-functions 2.6

[Carsten Gräser]

Hi,

Am 09.11.18 um 15:07 schrieb Felix Mueller:
[...]
> /auto lv = basis.localView();
> // loop over all leaf elements
> for (auto leafElement : leafElements) {
>   lv.bind(leafElement);
>   const auto& node = lv.tree();
>   const auto& fe = node.finiteElement();
>   auto feSize = fe.size();
> 
>   for (std::size_t i = 0; i < feSize; ++i) {
>     auto index = lv.index(node.localIndex(i);
>     // save vector[index] for later use
what exactly is saved here?

>     // ...
>   }
> }
> 
> lv.bind(coarseElement);
> const auto& node = lv.tree();
> const auto& fe = node.finiteElement();
> // do something with fe
> // .../
This is not allowed. Notice that a basis is always associated
to a grid view and using localView.bind(element) for any element
not contained in the grid view leads to undefined behavior.
> Note that /coarseElement/ does not access the basis's index set here.
The coarseElement itself does not, but lv.bind(coarseElement) does.
In fact this pre-computes and caches all indices such that lv.index()
is guaranteed to be fast. Hence it will in most cases access to grid
views index set with the given element, which again leads to undefined
behavior. Here you can see different variants of undefined behaviour:
UGGrid silently returns values which may or may not have any meaning,
while AlbertaGrid fails.
> Deactivating the assertion prevents the crash and makes my code run just
> fine.
It may run without any error message, but will not produce
what you intend to do.
> This tries to access the index set on /coarseElement/, which is not
> defined and not needed within my code either.
Yes, that's on purpose.

> *Dune-functions 2.5* had a seperate binding mechanism for the local view
> and index sets, providing only the method
This one would not lead to undefined behavior of the grid due
to the misuse of the indexset. However, the behavior of the
basis is still undefined. In it would be perfectly OK if
localView.bind(element) terminates with an error if the element
is not contained in the grid view. The only reason why it does
not do so is to avoid to possibly costly check.

The upshot is: Don't do this even if it does not produce
error messages. If you want to evaluate a shape function
on an element not contained in the grid you first have
to transform the coordinates to the proper element
using geometryInFather().

Best,
Carsten

> Thanks for your attention and helpful answers!

> 
> Felix Müller
> 
> 
> PS: You will find a small file reproducing the error attached. A
> filename for the AlbertaReader must be provided as argument to the main
> method. I have used *dune-grid/doc/grids/amc/grid-2-2.amc* as grid file
> for my tests.
> 
> 
> _______________________________________________
> dune-functions mailing list
> dune-functions at lists.dune-project.org
> https://lists.dune-project.org/mailman/listinfo/dune-functions
> 


-- 
Prof. Dr. Carsten Gräser
Freie Universität Berlin
Institut für Mathematik
Arnimallee 6
14195 Berlin, Germany
phone: +49 30 838 72637
fax  : +49 30 838 472637
email: graeser at mi.fu-berlin.de
URL  : http://page.mi.fu-berlin.de/graeser

-------------- next part --------------
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 488 bytes
Desc: OpenPGP digital signature
URL: <https://lists.dune-project.org/pipermail/dune-functions/attachments/20181120/913427dc/attachment.sig>