Re: [eigen] still the solve() API debate |
[ Thread Index |
Date Index
| More lists.tuxfamily.org/eigen Archives
]
- To: eigen@xxxxxxxxxxxxxxxxxxx
- Subject: Re: [eigen] still the solve() API debate
- From: Gael Guennebaud <gael.guennebaud@xxxxxxxxx>
- Date: Fri, 11 Sep 2009 10:25:23 +0200
- Dkim-signature: v=1; a=rsa-sha256; c=relaxed/relaxed; d=gmail.com; s=gamma; h=domainkey-signature:mime-version:received:in-reply-to:references :date:message-id:subject:from:to:content-type :content-transfer-encoding; bh=9YqQDQ9nErJsb0Utn8kkWGLBs0kyNorWE1qSCHeG2Jo=; b=pQAxCV4BCY9krpS/RKEVcvYhmWRHEYtz/A/0BDJBrhVWQl+5jFiMogp9o7btH/Tvxh D3AcCJVssY74SDVGRbKLSJnrNvP2AOF16R2Roq0FfUnUIOzlvuHXB77G42zsp7AkwrWm 9aSbthI+A8c7m46/fq1riBT0LoWd8/Oj+gVTI=
- Domainkey-signature: a=rsa-sha1; c=nofws; d=gmail.com; s=gamma; h=mime-version:in-reply-to:references:date:message-id:subject:from:to :content-type:content-transfer-encoding; b=JiksBFh5J1Km7rajfs5alr7xzQxqUeY032cuRJeMelGxNBE4hXgXZdVzHSEF0meV4g A59F8j7/1hSlulIcBE3cFlXr/3JDtRgzqWTBtpwI2ziavXKx4e3a2pgzPggqJA4Su06a rZMOe4Za8MB21rTQ9bAES73ym9t76vgV8PQ0o=
First of all I don't like the use of pointers here, because it is not
a standard c++ practice, and therefore users we'll start wondering: "
why do they use a pointer here ? maybe because this parameter is
optional ? hm does not make sense... unless if x==0, then the
computation is performed in place overwritten b ? or the solve
function creates the matrix object and then I have to delete it myself
? etc. ". OK, maybe I'm exaggerating a bit, but that's the idea.
ok, now let's summarize, ideally we want an API which satisfies the
following constraints:
a - concise (not too verbose)
b - readable
c - simple behavior
d - it should return the result so that we can directly use it in an expression
e - does not create unnecessary temporary
f - directly works with temporary proxy object
g - compliant with standard c++ practices
since I changed my mind about ReturnByValue inheriting MatrixBase, I
think opt 3 is the one which best satisfies all these constraints.
However, I (we) forgot to include in the discussion the API of the
in-place solve functions which should be as close as possible to the
generic solve(). Currently this only apply to the triangular solver,
LLT and LDLT, but I think this could also apply to PartialLU. Here the
problem is a bit different because the result is also and input and so
only option 1 and 4 could work:
opt 1: A.solveInPlace(&b_and_x);
opt 4: A.solveInPlace(b_and_x.ref());
Since opt 1 does not satisfy the constraint f and that this constraint
is particularly important for the in-place triangular solver, it
currently remains the opt 4 only.
In order to satisfy the constraint d, the solveInPlace functions could
return a reference to the inout argument. So at that point one could
do either:
Z = 2 * M * A.solveInPlace(b.ref()); // no temporary, and b is overwritten
Z = 2 * M * A.solve(b); // one temporary but b is preserved
So yes you got me right, I'm suggesting to rely on both the
ReturnByValue<> and Reference<> mechanisms ! and since they have not
been extensively tested yet I cannot promise that we won't discover
any issue with them (especially the Reference<> one).
Also, note that a variant of opt 4, is to replace .ref() by the (more
explicit ?) pair .out() and .inout(), in which case the "InPlace"
suffix may be removed, e.g.:
Z = 2 * M * A.solve(b.inout())
but note that how subtle the difference is with:
Z = 2 * M * A.solve(b);
so maybe it is better to stuck with .ref() and the InPlace suffix.
gael.
On Fri, Sep 11, 2009 at 1:10 AM, Benoit Jacob <jacob.benoit.1@xxxxxxxxx> wrote:
> 2009/9/10 Keir Mierle <mierle@xxxxxxxxx>:
>> On Thu, Sep 10, 2009 at 2:47 PM, Gael Guennebaud <gael.guennebaud@xxxxxxxxx>
>> wrote:
>>>
>>> oh, actually we can even write a single
>>> SolveReturnType<DecomositionType, DerivedB> class inheriting
>>> ReturnByValue which would be used for all x = solve(b) functions.
>>
>> I am not convinced the added complexity justifies the benefits. What about
>> using the result in an expression?
>> y = M*A.solve(b); // Parts of my code could use this today.
>
> My understanding is that Gael is now proposing a ReturnByValue class
> that inherits MatrixBase, so that one could use that in any
> expression. That's a departure from my assumption that he was
> proposing ReturnByValue not inheriting MatrixBase.
>
> Gael : is that right?
>
>> x.start(3) = A.solve(b) // Never used this.
>> Not saying that the 2nd use case isn't justified, just that I personally
>> haven't needed it. Can others speak to this?
>
> Same here, I don't think use cases are common. I can think of one:
> when implementing a blocked algorithm (e.g. block LU), sometimes some
> of the blocks can be obtained as a matrix solve.
>
> However, for such complex use cases, I think that it wasn't too bad
> (although not wonderful either) to force using named xpr objects
> instead of temporaries.
>
> As Keir notes, here we're giving a lot of weight to a quite unusual
> use case. I don't think it's worth going for option 4 just for that
> use case, as that would make the API heavier for everybody just to
> make it lighter in that unusual use case.
>
>>> On Thu, Sep 10, 2009 at 11:37 PM, Gael Guennebaud
>>> <gael.guennebaud@xxxxxxxxx> wrote:
>>> > you forgot the main advantage of sol 3 and 4 that is to be able to use
>>> > temporary proxy object without having to declare them, e.g.:
>>> >
>>> > opt 3: x.start(n) = A.solve(b); // works without temporary
>
> You're right, I was forgetting about that unusual use case, but again,
> i've never met it outside of the implementation of blocked algorithms.
>
> Now it's true that it has intrinsic value, and that it's a good goal
> for our API design to fully support the seamless implementation of
> blocked algorithms.
>
> My main concern, like Keir, was that this option 3 wouldn't allow the
> result of solve() to be used in expressions, because you earlier said
> that ReturnByValue shouldn't inherit MatrixBase. Now I didn't think as
> much as you did about this, but now you say,
>
>>> >
>>> > So why not giving a try to opt. 3 ? Indeed, actually the multiple
>>> > matrix products are now implemented using something very very similar
>>> > to opt 3 with a ReturnByValue inheriting MatrixBase, and they are not
>>> > that heavy.
>
> Ah OK, so now we're again discussing a ReturnByValue that inherits
> MatrixBase. So option 3 indeed looks pretty good, can be used in any
> expression.
>
>>> > What I mean is that all the limitations of opt.3 with a
>>> > ReturnByValue inheriting MatrixBase already exists with the matrix
>>> > products... (the limitation is tha the NestByValue flag *must* be
>>> > honored, otherwise you have a weird (currently) compilation error).
>
> OK, that's not too bad, especially as, as you said, this issue already
> exists with all the products (and there it really is the best solution
> anyway).
>
> At least we'll have some uniformity. What's needed then is a good dox
> page on writing functions that take MatrixBase as argument.
>
>>> > Now, honestly I don't really like to use pointers for return types,
>>> > unless they are optional.
>
> Do you mean: as opposed to natural syntax, result = function() ? As in option 3.
> Then I can understand. Me too, i prefer natural syntax!
>
> Or do you mean: as opposed to function(result.ref()) or passing by c++
> reference? Then I need to ask, why? What is bad about pointer, as
> opposed to references or custom structs? OK custom structs allow
> passing temporaries, but are you also saying that you prefer c++
> references over pointers? Then I really would like to know why.
>
> Benoit
>
>
>