|Re: [eigen] still the solve() API debate|
[ Thread 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: Thu, 10 Sep 2009 23:37:57 +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=0RKs8bB1+s//EfYIui7ySVc2xIIlW7IBUvzChOT3VtQ=; b=sOwrIr+WQ1It5ZpBdhVhDO9wRKgm4UpWbW95uQXwu+PJxrJ61NOLlpiDXZoSzyGISp RKSOSeZzOQIZ6HVKyLn9Aun27qeOQarXoLdkTHjaaFxKcqvxY+NnvOXWDYCd21ivlq0d qUQbzdq25x9mPV4BpcifChm+Dk01Lff+rmt8I=
- 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=eJ/EinrhkckV4pw/BjTjv1IaZjiPw8PmEGQDu8H7tfzNBWWdxDF46V9oHekoTMKP5f ey6yLMGsR3fnDCG5xLn7DKcCFc9gMSYyiIrT6nZNrKgKxprf9F0hLqaxqvbGCm3C4Xhb ELktde+ISDigwaxI+9FgwN8bEoSaI50haTWow=
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
opt 4: A.solve(b, x.start(n).ref());
VectorBlock<TypeOfVectorX, Dynamic> start_of_x = x.start(n);
The advantage of opt 4 versus opt 3, is that it is much simpler to implement.
Now, honestly I don't really like to use pointers for return types,
unless they are optional. Also I don't really like to have two
different API to do the same job.
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. 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).
Moreover, writing a solver using opt 3 should be less verbose than for
the products. Basically:
template<typename DerivedB> LU::solve(const MatrixBase<DerivedB>& b)
would return a LUSolver<LUType,DerivedB> class inheriting
NestByValue<LuSolver<LUType,DerivedB>, Matrix<Scalar,XXX,XXX> > where
the second template argument is the default evaluation matrix type.
Then you only have to implement 3 trivial functions:
- the ctor copying a ref to a LU and DerivedB objects,
- the functions rows() and cols()
and the main template<DerivedX> evalTo(DerivedX& x) function which
does all the job or call a function of LU, e.g.:
return m_lu.solve(m_b, &x);
On Thu, Sep 10, 2009 at 10:27 PM, Keir Mierle <mierle@xxxxxxxxx> wrote:
> Nice writeup.
> On Thu, Sep 10, 2009 at 12:04 PM, Benoit Jacob <jacob.benoit.1@xxxxxxxxx>
>> we really need to make a decision here to be able to move on...
>> *** The options ***
>> Here are the options on the table:
>> A.solve(b, &result)
>> Pro: simple, functional, predictable
>> Con: not very nice API; some object that passing ptrs isn't C++-ish.
>> result = A.solve(b) // returning a plain Matrix
>> Pro: most natural, simple, nice API
>> Con: the result may be copied redundantly which hurts performance.
>> result = A.solve(b) // returning a proxy object
>> Pro: nice, natural API in this case; no performance problem
>> Con: Unpredictable since A.solve(b) can't be used in an expression;
>> adds more code complexity
>> A.solve(b, result.ref())
>> Pro: combines the advantages of 1. while pleasing people who don't like
>> Con: ref() is more wordy than &; and It's still not as convenient as 2.
>> *** My proposal ***
>> I propose that we keep 1., where we actually put the implementation
>> (this is how we currently do), and alongside it, also offer 2. in the
>> API, as a trivial wrapper around 1.
> Simple. Avoids complicated code. Supports both use cases. RVO may be
> solveable; can you make a small case where RVO doesn't happen? Rather
> specific conditions must apply for RVO to happen.
>> Like this:
>> template<typename RightHandSideType, typename ResultType>
>> void Solver::solve(const RightHandSideType& b, ResultType *result)
>> /// bla bla, put the solver implementation there
>> template<typename RightHandSideType>
>> PlainMatrixType solve(const RightHandSideType& b)
>> PlainMatrixType m;
>> solve(b, &m);
>> return m;
>> Here are my arguments. Sorry if this contradicts opinions that i
>> expressed before.
>> I see 2 nontrivial points that I need to justify:
>> a) why plain pointer and not .ref()
>> b) why return by value and not ReturnByValue
>> Let's start with b). The advantage of returning a plain matrix is
>> obvious; the drawback is the potential performance problem. What
>> happens in practice? There are two cases, either the result is a
>> fixed-size matrix or it is dynamic-size. If the result is a fixed-size
>> matrix then (i did some experiments...) it seems that the RVO works
>> perfectly, so there's no performance issue at all. If the matrix is
>> dynamic-size (array on the heap) then indeed the RVO doesn't seem to
>> happen (by the way that's mysterious to me, why is that?) so we indeed
>> pay for a redundant malloc/copy/free. But does it matter? We're
>> talking about dynamic size, so we optimize for large sizes, where this
>> overhead is relatively small (the solving itself has cubic complexity
>> for matrix solve and quadratic complexity for a vector solve, so it's
>> always bigger!) Thus, I don't really believe anymore in the
>> performance argument against 2.
>> Plus, we still offer 1. in the API for the user who's concerned about
>> explicitly avoiding this potential performance issue.
>> Now let's discuss a). Really my only argument for &result instead of
>> result.ref() is that it's shorter to type and the result is more
>> obvious to everybody. What would be arguments for ref() in this case?
>> If I remember well the arguments for ref() were to honor constness in
>> swap() and have a single, uniform syntax for output args everywhere.
>> But there's a problem: it's heavy to write, especially in swap().
>> We're not going to force the user to write
>> So at least swap would have to escape that uniformization... so after
>> all, sorry, i'm afraid that the ref() idea is blocked by how heavy it
>> makes the API :( internally, I agree that it was a very good idea.
>> Anyway, in this particular case: A.solve(b,&result), there's no const
>> correctness issue, and the only intrinsic advantage of ref() over a
>> pointer is that many people don't like pointers in C++. But, here i
>> need to be honest: i don't really understand why. They're almost the
>> same thing as references (at least, _const_ pointers are. do you like
>> it better if we make it a const pointer)? I understand that they don't
>> have the same nice properties in some cases, but does any of that
>> matter in the present use case?
>> If you would prefer us to standardize on references instead of
>> pointers, speak up, i don't want to impose a decision onto the rest of
>> people here. I just thought that the advantage of using pointers here,
>> namely that in A.solve(b, &result) it is immediately clear which one
>> is the result, would outweigh the drawbacks (tell me if i'm missing