dolfin team mailing list archive

[Johan Hake <johan.hake@xxxxxxxxx>]

On Tuesday June 14 2011 03:33:59 Anders Logg wrote:
> On Tue, Jun 14, 2011 at 09:25:17AM +0100, Garth N. Wells wrote:
> > On 14/06/11 08:53, Anders Logg wrote:
> > > 14 jun 2011 kl. 09:18 skrev "Garth N. Wells" <gnw20@xxxxxxxxx>:
> > >> On 14/06/11 08:03, Marie E. Rognes wrote:
> > >>> On 06/13/2011 11:16 PM, Anders Logg wrote:
> > >>>>>> 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()
> > >> 
> > >> I favour this option, but I think that the name
> > >> 'LinearVariationalSolver' is misleading since it's not a
> > >> 'variational solver', but solves variational problems, nor should
> > >> it be confused with a LinearSolver that solves Ax = f.
> > >> LinearVariationalProblem is a better name. For total control, we
> > >> could have a LinearVariationalProblem constructor that accepts a
> > >> GenericLinearSolver as an argument (as the NewtonSolver does).
> > >> 
> > >>> For the eigensolvers, all arguments go in the call to solve.
> > >>> 
> > >>>> 2. simple option
> > >>>> 
> > >>>> u = solve(a, L, bc)
> > >> 
> > >> I think that saving one line of code and making the code less
> > >> explicit isn't worthwhile. I can foresee users trying to solve
> > >> nonlinear problems with this.
> > > 
> > > With the syntax suggested below it would be easy to check for errors.
> > > 
> > >>> Just for linears?
> > >>> 
> > >>>> 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
> > >> 
> > >> I don't like this on the same grounds that I don't like the
> > >> present design. Also, I don't follow the above syntax
> > > 
> > > I'm not surprised you don't like it. But don't understand why. It's
> > > very clear which is linear and which is nonlinear. And it would be
> > > easy to check for errors. And it would just be a thin layer on top of
> > > the very explicit linear/nonlinear solver classes. And it would
> > > follow the exact same design as for la with solver classes plus a
> > > quick access solve function.
> > 
> > Is not clear to me - possibly because, as I wrote above, I don't
> > understand the syntax. What does the '==' mean?
> 
> Here's how I see it:
> 
> 1. Linear problems
> 
>   solve(a == L, bc)
> 
>   solve the linear variational problem a = L subject to bc
> 
> 2. Nonlinear problems
> 
>   solve(F == 0, bc)
> 
>   solve the nonlinear variational problem F = 0 subject to bc
> 
> It would be easy to in the first case check that the first operand (a)
> is a bilinear form and the second (L) is a linear form.
> 
> And it would be easy to check in the second case that the first
> operand (F) is a linear form and the second is an integer that must be
> zero.
> 
> In both cases one can print an informative error message and catch any
> pitfalls.
> 
> The nonlinear case would in C++ accept an additional argument J for
> the Jacobian (and in Python an optional additional argument):
> 
>   solve(F == 0, J, bc);
> 
> The comparison operator == would for a == L return an object of class
> LinearVariationalProblem and in the second case
> NonlinearVariationalProblem. These two would just be simple classes
> holding shared pointers to the forms. Then we can overload solve() to
> take either of the two and pass the call on to either
> LinearVariationalSolver or NonlinearVariationalSolver.
> 
> I'm starting to think this would be an ideal solution. It's compact,
> fairly intuitive, and it's possible to catch errors.
> 
> The only problem I see is overloading operator== in Python if that
> has implications for UFL that Martin objects to... :-)

Wow, you really like magical syntaxes ;) 

The problem with this syntax is that who on earth would expect a 
VariationalProblem to be the result of an == operator... 

I see the distinction between FEniCS developers who have programming versus 
math in mind when designing the api ;) 

Also __eq__ is already used in ufl.Form to compare two forms.

Johan

> 
> --
> Anders
> 
> _______________________________________________
> Mailing list: https://launchpad.net/~dolfin
> Post to     : dolfin@xxxxxxxxxxxxxxxxxxx
> Unsubscribe : https://launchpad.net/~dolfin
> More help   : https://help.launchpad.net/ListHelp