dolfin team mailing list archive

[Martin Sandve Alnæs <martinal@xxxxxxxxx>]

On 15 June 2011 13:06, Anders Logg <logg@xxxxxxxxx> wrote:
> Would this work? What it does is to define a simple Equation class
> that holds a pair of forms (lhs, rhs) but can still be used for
> comparisons.

Short story:
I'm pretty sure this will work out fine, except that the
Equation.__eq__ function must be dropped or changed
as explained below.


The details:

> class Equation:
>
>    def __init__(self, lhs, rhs):
>        self.lhs = lhs
>        self.rhs = rhs


Allowing the syntax (a == L) == True has severe pitfalls,
and I see no use for it. It breaks the Python model.
In general, <expression> == True is not good python code,
as many boolean expressions are allowed to return non-bool
objects. That's why this works out in the first place, so we
should respect that :)

Anyway, we only need bool(a == L) to say "if a == L" and use
a and L in dicts and sets. So this would have to change:

(class Equation)
>    def __eq__(self, other):
>        if isinstance(other, bool):
>            return self.__nonzero__() == other
>        else:
>            return repr(self) == repr(other)

It can eventually be replaced with this:

   def __eq__(self, other):
       return isinstance(other, Equation) and\
           bool(self.lhs == other.lhs) and\
           bool(self.rhs == other.rhs)


As a generic implementation pattern, this should rather
use comparison functions from lhs and rhs:

>    def __nonzero__(self):
>        return repr(self.lhs) == repr(self.rhs)
    def __nonzero__(self):
        return self.lhs.eval_equals(self.rhs) # or another name, maybe
__cmp__ == 0

But for this particular case, it's probably ok since lhs and rhs are always
Form objects, and repr is used in that single __eq__ implementation today.

Note that I have only evaluated this for Form/Equation, my conclusion
that I guess this is ok does not extend to any other part of the universe...
E.g., if we were to use this trick for Expr as well, it could have severe
performance penalties, so that's not something to do in a hurry.

Martin


> class Form:
>
>    def __init__(self, a):
>        self.a = a
>
>    def __eq__(self, other):
>        return Equation(self, other)
>
>    def __repr__(self):
>        return str(self.a)
>
> f1 = Form(1)
> f2 = Form(2)
> f3 = Form(1)
>
> eq = f1 == f2
> print eq == False
> print eq == True
>
> if f1 == f2:
>    print "f1 == f2 is True"
>
> if f1 == f3:
>    print "f1 == f3 is True"
>
> --
> Anders
>
>
> On Tue, Jun 14, 2011 at 10:01:47PM +0200, Martin Sandve Alnæs wrote:
>> I agree it is pretty in a way, but here comes the showstopper... Sorry! ;)
>>
>> There can be no ufl.Form.__eq__ implementation
>> that does not conform to the Python conventions
>> of actually comparing objects and returning
>> True or False, because that will break usage
>> of Form objects in built in Python data structures.
>>
>> Martin
>>
>> On 14 June 2011 21:53, Anders Logg <logg@xxxxxxxxx> wrote:
>> > On Tue, Jun 14, 2011 at 09:03:31PM +0200, Marie E. Rognes wrote:
>> >> On 06/14/2011 08:35 PM, Garth N. Wells wrote:
>> >> >
>> >> >
>> >> >On 14/06/11 19:24, Anders Logg wrote:
>> >> >>On Tue, Jun 14, 2011 at 10:19:20AM -0700, Johan Hake wrote:
>> >> >>>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 ;)
>> >> >>
>> >> >>Yes, a pretty syntax has been a priority for me ever since we
>> >> >>started. I think it is worth a lot.
>> >> >>
>> >> >
>> >> >Magic and pretty are not the same thing.
>> >
>> > That's true, but some magic is usually required to make pretty.
>> >
>> > Being able to write dot(grad(u), grad(v))*dx is also a bit magic.
>> > The step from there to solve(a == L) is short.
>> >
>> >> >>>The problem with this syntax is that who on earth would expect a
>> >> >>>VariationalProblem to be the result of an == operator...
>> >> >>
>> >> >>I don't think that's an issue. Figuring out how to solve variational
>> >> >>problems is not something one picks up by reading the Programmer's
>> >> >>Reference. It's something that will be stated on the first page of any
>> >> >>FEniCS tutorial or user manual.
>> >> >>
>> >> >>I think the solve(a == L) is the one missing piece to make the form
>> >> >>language complete. We have all the nice syntax for expressing forms in
>> >> >>a declarative way, but then it ends with
>> >> >>
>> >> >>problem = VariationalProblem(a, L)
>> >> >>problem.solve()
>> >> >>
>> >> >>which I think looks ugly. It's not as extreme as this example taken
>> >> >>from cppunit, but it follows the same "create object, call method on
>> >> >>object" paradigm which I think is ugly:
>> >> >>
>> >> >>   TestResult result;
>> >> >>   TestResultCollector collected_results;
>> >> >>   result.addListener(&collected_results);
>> >> >>   TestRunner runner;
>> >> >>   runner.addTest(CppUnit::TestFactoryRegistry::getRegistry().makeTest());
>> >> >>   runner.run(result);
>> >> >>   CompilerOutputter outputter(&collected_results, std::cerr);
>> >> >>   outputter.write ();
>> >> >>
>> >> >>>I see the distinction between FEniCS developers who have programming versus
>> >> >>>math in mind when designing the api ;)
>> >> >>
>> >> >>It's always been one of the top priorities in our work on FEniCS to
>> >> >>build an API with the highest possible level of mathematical
>> >> >>expressiveness to the API. That sometimes leads to challenges, like
>> >> >>needing to develop a special form language, form compilers, JIT
>> >> >>compilation, the Expression class etc, but that's the sort of thing
>> >> >>we're pretty good at and one of the main selling points of FEniCS.
>> >> >>
>> >> >
>> >> >This is an exaggeration to me. The code
>> >> >
>> >> >   problem = [Linear]VariationalProblem(a, L)
>> >> >   u = problem.solve()
>> >> >
>> >> >is compact and explicit. It's a stretch to call it ugly.
>> >
>> > Yes, of course it's a stretch. It's not very ugly, but enough to
>> > bother me.
>> >
>> >> >>>Also __eq__ is already used in ufl.Form to compare two forms.
>> >> >>
>> >> >>I think it would be worth replacing the use of form0 == form1 by
>> >> >>repr(form0) == repr(form1) in UFL to be able to use __eq__ for this:
>> >> >>
>> >> >>class Equation:
>> >> >>   def __init__(self, lhs, rhs):
>> >> >>       self.lhs = lhs
>> >> >>       self.rhs = rhs
>> >> >>
>> >> >>class Form:
>> >> >>
>> >> >>   def __eq__(self, other):
>> >> >>       return Equation(self, other)
>> >> >>
>> >> >>I understand there are other priorities, and others don't care as much
>> >> >>as I do about how fancy we can make the DOLFIN Python and C++ interface,
>> >> >>but I think this would make a nice final touch to the interface.
>> >> >>
>> >> >
>> >> >I don't see value in it. In fact the opposite - it introduces complexity
>> >> >and a degree of ambiguity.
>> >
>> > Complexity yes (but not much, it would require say around 50-100
>> > additional lines of code that I will gladly contribute), but I don't
>> > think it's ambiguous. We could perform very rigorous and helpful
>> > checks on the input arguments.
>> >
>> >> Evidently, we all see things differently. I fully support Anders in
>> >> that mathematical expressiveness is one of the key features of
>> >> FEniCS, and I think that without pushing these types of boundaries
>> >> with regard to the language, it will end up as yet another finite
>> >> element library.
>> >>
>> >> Could we compromise on having the two versions, one explicit (based
>> >> on LinearVariational[Problem|Solver] or something of the kind) and
>> >> one terse (based on solve(x == y)) ?
>> >
>> > That's what I'm suggesting. The solve(x == y) would just rely on the
>> > more "explicit" version and do
>> >
>> >  <lots of checks>
>> >  LinearVariationalSolver solver(x, y, ...);
>> >  solver.solve(u);
>> >
>> > So in essence what I'm asking for is please let me add that tiny layer
>> > on top of what we already have + remove the Problem classes (replaced
>> > by the Solver classes).
>> >
>> >
>> > _______________________________________________
>> > Mailing list: https://launchpad.net/~dolfin
>> > Post to     : dolfin@xxxxxxxxxxxxxxxxxxx
>> > Unsubscribe : https://launchpad.net/~dolfin
>> > More help   : https://help.launchpad.net/ListHelp
>> >
>