dolfin team mailing list archive

[Ola Skavhaug <skavhaug@xxxxxxxxx>]

Anders Logg skrev den 01/04-2008 følgende:
> On Tue, Apr 01, 2008 at 11:23:29AM +0200, Martin Sandve Alnæs wrote:
> > 2008/4/1, Kent-Andre Mardal <kent-and@xxxxxxxxx>:
> > >
> > >  tir, 01.04.2008 kl. 11.08 +0200, skrev Martin Sandve Alnæs:
> > >
> > > > I don't see what's the point of that. The only reasons for
> > >  > Matrix/Vector are related to having a single LA backend in each
> > >  > application run, but a variety of solvers may be used in the same
> > >  > application.
> > >  >
> > >
> > >
> > > I thought the reason for having Matrix/Vector was to ensure that only
> > >  GenericMatrix/GenericVector functionality was used. I guess the PETSc
> > >  and the uBlas families now live seperate lives.
> > 
> > Then why don't you come up with an example of where this makes sense
> > for solvers.
> 
> Here's an example:
> 
> Choosing the type of solver based on some criterion and then
> reusing the solver:
> 
>    // At startup
>    LinearSolver* solver = 0;
>    if (solver_type == "direct")
>      solver = new LUSolver();
>    else
>      solver = new KrylovSolver()
> 
>    // Then reuse it
>    while (t < T)
>    {
>      ...
>      solver->solve(A, x, b);
>    }
> 
> This is used in the ODE solvers, where one may specify the type of
> solver used to solve the linearized equations in each iteration based
> on a parameter.
> 
> When calling LUSolver::solve(), it needs to choose a specific backend
> for the solution, depending on the backend used for A, x, b. Same for
> KrylovSolver.
> 
> Still not sure if it should be LinearSolver or GenericLinearSolver
> above.

What about a solution something like the following. It only uses the top level
generic types during solve:

GenericSolver* solver = 0;

if (solver_type == "direct")
    solver = A.factory()->directSolver();
else
    solver = A.factory()->iterativeSolver();

GenericPreconditioner* B = 0;
if (preconditioner_type == "multilevel") {
   GenericPreconditioner* B = A.factory()->multilevelPreconditioner(); 
   B->update(A); // Or
   //GenericPreconditioner* B = A.factory()->multilevelPreconditioner(A);
}

while (t < T)
{
    if (B != 0)
        solver->solve(B, A, x, b);
    else
        solver->solve(A, x, b);
}

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