TS_RK enhancements

[Kai Germaschewski]

On Thu, May 14, 2009 at 3:36 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:

  There is no one who is coordinating the TS interface at the moment. So you
> are welcome to propose changes/fixes/improvements and are in as good as a
> position as any of us.
>
>   I see there are four methods: prestep, poststep, update, postupdate; as
> you noticed the naming and usage of them is a bit inconsistent. Even the
> TSStep() function has a confusing name since it actually takes (potentially)
> many steps.  I suggest starting the following possible changes.
>
>   change TSStep to TSSolve(),


Yup, that makes a lot of sense. Actually, there is a TSSolve() already,
which however seems rarely used (all of the examples use TSStep())

  move the prestep() and poststep() methods to be the ones used AT the
> beginning and end of each timestep (not called in TSSolve())


Agreed. I don't really like changing the behavior of functions while keeping
the name the same, but I'll just assume there are few users of the current
PreStep / PostStep, since they don't really make much sense.

>
>   eliminate update() and postupdate() confusing names anyways
>   determine if there is a need for a pre and post called in TSSolve()
> could call them presolve and postsolve if you put them in.


I don't think there's a need for a user-settable presolve/postsolve, because
if they want to call functions before / after calling TSSolve(), they can
just call them directly. There may well be good reasons to introduce methods
like this into TS that can be set by a specific implementation (like
euler/beuler/...) to simplify the code there.


> do we really want prestep and poststep deciding to kick out of the
> iteration or do we want yet another method (corresponding to the converged()
> method in KSP and SNES) that is used to exit early from the solver?
>   I think using a flag to exit early from the solver is better than setting
> dt to zero.


A new method seems like a much cleaner way, and would avoid breaking the
existing method by changing the call sequence. Now the next question is for
a good name. In some sense, "SetConvergenceTest" is nice because it's the
same as used elsewhere, but it really doesn't quite describe what it really
does here. So how about "SetTerminationTest"?


>
>   having the steps and ptime output arguments from TSSolve() I think is
> wrong. They should be kept in the TS object and then calls like
> TSGetSolveTime() and TSGetNumberSteps()
>   to access them.


I'm not sure what you're referring to. For TSMonitor(), I don't think
there's a good reason to change the calling sequence. For TSStep() ->
TSSolve(), I definitely agree, TSSolve() shouldn't return those values. I
think the TSSolve(TS ts, Vec sol) interface as it already exists is just
fine.

There's already TSGetTime() and TSGetTimeStepNumber(), which can be used to
do the above. They don't have the most elegant names, but I think we both
believe that one shouldn't have multiple names for the same functionality,
and I personally don't like renaming functions just for aesthetics because
it's a pain to change applications and destroy backwards compatibility. (but
past experience makes me think that the petsc developers look at this
differently... ;)

My transition plan looks like this:
- Integrate the current TSStep() into TSSolve() and remove TSStep from the
public interface (which requires
changing current users of TSStep()
- Rename the TSStep_* methods -> TSSolve_*
- Add a generic TSSolve_TSStep() function, that provides a default timeloop
implementation so that simple timestepping methods only have to implement a
new TSStep_* method which does one (complete) step, but the actual looping,
bookkeeping etc is handled by generic code. This could end up being an
interface much like what Lisandro has to python. On the other hand, more
complicated time integrators could take full control of the entire
integration. (If it turns out there are no users left for such an interface,
there may be no point in keeping it)

So actually I came across one point I'm wondering:
All current time steppers seem to do

VecCopy(sol, update);
while (!done) {
  FindSolutionAtNewTime(update);
  VecCopy(update,sol);
}

Obviously for methods where steps may fail and will be retried, it's
necessary to preserve the old solution, however
many schemes just call KSPSolve(.., update) or SNESSolve(..., update). In a
linear method, I suppose by definition the operator does not depend on the
current "sol". In a nonlinear problem, it will, however I suppose that
SNESSolve() wouldn't update the solution vector until the very end, or am I
missing something?

Point being, all these VecCopy()'s aren't really necessary and actually end
up making the code a bit more complicated than needs be, and obviously
aren't helping performance. (For implicit methods, those copies will
generally be negligible, but for the explicit ones (euler eg), that's not so
clear...)

--Kai

-- 
Kai Germaschewski
Assistant Professor, Dept of Physics / Space Science Center
University of New Hampshire, Durham, NH 03824
office: Morse Hall 245F
phone:  +1-603-862-2912
fax: +1-603-862-2771
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20090516/5790d80d/attachment.html>