dolfin team mailing list archive

["Garth N. Wells" <g.n.wells@xxxxxxxxxx>]

On Thu, 2006-06-01 at 20:37 +0200, Anders Logg wrote:
> On Thu, Jun 01, 2006 at 11:08:40AM +0200, Garth N. Wells wrote:
> > On Thu, 2006-06-01 at 10:39 +0200, DOLFIN wrote:
> > > One or more new changesets pushed to the primary DOLFIN repository.
> > > A short summary of the last three changesets is included below.
> > > 
> > > changeset:   1981:5dd79b41a3d85cde459ff80e85edd968aeb90f78
> > > tag:         tip
> > > user:        "Garth N. Wells <g.n.wells@xxxxxxxxxx>"
> > > date:        Thu Jun  1 09:33:27 2006 +0200
> > > files:       src/kernel/la/uBlasKrylovSolver.cpp
> > > description:
> > > Add uBlasKrylovSolver.cpp
> > > 
> > > 
> > > changeset:   1980:f2770086f902eff6b19f85f6367d45054a998f09
> > > user:        "Garth N. Wells <g.n.wells@xxxxxxxxxx>"
> > > date:        Thu Jun  1 09:31:16 2006 +0200
> > > files:       src/kernel/la/dolfin/uBlasKrylovSolver.h src/kernel/la/dolfin/uBlasLUSolver.h
> > > description:
> > > Add uBlasKrylovSolver and uBlasLUSolver.
> > > 
> > > 
> > > changeset:   1979:6fef8c669a6a21635a3ff05a341d0b4fa9d2ae11
> > > user:        "Garth N. Wells <g.n.wells@xxxxxxxxxx>"
> > > date:        Thu Jun  1 09:15:15 2006 +0200
> > > files:       src/kernel/la/DenseVector.cpp src/kernel/la/GMRES.cpp src/kernel/la/KrylovSolver.cpp src/kernel/la/LU.cpp src/kernel/la/Makefile.am src/kernel/la/Makefile.in src/kernel/la/PETScKrylovSolver.cpp src/kernel/la/PETScLU.cpp src/kernel/la/dolfin/DenseVector.h src/kernel/la/dolfin/GMRES.h src/kernel/la/dolfin/KrylovSolver.h src/kernel/la/dolfin/LU.h src/kernel/la/dolfin/Makefile.am src/kernel/la/dolfin/Makefile.in src/kernel/la/dolfin/PETScKrylovSolver.h src/kernel/la/dolfin/PETScLU.h src/kernel/la/dolfin/Preconditioner.h src/kernel/la/dolfin/dolfin_la.h
> > > description:
> > > Add prefix "PETSc" to PETSc-based solvers and data structures. "KrylovSolver" and "LU" are typedefs to "PETScKrylovSolver" and "PETScLU" if PETSc is enabled, otherwise they are aliases to uBlasKrylovSolver and uBlasLUSolver.
> > > 
> > > This should make things simple throughout the code. If a PETSc-specific feature is required, use PETScFoo and add an "#ifdef", otherwise just use Foo.
> > > 
> > 
> > I'm having some trouble with the base class LinearSolver which has the
> > virtual function
> > 
> >  virtual uint solve(const SomeMatrix& A, SomeVector& x, const
> > SomeVector& b) = 0;
> >   
> > The problem is that I'd like to derive different solvers from it which
> > accept different matrix and vector types. For example
> > 
> >   PETScLU::solve(const PETScMatrix& A, PETScVector& x, const
> > PETScVector& b)
> >   uBlasLUSolver::solve(const uBlasMatrix& A, uBlasVector& x, const
> > uBlasVector& b)
> >  
> > Any ideas how to do this? I think that the Barton-Nackman trick (where
> > the base class is a template) which we tested for matrices might do the
> > job, at the expense of the complication introduced by the templates. 
> > 
> > Garth
> 
> How about the following:
> 
> In many places, we have adopted the practive that we place algorithms
> in separate classes to be able to separate data structures and
> algorithms and to not clutter the implementation. One example of this
> is the mesh library, where we have the class Mesh and a separate class
> MeshRefinement for mesh refinement. One may thus refine a mesh by
> 
>     MeshRefinement::refine(mesh);
> 
> In addition, we have also added a member function to Mesh, which just
> does MeshRefinement::refine(*this) so one may do
> 
>     mesh.refine();
> 
> We could follow the same principle for linear algebra, so that one may
> do either
> 
>     GMRES::solve(A, x, b);
> 
> or
> 
>     A.solveGMRES(x, b);
> 
> There can also be a default function solve() that calls GMRES for some
> matrix types and LU for other matrix types. What we would need to do
> would be to have a class PETScGMRES with a static function solve():
> 
>     static PETScGMRES::solve(const PETScSparseMatrix& A, ...);
> 
> We would also add a function solveGMRES() to PETScSparseMatrix that
> just calls this function with *this. Furthermore, we would have a
> class GMRES that takes a GenericMatrix:
> 
>     static GMRES::solve(GenericMatrix& A, ...);
> 
> This function would just call GenericMatrix::solveGMRES(), which would
> resolve to PETScSparseMatrix::solveGMRES() in the case of a PETSc
> sparse matrix, which in turn calls PETScGMRES::solve(*this).
>

Haven't had a chance to go through the above carefully yet, but it
sounds good. Static functions could fix things. The weakness is that the
parameter system isn't as neat as it isn't, for example, possible (in a
simple way) to associate different tolerances for different solver
objects of the same type. 

> When it comes to LinearSolver, we might just remove it. I grepped
> through the code and it seems we don't use it.
> 

I tried first not to derive the solver class from LinearSolver, but ran
into problems with NewtonSolver which has a pointer to a LinearSolver.
Hopefully the above will solve this.

Garth


> /Anders
> 
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev