[petsc-users] TS question 1: how to stop explicit methods because they do not use SNES(VI)?

[Barry Smith]

> On Feb 14, 2017, at 3:18 PM, Emil Constantinescu <emconsta at mcs.anl.gov> wrote:
> 
> 
> 
> On 2/14/17 3:04 PM, Barry Smith wrote:
>> 
>>> On Feb 14, 2017, at 2:55 PM, Emil Constantinescu <emconsta at mcs.anl.gov> wrote:
>>> 
>>> On 2/14/17 2:33 PM, Zhang, Hong wrote:
>>>> I think many users (including me) would like to start with academic examples, e.g. u_t=f(u)+g(u), when they try to learn PETSc TS solvers. This simple form allows for easy switch between all kinds of different integration methods.
>>> 
>>> Right, but then you can just write an if statement along the lines in ex31.c: (this was before RHSJacobian and it's not intended for IMEX)
>>> 
>>> TSGetType(ts,&time_scheme);
>>> if ((!strcmp(time_scheme,TSEULER)) || (!strcmp(time_scheme,TSRK)) || (!strcmp(time_scheme,TSSSP))) {
>>> /* Explicit time-integration -> specify right-hand side function ydot = f(y) */
>>> TSSetRHSFunction(ts,NULL,RHSFunction,&ptype[0]);
>>> } else if ((!strcmp(time_scheme,TSTHETA)) ||...
>>>             (!strcmp(time_scheme,TSARKIMEX))) {
>>>    /* Implicit time-integration -> specify left-hand side function ydot-f(y) = 0 */
>>>    /* and its Jacobian function                   */
>>>    TSSetIFunction(ts,NULL,IFunction,&ptype[0]);
>>>    TSSetIJacobian(ts,Jac,Jac,IJacobian,&ptype[0]);
>>>  }
>> 
>>  Multiple yucks! Switching methods should not require compiling code except when absolutely necessary and this is not a case of absolutely necessary!
> 
> Wait, no. You don't recompile. When you set the -ts_type then it defines RHSFunction (for Explicit and Implicit) or RHS Function & IFunction (for IMEX).

   Ok, you don't recompile but forcing that into user code is still disgusting. With my api the user code is

>>> TSSetRHSFunction(ts,NULL,RHSFunction,&ptype[0]);
>>> TSSetLHSFunction(ts,NULL,LHSFunction,&ptype[0]);
>>>    TSSetRHSJacobian(ts,Jac,Jac,RHSJacobian,&ptype[0]);
>>>   TSSetLHSJacobian(ts,Jac,Jac,LHSJacobian,&ptype[0]);

and -ts_type xxx works correctly for ALL methods, implicit, explicit and imex without requiring any special command line options for different methods.


> Adding all that logic to keep track of left sides and right sides for academic examples is likely not the best development.

  I don't think it is "just academic examples", it is all examples without a mass matrix. 

   Once the user has decided with ts_type to use for production if it is fully implicit or explicit then they can depending on the type selected, write just a left hand side, just a right hand side for higher efficiency (less update of ghost points, fewer iterations over loops etc).

   With a constant mass matrix we can have TSSetMassMatrix() and then TSSetIFunction() is reserved for when it is absolutely needed.

  Barry

> 
> Emil
> 
>> 
>>> 
>>> This should not be a performance bottleneck. IFunction can be a wrapper of the RHSFunction so not too much extra coding.
>>> 
>>> Emil
>>> 
>>>> Experienced users or experts who need to solve DAEs or complicate problems involving nontrivial mass matrix should be encouraged to try IFunction/IJacobian. Of course, we need better documentation to help them realize the switching is not always possible. Personally it took me a long time to get used to the IFunction business since I had been using u_t=f(u)+g(u) for a couple years before I jumped on PETSc. If we add support for this simple form, the learning curve would be less steep than it currently is.
>>>> 
>>>> Hong (Mr.)
>>>> 
>>>>> On Feb 14, 2017, at 11:55 AM, Jed Brown <jed at jedbrown.org> wrote:
>>>>> 
>>>>> Barry Smith <bsmith at mcs.anl.gov> writes:
>>>>>> Hence my suggestion to have TSSetLeftHandSideFunction() (or Jed's suggestion to have multiple Right Hand side functions) this will allow comparison of implicit, explicit, imex without recompiling (which we don't have currently) for the case with no mass matrix. With a mass matrix, unless we use mass lumping and have a TSSetMassMatrix() you cannot switch to full explicit, but that is due to mathematics, not the interface.
>>>>> 
>>>>> Who in reality has a problem that looks like
>>>>> 
>>>>> u_t = f(u) + g(u) + h(u)
>>>>> 
>>>>> where it may make sense to move some of this to the left (implicit)?
>>>>> I'm concerned that this isn't natural for many applications so adding an
>>>>> interface would be solving a problem that doesn't really exist outside
>>>>> perhaps a few academic examples.
>>