[dune-pdelab] VectorBackend, OnTheFlyOperator and storage

Peter Bastian peter.bastian at iwr.uni-heidelberg.de
Sat May 4 16:04:50 CEST 2013


Dear Christian,

Am 03.05.2013 um 22:33 schrieb Christian Engwer <christian.engwer at uni-muenster.de>:

> Dear all,
> 
> I recently stumbled upon some problems:
> 
> a) I wanted to reactivate the OnTheFlyOeprator...
>   Originally the template parameters are the
>   ISTLBlockVectorContainers. Since the backend doesn't inherit the
>   underlying container anymore (remember all the trouble with
>   inconsistent methods, some working on the blocks, some on the flat
>   matrix) this breaks the whole operator, as the solver passes a
>   BlockVector to the apply method of the operator, but the
>   operator expects an ISTLBlockVectorContainer. My first idea was to
>   change the template parameters and to the underlying BlockVector,
>   but this didn't work, because the operator calls go.jacobian_apply,
>   which in turn operates on the ISTLBlockVectorContainer.
You can (and should) run the solver with the ISTLBlockVectorContainer
as template parameter
(it has all the necessary methods). This is done in the parallel backends
in order to use the pdelab communication or things like set_constrained_dofs.

But I admit that although of course on the one hand having the
replacable backends is fine the mixture of container and underlying
object from the LA library is a pain in the ass. I have the feeling that
it currently has grown into something that works but it does not have
a clear concept. This becomes clear when you look into the
parallel backends.


> 
>   What can we do now?
> 
>   - first thing seems to me, that we have to seperate the container
>     from the storage, this means that ISTLBlockVectorContainer does
>     not store the BlockVector directly, but just a shared_ptr. This
>     way we can later update the storage object or we can create a
>     container, which wraps an existing BlockVector.
>   - then the next problem I faces was that this doesn't work either,
>     because shared_ptr has to hold mutable BlockVector, even if we
>     then only use a const reference of the container. On the other
>     hand the x parameter of the apply method is a const reference,
>     and we can't change this because we have to fulfill a virtual
>     interface.
>   - I'm not sure what would be a nice solution. We could try to add a
>     second kind of container, which is only for const objects, but
>     then we have to be able to convert between both...
> 
>  Perhaps someone has a completely different idea. Yeasterday in the
>  train I didn't come up with an other idea ;-)
> 
> b) I realized that, although we removed the inheritance for the vector
>   backends, we still use it for the matrix backend. I would propose
>   to change this for several reasons:
> 
>   - I prefer the clear seperation of backend and underlying LA lib
>   - we avoid potential problems with blocked vs flat view (as we had
>     them with the vector)
>   - we gain more flexibility regarding the implementation (see me
>     following point ©)
Now its not flat against block but we have the multi index stuff implemented
by Steffen. Clearly, we should omit the inheritance. Is it really used somewhere?
I can't think about a place.

> 
> c) Oliver asked we, why it is not possible to create an ISTL Matrix
>   and then later use PDELab to fill this matrix. While we completely
>   stay in the PDELab context one might argue that this is not
>   necessary and even the wrong approach, I think in situations where
>   we mix different codes Oli has a point. If we implement (b) and
>   then don't store the matgrix obect directly in the backend, but use
>   a shared_ptr, we could easily allow the user to allocate the matrix
>   before hand and then warp this matrix in a backend container.
Such a possibility is clearly something one should have. I just committed
a constructor on the matrix backend yesterday that allows you to
initialize a MatrixContainer with a given ISTL Matrix (instead of
a grid operator). In my case this 
is the triple product A_C = R A_DG A^T which is constructed using ISTL
but is then later to be used within PDELab (e.g. to insert essential boundary
conditions). 

Best,

Peter

------------------------------------------------------------
Peter Bastian
Interdisziplinäres Zentrum für Wissenschaftliches Rechnen
Universität Heidelberg
Im Neuenheimer Feld 368
D-69120 Heidelberg
Tel: 0049 (0) 6221 548261
Fax: 0049 (0) 6221 548884
email: peter.bastian at iwr.uni-heidelberg.de
web: http://conan.iwr.uni-heidelberg.de/people/peter/





More information about the dune-pdelab mailing list