Re: [ublas] [bindings][lapack] Initial high-level solve

Subject: Re: [ublas] [bindings][lapack] Initial high-level solve
From: Rutger ter Borg (rutger_at_[hidden])
Date: 2009-04-03 02:55:34

Thomas Klimpel wrote:

> Karl Meerbergen wrote:
>> In modern linear algebra notation (after Wilkinson, eg.), a vector is a
>> column vector, its transpose is a row vector. This is simple and I
>> propose to keep things simple. If a routine requires a matrix with one
>> row, simply use trans(v).
> 
> Good point.
> So how about as_matrix(b) for a Nx1 column major order matrix and
> as_c_matrix(b) for a Nx1 row major order matrix? But I can see for myself
> that it somehow won't work out. Even if the vector already has a vector
> bindings traits class, the stride of the vector might be different from 1,
> so that it's impossible to interpret an arbitrary vector as a Nx1 matrix.

If we would have a working is_bindable_vector metafunction, it becomes
possible to automatically catch vectors and translate them to matrices,
also for decorated types like trans(b). 

Strides are indeed a different issue, I think the leading dimension could be
used to correctly represent strides for row vectors. For column vectors 
I'm not so sure of a trick yet, this would require us to either assert
those strides to 1, or to rewrite the problem itself (e.g., by transposing
A like Jesse suggested).

> I also see a completely different problem, independent of using
> "column_major::as_1xN_matrix(b)" or "trans(b)": The solve template will
> capture it's argument by a non-const reference, but 
> "column_major::as_1xN_matrix(b)" or "trans(b)" will return temporary proxy
> objects. If I remember correctly, it is not allowed to pass temporary
> objects to an argument captured by a non-const reference.

This problem is already in effect for the pivot vector in the solve
function. Some of the low-level drivers have two versions of their static
member functions solve, overloaded on the constness of their arguments. 

Cheers, 

Rutger