> On Feb 14, 2017, at 3:18 PM, Emil Constantinescu <[email protected]> wrote: > > > > On 2/14/17 3:04 PM, Barry Smith wrote: >> >>> On Feb 14, 2017, at 2:55 PM, Emil Constantinescu <[email protected]> >>> 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 <[email protected]> wrote: >>>>> >>>>> Barry Smith <[email protected]> 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. >>
