← Back to team overview

dolfin team mailing list archive

Re: [HG] Add uBlasKrylovSolver.cpp

 

On Tue, 2006-06-06 at 20:10 +0200, Anders Logg wrote:
> On Tue, Jun 06, 2006 at 04:41:31PM +0200, Garth N. Wells wrote:
> > 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.
> 
> But that seems to make the class LinearSolver unusable? Then you can't
> do
> 
> LinearSolver* solver = new LU();
> solver.solve(A, x, b);
> 
> since LU does not need to implement solve() with the matrix type of A.
> 

It does work, so long the function "solve" is not a pure virtual
function. I've just checked it in so try it out.


> > > 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.  
> 
> Yes, I think they should be separate. But is it possible without
> working putting in a small function in the matrix interface that
> chooses the correct solver (or by templating)?
> 

Sorry, I wasn't thinking too straight on this point before. It would be
nice and clean to keep the solvers completely separate from the matrix
classes though. 

Garth

> Another option is to overload on all possible matrix types in
> LinearSolver and LU, and choose the specific solver there.
> 
> /Anders
> 
> 
> > 
> > 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
> > 
> > 
> > _______________________________________________
> > DOLFIN-dev mailing list
> > DOLFIN-dev@xxxxxxxxxx
> > http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev
> 
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev




Follow ups

References