← Back to team overview

dolfin team mailing list archive

Re: Function.vector() and solve()

 

2008/3/31, Anders Logg <logg@xxxxxxxxx>:
> On Mon, Mar 31, 2008 at 03:28:16PM +0200, Ola Skavhaug wrote:
>  > Anders Logg skrev den 31/03-2008 følgende:
>  > > On Mon, Mar 31, 2008 at 03:00:27PM +0200, Ola Skavhaug wrote:
>  > > > Anders Logg skrev den 31/03-2008 følgende:
>  > > > > On Mon, Mar 31, 2008 at 02:47:30PM +0200, Ola Skavhaug wrote:
>  > > > > > Anders Logg skrev den 31/03-2008 følgende:
>  > > > > > > On Mon, Mar 31, 2008 at 02:05:25PM +0200, Ola Skavhaug wrote:
>  > > > > > > > Anders Logg skrev den 31/03-2008 følgende:
>  > > > > > > > > On Mon, Mar 31, 2008 at 01:11:42PM +0200, Martin Sandve Alnæs wrote:
>  > > > > > > > > > 2008/3/31, Anders Logg <logg@xxxxxxxxx>:
>  > > > > > > > > > > On Sun, Mar 30, 2008 at 11:04:57PM +0200, Martin Sandve Alnæs wrote:
>  > > > > > > > > > >  > 2008/3/30, Anders Logg <logg@xxxxxxxxx>:
>  > > > > > > > > > >  > > On Sun, Mar 30, 2008 at 10:25:06PM +0200, Martin Sandve Alnæs wrote:
>  > > > > > > > > > >  > >  > Then solve should be fixed to expect a GenericVector, what's the problem?
>  > > > > > > > > > >  > >
>  > > > > > > > > > >  > >
>  > > > > > > > > > >  > > Yes, but that needs some thinking. The solve() function needs to check
>  > > > > > > > > > >  > >  what kind of arguments it gets and redirect to the correct backend.
>  > > > > > > > > > >  > >  For example, if (A, x, b) are PETSc objects, then it needs to call a
>  > > > > > > > > > >  > >  PETSc solver (which does not work through the GenericFoo interface).
>  > > > > > > > > > >  > >
>  > > > > > > > > > >  > >  I guess you can work your magic try { dynamic_cast<...> } here?
>  > > > > > > > > > >  >
>  > > > > > > > > > >  > Sure, that should be easy. I'll look at it tomorrow.
>  > > > > > > > > > >
>  > > > > > > > > > >
>  > > > > > > > > > > Fint!
>  > > > > > > > > > >
>  > > > > > > > > > >  solve() ligger i dolfin/la/solve.{h,cpp}
>  > > > > > > > > >
>  > > > > > > > > > Vector creates problems when testing the type of a GenericVector.
>  > > > > > > > > >
>  > > > > > > > > > All code that wants to test the type of a GenericVector
>  > > > > > > > > > will depend on Vector.
>  > > > > > > > > >
>  > > > > > > > > > Do you really want all this type trouble just to be able to write
>  > > > > > > > > >
>  > > > > > > > > >   Vector v;
>  > > > > > > > > >
>  > > > > > > > > > instead of
>  > > > > > > > > >
>  > > > > > > > > >   typedef FooVec Vector;
>  > > > > > > > > >   ...
>  > > > > > > > > >   GenericVector *v = new Vector()
>  > > > > > > > >
>  > > > > > > > > Yes! "Vector v" looks almost infinitely better than having to write
>  > > > > > > > > typedefs and calling new.
>  > > > > > > > >
>  > > > > > > > > But we don't really need to fix solve() if it's too much trouble.
>  > > > > > > > > If we don't want to dynamically resolve the representation of matrices
>  > > > > > > > > and vectors, then there will be some functions such as assemble() that
>  > > > > > > > > work for any type of vector (GenericVector) and there will be others
>  > > > > > > > > that work only for a specific implementation.
>  > > > > > > > >
>  > > > > > > > > So, if you have a Vector, then you can use it with everything:
>  > > > > > > > >
>  > > > > > > > >    assemble()
>  > > > > > > > >    solve()
>  > > > > > > > >    LUSolver
>  > > > > > > > >    KrylovSolver
>  > > > > > > > >
>  > > > > > > > > but if you have a FooVector, then you must use it with the
>  > > > > > > > > Foo solvers:
>  > > > > > > > >
>  > > > > > > > >    FooLUSolver
>  > > > > > > > >    FooKrylovSolver
>  > > > > > > > >
>  > > > > > > > > So, if it's too much trouble to resolve the backend dynamically, then
>  > > > > > > > > we don't really need to do it.
>  > > > > > > > >
>  > > > > > > > > The basic idea here is that if you don't care about the backend,
>  > > > > > > > > just use Matrix and Vector and everything is simple (assemble, solve).
>  > > > > > > > >
>  > > > > > > > > But if you want to use some special backend, then use FooMatrix,
>  > > > > > > > > FooVector and FooSolver.
>  > > > > > > >
>  > > > > > > > That is the problem in having a simpler design, omitting the Matrix-wrapper
>  > > > > > > > construction, and using only a typedef? F.ex.
>  > > > > > > >
>  > > > > > > > typedef Matrix uBlasMatrix;
>  > > > > > > >
>  > > > > > > > If we are consistent in the DOLFIN library, only communicating GenericFoo
>  > > > > > > > across the interfaces, the chances of non-compatible code being propagated
>  > > > > > > > into the library is minimal. The only place the additional methods of a
>  > > > > > > > particular backend can be used is in user space (and hence not our
>  > > > > > > > concern).
>  > > > > > > >
>  > > > > > > > This will give the nice syntax Anders wants, and conforms with standard
>  > > > > > > > oo-design. The cost is that methods added to a backend will be accessible
>  > > > > > > > through the Matrix typedef, and this might lead to non-portable _user_ code.
>  > > > > > > >
>  > > > > > > > ??
>  > > > > > > >
>  > > > > > > > Ola
>  > > > > > >
>  > > > > > > Yes, the problem with what we had before, namely
>  > > > > > >
>  > > > > > >   typedef FooMatrix Matrix;
>  > > > > > >   typedef FooVector Vector;
>  > > > > > >
>  > > > > > > was that it allowed writing code depending on Matrix and Vector which
>  > > > > > > was not portable (when changing the typedef).
>  > > > > > >
>  > > > > > > For example, it was fully possible to do
>  > > > > > >
>  > > > > > >   Matrix A;
>  > > > > > >   MatSetValues(A.mat(), ...);
>  > > > > > >
>  > > > > > > but this wouldn't work unless DOLFIN was compiled with
>  > > > > > >
>  > > > > > >   typedef PETScMatrix Matrix;
>  > > > > > >
>  > > > > > > So, the purpose of the current design with a class Matrix inheriting
>  > > > > > > from GenericMatrix and having a FooMatrix pointer is to make sure that
>  > > > > > > Matrix can only do the things defined in GenericMatrix, nothing less
>  > > > > > > (which is handled by the compiler) and nothing more (which we must
>  > > > > > > take care of manually).
>  > > > > >
>  > > > > > I understand this. My point is that if we use GenericFoo consistently though
>  > > > > > the interface of DOLFIN, the (mis)use of Matrix will only live in user code.
>  > > > > >
>  > > > > > Ola
>  > > > >
>  > > > > Do mean that we don't need the Matrix and Vector classes, that we can
>  > > > > use typedefs instead since we will only break someone else's code, not
>  > > > > our own?
>  > > >
>  > > > Almost. They will break their own code, and will probably be able to fix it.
>  > > > But your point is valid; we wont need the Matrix and Vector classes. Typedefs
>  > > > should suffice.
>  > >
>  > > It wasn't my point, it was yours. :-) I'm not completely convinced we
>  > > should replace the classes with typedefs yet.
>  > >
>  > > Having the classes there instead of typedefs is just a safety measure
>  > > and the question is how much we want to pay for it (in terms of
>  > > maintenance). It seem the cost is pretty small (since they are already
>  > > there). We just need to add a few more functions at some point in the
>  > > near future.
>  >
>  > I don't think the cost is small. If we need to resolve the actual type of a
>  > Matrix everywhere we get a GenericMatrix, this will give circular dependencies.
>  > Having a clean cut object oriented design will enable us to simply
>  > dynamic_cast the generic base class pointer to the particular backend instance
>  > and do the work. It is very ambitions to not use dynamic casts, and it will
>  > for sure lead to much work. On the other hand, if someone implements a backend
>  > and is able to implement all functionality by only using the base class
>  > interface, the design we probose will not stop that.
>  >
>  > Ola
>
>
> You mean that if the Matrix that comes in to some function that wants
>  a GenericMatrix is actually a PETScMatrix, then it will be easy to do
>  the dynamic cast, but if it's a Matrix, then it will be very
>  difficult? ok, I can buy that.

It's not very difficult, but cumbersome, and it can't be done in f.ex.
uBlasVector since uBlasVector can't depend on Vector, since that would
cause a circular dependency.

>  But in that case, do we need to do dynamic cast at all? Instead of
>
>   void foo(GenericMatrix& A)
>   {
>      try { PETScMatrix }
>      try { uBlasMatrix }
>   }

(Sidenote: there's no try involved in dynamic_cast. It just returns NULL.)

>  we can just overload:
>
>   void foo(PETScMatrix& A) {}
>   void foo(uBlasMatrix& A) {}
>
>
>  ?

The compiler can only select between these foo's at compile-time.

A simple use case which requires type testing:
I have a DiscreteFunction, and want to pass the underlying vector to
some linear algebra code using petsc, ublas or trilinos directly. Then
I need to get the underlying object and cast it. This is not a very
exotic problem!

-- 
Martin


Follow ups

References