I think I'm going to go ahead and start on a Preconditioner base
class... that doesn't inherit from LinearSolver for now. I'll get it
up and running with a PetscPreconditioner class.
This will pretty much meet my requirements for now.
Derek
On Jan 23, 2009, at 2:19 PM, Derek Gaston wrote:
> Ok - all good arguments... and I'm seeing some of your reasoning
> now.... but let me show you the use case I'm thinking of which might
> better explain why I think the LinearSolver link is useful...
>
> I'm needing to solve individual systems _just_ for preconditioning
> purposes... and they use a preconditioner _as a solver_. If a
> Preconditioner was _not_ a LinearSolver this would be what I would
> have to do (and is close to what I'm doing now... although there is
> even more manual manipulation).
>
> ApplyPreconditioner(Vector x, Vector & y)
> {
> LinearImplicitSystem sys;
> sys.add_variable("u");
>
> AssemblePrecMatrix(sys.matrix);
>
> copyCorrectPiecesOfX(x, sys.rhs);
>
> Preconditioner pre(sys.matrix);
> pre.apply(sys.rhs, sys.solution);
>
> copyCorrectPiecesOfSolution(sys.solution, y);
> }
>
> This really doesn't look too bad... the problem is that there's
> datastructures that aren't used... or are redundant. For instance
> there's a KSP object inside of that LinearImplicitSystem that never
> gets used... and it has it's own PC object... again that never gets
> used. Then there's the fact that you have to call apply() and pass
> in two things that are both held within the system.
>
> If you could just create the LinearSystem, passing in the
> Preconditioner then just call sys.solve().... it would simplify
> things. I suppose that there's nothing stopping us from adding a
> precondition() method to LinearSystem.... that just calls the
> Preconditioner object (that you can pass in at construction time or
> attach).
>
> What do you guys think?
>
> Derek
>
> On Jan 23, 2009, at 12:18 PM, John Peterson wrote:
>
>> On Fri, Jan 23, 2009 at 12:01 PM, Derek Gaston <[email protected]>
>> wrote:
>>> Alright! Now we have a debate! Let me see if I can persuade
>>> you... and if
>>> I can't... I'm probably going to code it up anyway to give you a
>>> chance to
>>> see what I mean.
>>>
>>> First, let's look at the two different forms of preconditioning
>>> (I'll be
>>> using GMRES the whole time because it's easy to think about)...
>>>
>>> Solving: Ax=b
>>>
>>> Left Preconditioning: Find P^-1 Ax = P^-1b
>>> Right Preconditioning: Find A P^-1 P x = b
>>
>> You're right, I wasn't taking right-preconditioning into account in
>> my
>> previous email.
>>
>>
>>> Further, from Axelsson: Iterative Solution Methods (2000) pg. 472
>>> when
>>> talking about left preconditioning GMRES:
>>>
>>> "The statement h:=C^-1 r is to be interpreted as solving the
>>> system Ch = r."
>>>
>>> Where h = y and C = P.
>>
>> What you've perfectly described is a Preconditioner object that
>> *uses*
>> a LinearSolver during its application, not something which "is a"
>> LinearSolver. Of course, this doesn't preclude you from also
>> having a
>> Preconditioner::solve() member, which just applies the
>> preconditioner.
>> It doesn't, in my opinion, make it a LinearSolver, in the sense in
>> which we use that name in Libmesh.
>>
>>
>>> You say that a Preconditioner doesn't satisfy the "is a"
>>> relationship for
>>> inheritance... but I hope I've shown you that in fact it does. A
>>> Preconditioner "is a" LinearSolver that solves a specific Py=q
>>> system... the
>>> solution of which is used internally by another LinearSolver.
>>
>> Perhaps, but in LibMesh a LinearSolver isn't a just thing that
>> "solves
>> a specific Py=q system," although it does do that. It's an abstract
>> interface to 3rd-party concrete linear algebra libraries as I
>> mentioned before. So, it's a little unfortunate that we have already
>> taken the name that you would like to use but I still argue strongly
>> against adding your proposed Preconditioner to this hierarchy. We
>> might need to have a discussion about expanding the meaning of
>> LinearSolver in libmesh, but I think the current system is elegant
>> and
>> relatively simple to understand for new users.
>>
>> --
>> John
>
------------------------------------------------------------------------------
This SF.net email is sponsored by:
SourcForge Community
SourceForge wants to tell your story.
http://p.sf.net/sfu/sf-spreadtheword
_______________________________________________
Libmesh-devel mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/libmesh-devel
