From: Derek G. <fri...@gm...> - 2009-01-20 18:50:23
|
I know there was a lot of discussion that led up to adding the ShellMatrix object to libMesh.... and looking back through it I see that preconditioning was mentioned... but now that I'm looking at the code I don't see a straight forward way to add a ShellMatrix as a preconditioner.... but maybe I'm missing it. What I mean is... there needs to be a call to PCSetType(pc,PCSHELL) (which there is) and a call to PCShellSetApply.... which there isn't... or am I just missing it? Currently I'm just using my own hacked up code that doesn't even use ShellMatrix... but it would be nice to use the interfaces (for easier transition to Trilinos). Also... on a related question.... what do you guys think about providing a Preconditioner base class in libMesh? Concrete implementations of this class would be HyprePreconditioner, MLPreconditioner, LinearSolvePreconditioner, etc... I guess they are similar to linear solvers... but not exactly the same.... Derek |
From: John P. <jwp...@gm...> - 2009-01-20 18:58:57
|
On Tue, Jan 20, 2009 at 12:49 PM, Derek Gaston <fri...@gm...> wrote: > I know there was a lot of discussion that led up to adding the > ShellMatrix object to libMesh.... and looking back through it I see > that preconditioning was mentioned... but now that I'm looking at the > code I don't see a straight forward way to add a ShellMatrix as a > preconditioner.... but maybe I'm missing it. All I was able to finish was checking in Tim's patch, I didn't have a chance to redo the whole NumericMatrix hierarchy yet. It's on my todo list but there are day-job things taking up time currently. > What I mean is... there needs to be a call to PCSetType(pc,PCSHELL) > (which there is) and a call to PCShellSetApply.... which there > isn't... or am I just missing it? > > Currently I'm just using my own hacked up code that doesn't even use > ShellMatrix... but it would be nice to use the interfaces (for easier > transition to Trilinos). > > Also... on a related question.... what do you guys think about > providing a Preconditioner base class in libMesh? Concrete > implementations of this class would be HyprePreconditioner, > MLPreconditioner, LinearSolvePreconditioner, etc... I guess they are > similar to linear solvers... but not exactly the same.... I guess I don't see the need for this, but that's not to say I can't be convinced. There's nothing like that in PETSc AFAIK. -- John |
From: Jed B. <je...@59...> - 2009-01-20 19:35:36
|
On Tue, Jan 20, 2009 at 09:58, John Peterson <jwp...@gm...> wrote: >> Also... on a related question.... what do you guys think about >> providing a Preconditioner base class in libMesh? Concrete >> implementations of this class would be HyprePreconditioner, >> MLPreconditioner, LinearSolvePreconditioner, etc... I guess they are >> similar to linear solvers... but not exactly the same.... > > I guess I don't see the need for this, but that's not to say I can't > be convinced. There's nothing like that in PETSc AFAIK. There sure is, but the linear solve is actually a method for Mat. In OO language, all PETSc objects come from a factory which can be influenced by command line parameters. When you use -pc_type hypre, you make the factory create an instance of the derived class PC_HYPRE. At least that is how a similar API would be implemented in C++ so that the user never hard-codes concrete types. In the current interface, the factories are not bona-fide objects, but it's a reasonable way to think of PETSc objects before XXXSetFromOptions or whatever is called. The trouble with implementing this on the Libmesh side is that you duplicate a lot of work and lose significant functionality unless you reimpliment a lot of PETSc. PCSHELL is always there for when you have a clever way to form a preconditioner. I don't know if it is useful for anyone on this list, but I like to think of forming a preconditioner as PCApply = PC(A,Ap) where A is an implementation of MatMult which is frequently ignored and Ap is ``extra information'' which is usually an assembled sparse matrix, but could be any information allowing you to create a (nominally) linear operation PCApply so that PCApply \circ A or A \circ PCApply has favorable properties for the Krylov solver being used ($\circ$ is function composition). This view makes it clear that the preconditioner is intimately coupled to the format of Ap. PETSc handles this by implementing parts of many preconditioners as methods of Mat (factor, solve, relax, etc). Many third party preconditioners require their own matrix type, usually derived from AIJ. Jed |
From: Derek G. <fri...@gm...> - 2009-01-20 19:37:49
|
On Jan 20, 2009, at 11:58 AM, John Peterson wrote: > All I was able to finish was checking in Tim's patch, I didn't have a > chance to redo the whole NumericMatrix hierarchy yet. It's on my todo > list but there are day-job things taking up time currently. Do you want me to stab an interface into LinearImplicitSystem? I would just basically follow what's there for shell matrices... with an attach_shell_preconditioner() function and all the rest. Or should we wait until the whole hierarchy gets redone? > I guess I don't see the need for this, but that's not to say I can't > be convinced. There's nothing like that in PETSc AFAIK. Right... part of the problem here is that I'm trying to write specialized preconditioner code.... but do it in a way that's agnostic of Petsc/Trilinos. In particular I need to be able to partially solve linear systems using multigrid solvers / preconditioners... it would be nice just to be able to instantiate a MultiGridPreconditioner object... having the correct one chose depending on how libMesh was configured and allowing you to set a lot of multigrid only options (for instance smoothing, etc.). Right now it's looking like I'm basically going to have two completely different code paths in my code that do the same thing (but call different solvers in the end).... and I can definitely just add abstractions into my own code for some of this stuff... but some of those abstractions might be useful to others in the community. I honestly don't know right now. Just thought I would get some opinions on this. I suppose I could just add a MultigridLinearSolver object... and concrete implementations of it... but that kind of goes against what we have setup with Petsc where there are a lot of actual solvers to choose from underneath (including multigrid options). The trouble is how to tell a specific PetscLinearSolver object to solve using Hypre... in a completely solver library agnostic way. Currently, when you are just solving one system it's ok to set options on the command line (like -pc_type hypre)... but that doesn't work for the stuff we're doing now... where we have _lots_ of small linear systems and some (not so small) nonlinear systems all being solved simultaneously and all preconditioned / solved differently. Currently we have to get ahold of the KSP object for each system and manually twiddle the options with Petsc calls.... but again... that means we have to do the same with Trilinos........ and that's a lot of duplicated logic. Anyway... at this point I'm rambling. If anyone has ideas, let me know. Derek |
From: Roy S. <roy...@ic...> - 2009-01-20 19:44:20
|
On Tue, 20 Jan 2009, Derek Gaston wrote: > I suppose I could just add a MultigridLinearSolver object... and > concrete implementations of it... but that kind of goes against what > we have setup with Petsc where there are a lot of actual solvers to > choose from underneath (including multigrid options). I'm not yet certain if the Preconditioner classes you mentioned earlier are the right way to solve this, but I agree 100% (and for the same reasons) that adding MultigridLinearSolver (+ ILULinearSolver, ...) classes would be the wrong way. --- Roy |
From: Jed B. <je...@59...> - 2009-01-20 19:57:17
|
On Tue, Jan 20, 2009 at 10:36, Derek Gaston <fri...@gm...> wrote: > Currently we have to get > ahold of the KSP object for each system and manually twiddle the > options with Petsc calls Set a prefix any time you have more than 1 KSP or SNES object. Then you can put this twiddling in a config file (see -options_file) and still override it from the command line. Presumably something similar is possible in Trilinos, but I think hard-coding this stuff is really bad unless your code is actively making solver decisions based on how the analysis is going. Jed |