Re: [boost] ublas interoperability

From: Joerg Walter (jhr.walter_at_[hidden])
Date: 2002-09-30 15:17:07

----- Original Message -----
From: "Kresimir Fresl" <fresl_at_[hidden]>
To: <boost_at_[hidden]>
Sent: Thursday, September 26, 2002 6:47 PM
Subject: Re: [boost] ublas interoperability

>
> Michael Stevens wrote:
>
> >> From:
> >> nbecker_at_[hidden] (Neal D. Becker)
>
> >> Can ublas matrixes interconvert to/from plain old C pointers (to
> >> interop with legacy code)?  If so, how can a Matrix be converted to a
> >> pointer to the start of storage?  Is conversion in the other direction
> >> possible?
>
> > It is very easy to get a pointer to the start of storage.
>
> Unfortunately, not always. Or, to be more exact, it seems that there
> is no simple `universal' solution.
>
> > The uBLAS
> > contianer (matrix,vector etc) all wrap an underling storage type. This
> > is usualy the 'unbounded_array' type. You can easily get a reference to
> > this array using the accessor functions '.data()'.
>
> > From then it is easy to get a pointer to the begin of the storage
> > using '.begin()'
>
> > Therefore to get a pointer to the first element of matrix A used
> >     A.data().begin()
>
> Yes, ublas::unbounded_array<>::begin() returns
> ublas::unbounded_array<>::iterator (or const_iterator),
> which *happens* to be a pointer. (The same is true
> for ublas::bounded_array<>.)
>
> But, for other storage types this need not be true; e.g. storage
> type can be std::vector<> and its iterator can be `proper' class
> (and it is in gnu library).
>
> One can try
>
>     &A.data()[0]
>
> that is, take the address of the first element of storage array.
>
> But this does not work when storage type is
> ublas::(un)bounded_array<> and A is const. Namely,
> return type of `operator[]() const' is
> ublas::type_traits<>::const_reference, which for
> built-in types (e.g. float and double) is not reference,
> but value (e.g. not `double const&' but `double').

OK, here we're discussing at the border of genericity vs. optimization. If
genericity would be of greater importance, we could disable this
optimization.

> Maybe:
>
>     &*A.data().begin()

i.e. &A.data().begin()[0] ?

>
> This works both with ublas::(un)bounded_array<>
> and std::vector<>. But if storage type defines
> const_iterator type whose const_reference is
> ublas::type_traits<>::const_reference, there's again
> a problem.

OK, but the set of supported storage containers currently is limited to
ublas::(un)bounded_array and std::vector. Maybe we should consider
boost::array, too (if your patch is accepted eventually ;-) For that
scenario someone would have to implement a new container.

> Maybe (let's suppose that A is const boost::vector<double>):
>
>    &const_cast<vector<double>&> (A).data()[0]
>
> i.e. remove vector's const; or:
>
>    &const_cast<unbounded_array<double>&> (A.data())[0]
>
> or (more general):
>
>    &const_cast<vector<double>::array_type&> (A.data())[0]
>
> i.e. remove const from storage.

You've lost me :-)

> I think that the proper solution is some kind of a traits class;
> something like:
>
>    template <typename V>
>    struct storage_traits {
>      typedef typename V::pointer pointer;
>      static pointer storage (V& v) {
>        return &v.data()[0];
>      }
>    };
>
>    template <typename V>
>    struct storage_traits<V const> {
>      typedef typename V::const_pointer pointer;
>      static pointer storage (V const& v) {
>        V& vv = const_cast<V&> (v);
>        return &vv.data()[0];
>      }
>    };
>
>    template <typename V>
>    inline typename storage_traits<V>::pointer storage (V& v) {
>      return storage_traits<V>::storage (v);
>    }
>
> Of course, all these work if storage type has continuous
> storage.

This sounds like a reasonable alternative.

> [...]
>
> > I don't think it is possible at present to go the other way.
>
> Oh, it is. There is (undocumented) class array_adaptor<>
> (in storage.hpp):
>
>    double b[] = { 1., 2. };
>    array_adaptor<double> aa1 (2, b);    // `2' is array size
>    vector<double, array_adaptor<double> > v1 (2, aa1);
>
>    double *c = new double[2];
>    array_adaptor<double> aa2 (2, c);
>    vector<double, array_adaptor<double> > v2 (2, aa2);
>    c[0] = 1.; c[1] = 2.;
>    cout << v2 (0) << "  " << v2 (1) << endl;
>
> [...]

Good to know, that someone else knows this stuff :-) BTW, array_adaptor is
currently undocumented and in a certain sense experimental, because it
encapsulates reference counting, which is the base of other implementations
(still not sure, if this is correct ;-).

Regards

Joerg