← Back to team overview

dolfin team mailing list archive

Re: [HG] Add uBlasKrylovSolver.cpp

 

On Tue, 2006-06-06 at 13:50 +0200, Anders Logg wrote:
> On Tue, Jun 06, 2006 at 12:51:58PM +0200, Garth N. Wells wrote:
> > On Sat, 2006-06-03 at 10:26 +0200, Anders Logg wrote:
> > 
> > > > > > 
> > > > > > Just done some test and this works fine, but the function cannot be a
> > > > > > pure virtual function. 
> > > > > 
> > > > > Strange. The following compiles for me:
> > > > > 
> > > > > class Abstract
> > > > > {
> > > > > public:
> > > > > 
> > > > >   virtual void foo() = 0;
> > > > > 
> > > > > };
> > > > > 
> > > > > class Foo
> > > > > {
> > > > > public:
> > > > > 
> > > > >   virtual void bar(Abstract& abstract) = 0;
> > > > > 
> > > > > };
> > > > > 
> > > > > class Bar : public Foo
> > > > > {
> > > > > public:
> > > > > 
> > > > >   void bar(Abstract& abstract) {}
> > > > > 
> > > > > };
> > > > > 
> > > > > int main()
> > > > > {
> > > > >   Bar bar;
> > > > > 
> > > > >   return 0;
> > > > > }
> > > > > 
> > > > > Is this not what you mean?
> > > > > 
> > > > 
> > > > This doesn't work for me. Can you send the actual code? I'll try it out.
> > > > 
> > > > Garth
> > > 
> > > The code above is the actual code. I've attached it. Just compile with
> > > with
> > > 
> > >     g++ main.cpp
> > > 
> > > Works for me with g++ 4.0.3.
> > > 
> > > /Anders
> > > 
> > > 
> > 
> > This works for me, but the problem I have is slightly different. The
> > problem I have is that the equivalent of Bar::bar (a solver) does't want
> > "Abstract" as an argument, rather "SomeAbstract" type. I've pasted some
> > code below. The solution seems be not to make the functions in
> > LinearSolver pure virtual functons.
> > 
> > Garth
> > 
> > 
> > class GenericMatrix
> > {
> > public:
> >   
> >   virtual void foo() = 0;
> >   
> > };
> > 
> > 
> > class SomeMatrix : public GenericMatrix
> > {
> > public:
> >   
> >   
> > };
> > 
> > class LinearSolver
> > {
> > public: 
> >   
> >   virtual void bar(GenericMatrix& matrix) = 0;
> >   
> > };
> > 
> > class LU : public LinearSolver
> > {
> > public:
> >   
> >   // This doesn't work
> >   void bar(SomeMatrix& matrix) {}
> > 
> >   // This works
> > //  void bar(GenericMatrix& matrix) {}
> >   
> > };
> > 
> > int main()
> > {
> >   LinearSolver* solver;
> > 
> >   solver = new LU;
> >   
> >   return 0;
> > } 
> 
> ok, I see the problem.
> 
> We will need to do something like the following:
> 
>   class LinearSolver
>   {
>   public: 
>     virtual void solve(GenericMatrix& A, ...) = 0;
>   };
> 
>   class GenericMatrix
>   {
>   public:
>     virtual void solveLU(Vector& x, const Vector& b)
>     {
>       dolfin_info("Informative error message, not implemented.");
>     }
>   }
> 
>   class LU : public LinearSolver
>   {
>   public: 
>     void solve(GenericMatrix& A, ...)
>     {
>       A->solveLU(x, b);
>     }
>   };
> 
>   class PETScSparseMatrix : public GenericMatrix
>   {
>   public:
>     void solveLU(Vector& x, const Vector& b)
>     {
>       PETScLUSolver solver;
>       solver.solve(this, x, b);
>     }
>   };
> 
> PETScLUSolver does not inherit from LinearSolver. It could inherit
> from PETScLinearSolver if we want.
> 

A simpler solution is to not make the virtual function in LinearSolver
pure. This generates a compile-time error if the declared solver doesn't
accept the matrix/vector types given to it.



> Thus, all the work takes place in the specialized classes. We have
> PETScSparseMatrix and we have PETScLUSolver (or PETScSparseLUSolver?).
> 
> The classes LinearSolver and LU, GMRES etc are just very small classes
> that we add for convenience.

The matrix classes will become very large. I think that we should keep
the solvers separate from the matrix classes because the solvers are
independent of whether a matrix is dense or sparse in the case of PETSc
and uBlas.  


Garth

> 
> When it comes to parameters, we could maybe add a ParameterList as an
> extra parameter to the solve function in GenericMatrix. Thus,
> LinearSolver would inherit from Parametrized so one can do
> 
>     Matrix A; // defaults to PETScSparseMatrix
>     LU solver;
>     solver.set("...", ...);
>     solver.solve(A, x, b);
> 
> The last call would call LU::solve() which calls
> PETScSparseMatrix::solve() which calls PETScLUSolver::solve().
> 
> /Anders
> 
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev




Follow ups

References