dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

On Mon, Jun 13, 2011 at 10:25:27PM +0200, Marie E. Rognes wrote:
> On 06/13/2011 09:04 PM, Anders Logg wrote:
> >On Mon, Jun 13, 2011 at 09:24:43AM -0700, Johan Hake wrote:
> >>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<logg@xxxxxxxxx>  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!
> >
> >Further comments? What does Marie say?
> >
>
> Marie thinks this discussion went from something fairly fundamental to
> something rather extreme quite quickly ;-)
>
> I absolutely see that [Linear|Nonlinear]Problem is more
> explicit. However, I still think that just splitting
> VariationalProblem into LinearProblem and NonlinearProblem is a change
> that implies a significant amount of work (which I am sure the editors
> of a certain book are very aware of) with little gain in terms of
> functionality/simplicity.

I don't think it involves that much work. For linear problems, it will
just be a matter of replacing Variational --> Linear, but see below.

> But while we are heading in that direction, how about abolishing the
> *Problem class(es) altogether, and just use LinearVariationalSolver
> and NonlinearVariationalSolver/NewtonSolver taking as input (a, L, bc)
> and (F, dF, bcs), respectively.

This will be in line with an old blueprint. We noted some time ago
that problems/solvers are designed differently for linear systems
Ax = b than for variational problems a(u, v) = L(v). For linear
systems, we have solvers while for variational problems we have both
problem and solver classes.

> I mean, the main difference lies in
> how to solve the problems, right?

It looks like the only property a VariationalProblem has in addition
to (forms, bc) + solver parameters is the parameter
symmetric=true/false.

If we go this route, we could mimic the design of the linear algebra
solvers and provide two different options, one that offers more
control, solver = KrylovSolver() + solver.solve(), and one quick
option that just calls solve:

1. complex option

solver = LinearVariationalSolver() # which arguments to constructor?
solver.parameters["foo"] = ...
u = solver.solve()

2. simple option

u = solve(a, L, bc)

3. very tempting option (simple to implement in both C++ and Python)

u = solve(a == L, bc)    # linear
u = solve(F == 0, J, bc) # nonlinear

--
Anders