On Monday June 13 2011 08:09:36 Anders Logg wrote: > On Mon, Jun 13, 2011 at 05:04:53PM +0200, Kristian Ølgaard wrote: > > On 13 June 2011 16:37, Anders Logg <l...@simula.no> wrote: > > > On Mon, Jun 13, 2011 at 03:30:22PM +0100, Garth N. Wells wrote: > > >> On 13/06/11 14:54, Anders Logg wrote: > > >> > On Mon, Jun 13, 2011 at 02:42:57PM +0100, Garth N. Wells wrote: > > >> >> On 13/06/11 14:10, Anders Logg wrote: > > >> >>> On Mon, Jun 13, 2011 at 01:09:20PM +0100, Garth N. Wells wrote: > > >> >>>> On 13/06/11 12:51, Marie E. Rognes wrote: > > >> >>>>> On 06/13/2011 10:03 AM, Garth N. Wells wrote: > > >> >>>>>> On 12/06/11 22:54, Anders Logg wrote: > > >> >>>>>>> On Wed, Jun 08, 2011 at 11:47:10PM +0200, Anders Logg wrote: > > >> >>>>>>>> I'm on a cell phone and can't engage in any discussions but I > > >> >>>>>>>> just want to throw something in. Marie and I discussed at > > >> >>>>>>>> some point another option which is to write everything as F > > >> >>>>>>>> = 0 and let UFL compute J even for linear problems. J would > > >> >>>>>>>> be a optional variable. Then we could have a common > > >> >>>>>>>> interface and also a common solver. > > >> >>>>>> > > >> >>>>>> I think that this would be a nice option, but I don't think > > >> >>>>>> that it can work without a more elaborate interface than just > > >> >>>>>> > > >> >>>>>> pde = VariationalProblem(F, bcs) > > >> >>>>>> > > >> >>>>>> because UFL cannot know which coefficient in a form is unknown, > > >> >>>>>> and therefore which function to compute the linearisation with > > >> >>>>>> respect to. Moreover, > > >> >>>>>> > > >> >>>>>> - UFL cannot compute the linearisation for all forms > > >> >>>>>> - In some cases it's not desirable to let UFL do the > > >> >>>>>> linearisastion > > >> >>>>> > > >> >>>>> Some specific responses: First, the moment solve is called (with > > >> >>>>> a coefficient), which function to compute the linearization > > >> >>>>> with respect to, is known. > > >> >>>> > > >> >>>> Which means that the simple call > > >> >>>> > > >> >>>> u = pde.solve() > > >> >>>> > > >> >>>> will not be possible. > > >> >>> > > >> >>> Why do you call it "pde" and not "problem"? > > >> >> > > >> >> The name has nothing to do with the discussion. > > >> > > > >> > This discussion is about just that: which names to use. > > >> > > >> Yes, the class name, but not what I dream up for an object: > > >> > > >> xxx = FooVariationalFoo(....) > > >> > > >> u = xxx.solve() > > >> > > >> >> I could call it xxx. > > >> >> > > >> >>> Once upon a time it was > > >> >>> called that, so yet another option would be > > >> >>> > > >> >>> LinearPDE > > >> >>> NonlinearPDE > > >> >>> > > >> >>> which is short enough. > > >> >>> > > >> >>> What speaks against this, and likely one of the reasons we > > >> >>> switched to VariationProblem is that the above looks weird when > > >> >>> used to define a projection which does not involve any > > >> >>> derivatives. > > >> >>> > > >> >>> And here's yet another option: > > >> >>> > > >> >>> LinearProblem > > >> >>> NonlinearProblem > > >> >> > > >> >> This looks ok, but it's a secondary issue for if we agree that we > > >> >> should have separate classes for linear and nonlinear equations > > >> >> (or have we agreed this?). > > >> > > > >> > Not yet. I think we should settle first one whether we want > > >> > > > >> > 1. a common class for (non)linear problems > > >> > 2. two classes > > >> > > >> Agree. I strongly support 2. > > >> > > >> What I oppose vehemently is distinguishing between linear and > > >> nonlinear by the order of the lhs/rhs arguments. As pointed out by > > >> Martin, the line > > >> > > >> VariationalProblem(C, D) > > >> > > >> says nothing to the reader about the type of problem, and it has led > > >> to confusion. > > >> > > >> The code > > >> > > >> LinearVariationalProblem(C, D) > > >> > > >> or > > >> > > >> NonlinearVariationalProblem(C, D) > > >> > > >> (class names subject to change) are clear, and the code can test for > > >> the arity of C and D and throw an error if the order is mixed up. > > > > > > I agree with this. It's important that we can check input arguments. > > > > > > What are the votes in favor of either of these two options? > > > > > > 1. One class VariationalProblem for all (non)linear problems. > > > > > > 2. Two classes LinearProblem and NonlinearProblem.
2 gets my vote. > > I vote for 2 with the arguments of readability/debugging as Martin and > > Garth have promoted, > > and like Garth I also don't think that one class should try to do > > everything. > > ok. Note the new choice of names LinearProblem/NonlinearProblem > (without Variational). > > This would also let us keep the VariationalProblem class and issue a > suitable error message that it has been replaced by those other > classes. Sounds good! Johan > -- > Anders > > > Kristian > > > > >> >>>>> Second, the Jacobian could definitely be provided as an > > >> >>>>> optional argument for those cases where it is not > > >> >>>>> feasible/desirable for UFL to do the linearization. See below > > >> >>>>> for longer, more general rant. > > >> >>>>> > > >> >>>>>>> Any further thoughts on this? It's important we get this > > >> >>>>>>> right. > > >> >>>>>>> > > >> >>>>>>> I can live with > > >> >>>>>>> > > >> >>>>>>> LinearVariationalProblem > > >> >>>>>>> NonlinearVariationalProblem > > >> >>>>>>> > > >> >>>>>>> if everyone else thinks that's the best option, > > >> >>>>> > > >> >>>>> I don't think it is the best option. > > >> >>>>> > > >> >>>>>>> but don't really like > > >> >>>>>>> it since it's very long. > > >> >>>>> > > >> >>>>> My reasons do not really relate to the length of the name, but > > >> >>>>> rather that I do not see why we should make the difference > > >> >>>>> between linear and nonlinear variational problems greater in > > >> >>>>> the _highest level interface_ of DOLFIN. I think this just > > >> >>>>> makes things more cumbersome (a) in terms of debugging and (b) > > >> >>>>> in terms of explaining how to solve variational problems to > > >> >>>>> FEniCS beginners. > > >> >>>> > > >> >>>> My experience is the complete opposite of this. > > >> >>> > > >> >>> I agree with Marie here. It's a bit weird when we do > > >> >>> > > >> >>> u = TrialFunction(V) > > >> >>> ... > > >> >>> u = Function(V) > > >> >>> problem.solve(u) # or u = problem.solve() > > >> >> > > >> >> I agree that this is not elegant, but I think that it's a secondary > > >> >> issue (but still a something to improve). > > >> > > > >> > I don't think this is an issue we can fix in small increments. We > > >> > need to find out what kind of interface we want, then implement it > > >> > and then not touch it ever again. > > >> > > >> The above is part of the Python interface. If the C++ interface is > > >> resolved, we can turn our attention to prettifying Python issues, like > > >> the above. > > >> > > >> >>> It's almost as confusing as the V=V thing we had for a while. > > >> >>> > > >> >>>>> My typical pedagogical use case looks something like this (in > > >> >>>>> Python): Say you have explained to someone how to solve > > >> >>>>> stationary Stokes with Taylor-Hood, and then want to extend to > > >> >>>>> Navier-Stokes. Then ideally all that should be necessary is to > > >> >>>>> add the nonlinear term to the form(s). But no. In addition, you > > >> >>>>> at least have to replace the TrialFunction with a Function > > >> >>>>> (which typically leads to a longer discussion on what a > > >> >>>>> TrialFunction represents and why that is used for linear > > >> >>>>> problems when the same ansatz is made for the two in most text > > >> >>>>> books), and then you have to call solve with the same Function > > >> >>>>> instead of just returning it (which is also on my top 3 list of > > >> >>>>> pitfalls). I don't see how adding > > >> >>>>> LinearVariationalProblem/NonlinearVariationalProblem to the list > > >> >>>>> makes it better. The same (but the other way) goes when > > >> >>>>> wanting to debug Navier-Stokes by reducing to Stokes. > > >> >>>>> > > >> >>>>> In essence, I've started treating all variational problems as > > >> >>>>> nonlinear. > > >> >>>> > > >> >>>> This might be fine for toy problems, but not for real problems. > > >> >>>> There is overhead in a nonlinear solve to check for convergence, > > >> >>>> which is not required if the eqn is linear. > > >> >>> > > >> >>> I don't see why this should be restricted to toy problems. > > >> >>> Defining the Jacobian is always an option and then it would be > > >> >>> just like today. > > >> >> > > >> >> Because if I'm solving a very large linear problem, I'm going to be > > >> >> annoyed if the code insists on unnecessarily assembling the > > >> >> residual. > > >> > > > >> > Aren't all "real world problems" nonlinear? ;-) > > >> > > >> Fortunately no. > > >> > > >> >>> Perhaps it would be possible for UFL to check that a problem is > > >> >>> linear (by differentiating twice and getting zero) and then we > > >> >>> could remove the overhead of the extra residual evaluation. > > >> >> > > >> >> Too much magic over simplicity for my liking. We also need to bear > > >> >> in mind that we're talking about the C++ interface too. I think a > > >> >> better starting point is the C++ interface. The Python interface > > >> >> will follow naturally, and sugar can be added later. > > >> > > > >> > I agree with that. Say for a moment that we choose to implement the > > >> > common (non)linear interface. Then what should the interface look > > >> > like in C++? > > >> > > > >> > 1. Poisson::BilinearForm a(V, V); > > >> > Poisson::LinearForm L(V); > > >> > > > >> > 2. Poisson::Residual F(V); > > >> > Poisson::Jacobian J(V, V); > > >> > > > >> > 3. Poisson::VariationalProblem F(V, V); > > >> > > > >> >> What if the equation is linear, but I want to solve it as a > > >> >> nonlinear equation? I do this all the time as a first step for > > >> >> testing when testing a nonlinear solver? I know we can always 'add > > >> >> another option' to handle various case, but it all add complexity. > > >> > > > >> > Isn't that an argument for having a common nonlinear solver for all > > >> > problems? > > >> > > >> No. It's an argument that we shouldn't try to make one class do > > >> everything, and that the code should be explicit. > > >> > > >> >> Another factor is that there is more to solving nonlinear equations > > >> >> than vanilla versions of Newton's method. It would nice in the > > >> >> future to have a NonlinearFoo that does more than what we now, and > > >> >> without increasing the code complexity for a straightforward > > >> >> linear eqn solve. > > >> > > > >> > Perhaps, but one could also argue that is yet another argument in > > >> > favor of a single nonlinear solver. There is some code duplication > > >> > already in the two linear/nonlinear solvers, setting options etc. > > >> > > >> Ditto what I wrote above. Trying to make one class do too much leads > > >> to complications. > > > > > > _______________________________________________ > > > Mailing list: https://launchpad.net/~dolfin > > > Post to : dolfin@lists.launchpad.net > > > Unsubscribe : https://launchpad.net/~dolfin > > > More help : https://help.launchpad.net/ListHelp > > _______________________________________________ > Mailing list: https://launchpad.net/~dolfin > Post to : dolfin@lists.launchpad.net > Unsubscribe : https://launchpad.net/~dolfin > More help : https://help.launchpad.net/ListHelp _______________________________________________ Mailing list: https://launchpad.net/~dolfin Post to : dolfin@lists.launchpad.net Unsubscribe : https://launchpad.net/~dolfin More help : https://help.launchpad.net/ListHelp