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. > > 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. -- 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
