Re: [ublas] [bindings] [lapack] geev and the type of the input matrix

Subject: Re: [ublas] [bindings] [lapack] geev and the type of the input matrix
From: Marco Guazzone (marco.guazzone_at_[hidden])
Date: 2010-07-12 07:22:50

On Tue, Jul 6, 2010 at 11:44 AM, Marco Guazzone
<marco.guazzone_at_[hidden]> wrote:
> On Tue, Jul 6, 2010 at 11:34 AM, Thomas Klimpel
> <Thomas.Klimpel_at_[hidden]> wrote:
>> Marco Guazzone wrote:
>>> Just another question (sorry for abusing of your patience)
>>>
>>> I'd like to create an "eigen" function (as a wrapper of lapack::geev)
>>> which accepts as input matrix both:
>>> * a const& matrix (since it is not modified anymore), and
>>
>> What do you mean by "(since it is not modified anymore)"? lapack will definitively overwrite the original matrix during the computation and not restore it at the end.
>
> Sorry,
> This refer to the point of view of users of my lib...
> but as you have pointed out I've just realized that LAPACK overwrites it
>
> A       (input/output) DOUBLE PRECISION array, dimension (LDA,N)
>          On entry, the N-by-N matrix A.
>          On exit, A has been overwritten.
>
> So to keep it consistent with the point of view of my lib, I need to
> create a temp matrix and pass it to lapack::geev
>
>
>>
>>> * a matrix_expression (in order to allow the user to call say
>>> eigen(trans(A),...) )
>>
>> If you want that, you probably have to create a temporary matrix that lapack can modify during the computation.
>
> OK
>
> Thank you so much!!
>

Hi,

I've found in geev.hpp that some definition of geev takes a "const&"
for the input matrix.
For instance (the input matrix is "a"):

--- [code] ---
template< typename MatrixA, typename VectorW, typename MatrixVL,
        typename MatrixVR, typename Workspace >
inline typename boost::enable_if< detail::is_workspace< Workspace >,
        std::ptrdiff_t >::type
geev( const char jobvl, const char jobvr, const MatrixA& a, VectorW& w,
        const MatrixVL& vl, const MatrixVR& vr, Workspace work ) {
    return geev_impl< typename bindings::value_type<
            MatrixA >::type >::invoke( jobvl, jobvr, a, w, vl, vr, work );
}
--- [/code] ---

However, inside of these functions "a" is passed as-is, whithout
making a temporary copy.
Moreover, none of the "invoke" methods seems to accept a "const&" for
the input matrix.

Is it correct?

PS: I've found similar pattern in other bindings::lapack files.

Thank you very much!

Best,

-- Marco