|
|
Boost : |
From: Larry Evans (jcampbell3_at_[hidden])
Date: 2001-07-01 11:27:36
Hubert HOLIN wrote:
> Paris (U.E.), le 28/06/2001
>
> --- In boost_at_y..., Larry Evans <jcampbell3_at_p...> wrote:
>
> [SNIP]
>
> > I'm confused by "accessing each index". Do you mean accessing the
> > value corresponding to that index. For example, if the matrix is 'float A[N][M]`
> > and if 'i' is such an iterator, then starting from the initial value of i:
> > *i == A[0][0]
> > *(++i) == A[0][1]
> > ...
> > *(++i) == A[N][M]
> > ?
> > In other words, the iterator produces a 'flattened' version or 'ravel' (apl term)
> > of A, and then access the elements just as a regular stl iterator. Is that
> > what you mean?
>
> Yes.
Attached is code demonstrating an implementation of this. It allows
permuting of the axes, as shown by the last test. Comments welcome.
>
>
> The quirk is that on the one hand the traversal can be somewhat
> more complicated than this (if required by the problem these iterators
> are used for), and on the other hand not every index need be accessed
> but only representatives of certain classes (such as in the case of a
> row index), in which case a further level of iteration is needed after
> dereferencing.
Could you clarify please? I'm guessing that maybe you mean, for rank3 array, A,
A[all][3][all].begin() would be an iterator for a rank 2 array? And by
"further level of iteration" do you mean something like (*i).begin()?
//Purpose:
// Demonstrate section 6.3 "Sequential Access" of
// Budd, Timothy
// _Apl Compiler_
// Springer-Verlag, 1988
//
#include <iostream>
#include <iomanip>
#include <numeric>
#include <vector>
using namespace std;
typedef
int
offset_type
;
template
< offset_type Rank
, typename Elem
>
class
vec_shape
: public vector<Elem>
{
public:
typedef
vector<Elem>
super_type
;
vec_shape(void)
: super_type(Rank)
{}
vec_shape(const Elem x[Rank])
: super_type(Rank)
{ copy(x,x+Rank,begin())
;}
}
;
offset_type const
rank
= 3
;
typedef
vec_shape<rank,offset_type>
index_type
;
typedef
vec_shape<rank+1,offset_type>
evec_type
;
typedef
vector<evec_type>
array_type
;
template
< typename Index
>
ostream&
operator<<(ostream&sout, vector<Index> const& x)
{
; int i=0
; int n=x.size()
; for(;i<n;++i) sout<<x[i]<<endl
; return sout
;}
template
< offset_type Shape
, typename Index
>
ostream&
operator<<(ostream&sout, vec_shape<Shape,Index> const& x)
{
; int i=0
; int n=Shape
; for(;i<n;++i) sout<<setw(2)<<x[i]<<" "
; return sout
;}
class
index_permut
: public index_type
{
public:
index_permut(void)
{ for(unsigned i=0; i<size(); ++i){at(i)=i;}
;}
void
swap(unsigned i, unsigned j)
{ std::swap(at(i),at(j))
;}
}
;//end index_permut class
void
shape2evec
( const index_type&a_shape
, vector<offset_type>::reverse_iterator a_evec
)
{
; partial_sum(a_shape.rbegin(),a_shape.rend(),a_evec,multiplies<offset_type>())
;}
void
shape2evec
( const index_type&a_shape
, evec_type&a_evec
)
{
; a_evec[rank]=1
; shape2evec(a_shape,++a_evec.rbegin())
;}
class
ravel_iter
{
public:
ravel_iter
( const index_permut& a_permut
, const index_type& a_shape
, const array_type& a_array
)
: m_offset_inp(0)
, m_offset_out(0)
, m_array(a_array)
{
; cout<<"ravel_iter:a_permut="<<a_permut<<endl
; evec_type l_evec_inp
; shape2evec(a_shape,l_evec_inp)
; cout<<"ravel_iter:l_evec_inp="<<l_evec_inp<<endl
; index_type l_shape_out
; index_type l_evec_out
; unsigned n = rank()
; int i
; for(i=0; i<n; ++i)
{
; l_shape_out[a_permut[i]] = a_shape[i]
; l_evec_out[a_permut[i]] = l_evec_inp[i+1]
;}
//At this point:
// here [Budd88,sec6.3]
// ---- ---------------
// l_shape_out[i] s[i+1]
// l_evec_out[i] gamma[i+1]
; cout<<"ravel_iter:l_shape_out="<<l_shape_out<<endl
; cout<<"ravel_iter:l_evec_out="<<l_evec_out<<endl
; i = n-1
; m_delta[i] = l_evec_out[i]
; for(--i; -1<i; --i)
{
; m_delta[i] = l_evec_out[i] - l_evec_out[i+1] * l_shape_out[i+1]
;}
; cout<<"ravel_iter:m_delta="<<m_delta<<endl
; shape2evec(l_shape_out,m_evec_out.rbegin())
; cout<<"ravel_iter:m_evec_out="<<m_evec_out<<endl
//At this point:
// here [Budd88,sec6.3]
// ---- ---------------
// m_delta[i] delta[i+1]
// m_evec_out[i] *\s[n-i...1] ( where / here is apl scan op == stl partial_sum )
;}
offset_type
offset(void)const
{ return m_offset_inp
;}
offset_type
inc(void)
//This implements the "Value Phase" algorithm
//in [Budd88,sec6.3], except m_evec_out[i] is used
//instead of s[i] and t to avoid mulitiplications.
{
; ++m_offset_out
; unsigned n = rank()
; int i = n - 1
; for(; -1<i; --i)
{
; m_offset_inp += m_delta[i]
; if(m_offset_out % m_evec_out[i] != 0 ) break
;}
; return m_offset_inp
;}
offset_type
operator++(void)
{ return inc()
;}
evec_type const&
operator()(void)const
{ return m_array[offset()]
;}
unsigned
size(void)const
{ return m_array.size()
;}
unsigned
rank(void)const
{ return m_delta.size()
;}
index_type
indices(void)const
//The following implements [Budd88,p62]
{
; index_type ndx
; unsigned n = rank()
; unsigned i = 0
; unsigned s = 0
; if(1<n)
{ ndx[i] = m_offset_out/m_evec_out[i+1] //a[1] in [Budd88]
; for(++i; i<n-1; ++i)
{ s = m_evec_out[i]/m_evec_out[i+1]
; ndx[i] = ( m_offset_out/m_evec_out[i+1] ) % s //eq(5.7)
;}
;}
; s = m_evec_out[i]
; ndx[i] = m_offset_out % s //a[n] in [Budd88]
; return ndx
;}
private:
index_type
m_evec_out
//This is used in place of s[i] and t in [Budd88,sec6.3]
;
index_type
m_delta
//This is delta in [Budd88,sec6.3]
;
offset_type
m_offset_out
//This is offset (i.e. "output" offset) in [Budd88,sec6.3]
;
offset_type
m_offset_inp
//This is offset' (i.e. "input" offset) in [Budd88,sec6.3]
;
const
array_type&
m_array
;
}
;//end ravel_iter class
ostream&
operator<<
( ostream&sout
, ravel_iter&r_iter
)
{ unsigned n = r_iter.size()
; unsigned i = 0
; for(; i<n; ++i, ++r_iter) sout<<"["<<r_iter.indices()<<"]="<<r_iter()<<endl
; return sout
;}
int
main(void)
{ offset_type ishape[rank]={3,4,2}
; index_type shape(ishape)
; evec_type evec
; shape2evec(shape,evec)
; cout<<"array shape="<<shape<<endl
; cout<<"array evec="<<evec<<endl
; index_type ndx
; int ir=-1
; array_type arr(evec[0])
; array_type::iterator iter(arr.begin())
; for(++ir, ndx[ir]=0; ndx[ir]<shape[ir]; ++ndx[ir])
{ for(++ir, ndx[ir]=0; ndx[ir]<shape[ir]; ++ndx[ir])
{ for(++ir, ndx[ir]=0; ndx[ir]<shape[ir]; ++ndx[ir])
{ evec_type&ev=*iter
; copy(ndx.rbegin(),ndx.rend(),ev.rbegin())
; offset_type off_set =
inner_product
( ndx.rbegin()
, ndx.rend()
, evec.rbegin()
, (offset_type)0
, plus<offset_type>()
, multiplies<offset_type>()
)
; ev[0] = off_set
; ++iter
;}
; --ir
;}
; --ir
;}
; cout<<"arr="<<endl
; cout<<arr<<endl
; {
; cout<<"r_iter w id permut"<<endl
; index_permut ip
; ravel_iter r_iter(ip,shape,arr)
; cout<<"r_iter print:"<<endl
; cout<<r_iter
;}
; {
; cout<<"r_iter w trans permut"<<endl
; index_permut ip
; ip.swap(0,2)
; ravel_iter r_iter(ip,shape,arr)
; cout<<"r_iter print:"<<endl
; cout<<r_iter
;}
; return 0
;}
Boost list run by bdawes at acm.org, gregod at cs.rpi.edu, cpdaniel at pacbell.net, john at johnmaddock.co.uk