dolfin team mailing list archive
-
dolfin team
-
Mailing list archive
-
Message #02819
Re: uBlas Krylov solver and preconditioner
On Tue, Jul 04, 2006 at 08:26:47AM +0200, Garth N. Wells wrote:
> On Mon, 2006-07-03 at 17:31 +0200, Anders Logg wrote:
> > I've extended the uBlas Krylov solver so that it now works with
> > matrices of type uBlasKrylovMatrix, only defined by having two
> > functions size() and mult(). I had to make some minor changes to avoid
> > using axpy_prod() but the tests I made showed now reduction in
> > speed. (There is one extra copy of values into a vector compared to
> > the original code.)
> >
>
> Make sure that you use noalias() whenever possible when working with
> uBlas data types. In the code, r += b should be noalias(r) += b or
> r.plus_assign(b). Otherwise, a copy of r is make internally which gets
> expensive for large vectors.
ok, I'll fix.
> > I have temporarily disabled the preconditioning. See all comments
> >
> > // FIXME: Preconditioner temporarily disabled
> >
> > This was necessary since uBlasPreconditioner assumes uBlasSparseMatrix.
> > We need to make it so that it only assumes uBlasKrylovMatrix. Then we
> > can have a special subclass, like
> >
> > class uBlasILUPreconditioner : public uBlasPreconditioner
> >
> > or similarly. This class should then take in the constructor as
> > argument a uBlasSparseMatrix. Then we can overload solve()
> > uBlasKrylovSolver so we have two versions: one with and one without an
> > extra preconditioner argument.
> >
> > I don't think it's possible to let the Krylov solver create a
> > preconditioner on its own since it only knows how to do mult(), that
> > is, it only knows that the matrix is uBlasKrylovMatrix. So it can't
> > create a preconditioner that assumes it can do M.assign(A).
> >
> > Does that make sense?
> >
>
> Sounds ok. The preconditioner interface does need to be sorted out. My
> first implementation was very simple just to get something working.
>
> Garth
ok, I'll make an attempt to get the preconditioning going again and
then you can go ahead with the templating. :-)
/Anders
References