This package provides classes for various arrays. All of the classes and functions are defined in the array namespace. To use this package you can either include the header for the class you are using
#include "array/MultiArray.h"
or you can include the convenience header.
#include "array/all.h"
(Here I assume that stlib/packages is in your include path, i.e. g++ -Istlib/packages ....)
There are five array classes:
There are five multidimensional array classes:
Each of these classes is templated on the value type and the dimension (rank). For example MultiArray has the following declaration.
template<typename _T, std::size_t _Dimension> class MultiArray;
The value type can be any class or built-in type. Below we construct a 3-dimensional array of integers.
array::MultiArray<int, 3> a;
Each of the multidimensional array classes have different constructors. Consult the class documentation for details.
We can construct a multidimensional array by specifying the index extents. We use an array of std::size_t for the extents. Below we define some types and construct a 3x4 2-D array.
typedef array::MultiArray<double, 2> MultiArray; typedef MultiArray::SizeList SizeList; MultiArray a(SizeList(3, 4));
The analagous C multidimensional array would be declared as:
double ca[3][4];
We can treat 2-D array as simply a container of 12 elements. The multidimensional array classes support the standard STL-style interface. For those of you who are familiar with the C++ Standard Template Library (STL), each of the multidimensional arrays classes fulfill the requirements of a random access container. For the rest of you, I suggest you read "Generic Programming and the STL" by Matthew H. Austern.
assert(! a.empty());
assert(a.size() == 12);
assert(a.max_size() == 12);
std::fill(a.begin(), a.end(), 0);
for (std::size_t i = 0; i != a.size(); ++i) {
a[i] = i;
} Note that operator[] performs container indexing and not multidimensional array indexing.
The array classes use operator() to perform multidimensional indexing. The argument is an array of integers. The index range for the array we have declared is [0..2]x[0..3]. Below we define the multidimensional index type and iterate over the elements and assign each the sum of the indices.
typedef MultiArray::IndexList IndexList;
IndexList i;
for (i[0] = 0; i[0] != 3; ++i[0]) {
for (i[1] = 0; i[1] != 4; ++i[1]) {
a(i) = sum(i);
}
} The array we declared has indices that are zero-offset. That is, the lower bounds for the index ranges are zero. We can also declare an array with different index bases. The following array has index ranges [1..3]x[1..4].
MultiArray b(SizeList(3, 4), IndexList(1, 1));
The extents() and bases() accessors give the array extents and the index bases.
const SizeList extents = b.extents();
const IndexList bases = b.bases();
const IndexList upper = bases + extents;
IndexList i;
for (i[0] = bases[0]; i[0] != upper[0]; ++i[0]) {
for (i[1] = bases[0]; i[1] != upper[1]; ++i[1]) {
b(i) = sum(i);
}
} Each of the multidimensional array classes is an STL-compliant random access containers. The following types are related to this functionality. First the types defined in both the mutable and constant array classes.
value_type is the element type of the array.const_pointer is a pointer to a constant array element.const_iterator is an iterator on constant elements in the array.const_reverse_iterator is a reverse iterator on constant elements in the array.const_reference is a reference to a constant array element.size_type is the size type.difference_type is the pointer difference type.The mutable array classes define mutable versions of the pointers, iterators, and references.
pointer is a pointer to an array element.iterator is an iterator on elements in the array.reverse_iterator is a reverse iterator on elements in the array.reference is a reference to an array element.The remainder of the types start with a capital letter. They support array indexing and working with different views of an array.
Parameter is the parameter type. This is used for passing the value type as an argument. If the value type is a built-in type then it is value_type, otherwise it is const value_type&.Index is a single array index, which is a signed integer.IndexList is an Array of indices: Array<int, N>. This type is used for multidimensional array indexing.SizeList is an Array of sizes: Array<std:size_t, N>. This type is used to describe the index extents.Storage is a class that specifies the storage order.Range is a class that represents the index range for the array.ConstView is the class for a constant view of the array: MultiArrayConstView<value_type, N>.View is the class for a mutable view of the array: MultiArrayView<value_type, N>.The multidimensional array classes provide a number of member functions for using the array as an STL random access container.
empty() returns true if the array has zero elements.size() returns the number of elements.max_size() returns the number of elements as well. (Dynamically-sized containers like std::vector return the maximum number of elements that the container could hold, which is determined by the integer precision. For statically-sized containers, the max_size() is the same as the size().operator[]() returns the specified element in the array. Below we sum the elements of an array and check it against the sum() function.
value_type s = 0;
for (size_type i = 0; i != a.size(); ++i) {
s += a[i];
}
assert(s == sum(a)); operator[] performs container indexing and not array indexing. In the following example we create a 1-D array with index range [-5..5] but use container indexing to initialize the elements.
typedef array::MultiArray<double, 1> MultiArray;
typedef MultiArray::SizeList SizeList;
typedef MultiArray::IndexList IndexList;
typedef MultiArray::size_type size_type;
MultiArray a(SizeList(11), IndexList(-5));
for (size_type i = 0; i != a.size(); ++i) {
a[i] = i;
} begin() returns a random access iterator to the first element.end() returns a random access iterator to one past the last element. Below we copy the elements of an array a to a buffer. std::copy(a.begin(), a.end(), buffer);
rbegin() returns a random access reverse iterator to the last element.rend()returns a random access reverse iterator to one past the first element. Below we copy the elements from a to b in reverse order. assert(a.size() == b.size()); std::copy(a.rbegin(), a.rend(), b.begin());
fill() fills the array with the specified value. The following two lines are equivalent. std::fill(a.begin(), a.end(), 1); a.fill(1);
The array classes offer only one way to perform multidimensional array indexing: use operator() with an IndexList as an argument. This helps keep the interface simple and makes it unlikely that one will confuse array indexing and container indexing. Below we create a 3-D array of integers with index range [-5..5]x[-5..5]x[-5..5] and set the element values to the product of the indices.
typedef array::MultiArray<int, 3> MultiArray;
typedef MultiArray::SizeList SizeList;
typedef MultiArray::IndexList IndexList;
MultiArray a(SizeList(11, 11, 11), IndexList(-5, -5, -5));
const IndexList lower = a.bases();
const IndexList upper = a.bases() + a.extents();
IndexList i;
for (i[0] = lower[0]; i[0] != upper[0]; ++i[0]) {
for (i[1] = lower[1]; i[1] != upper[1]; ++i[1]) {
for (i[2] = lower[2]; i[2] != upper[2]; ++i[2]) {
a(i) = product(i);
}
}
} The MultiIndexRange class, as the name suggests, is used to represent index ranges. We can describe a continuous index range by its extents and bases. Consider an index range in 2-D. Below we construct the range [0..2]x[0..3]
typedef array::MultiIndexRange<2> Range; typedef Range::SizeList SizeList; Range range(SizeList(3, 4));
If the index range is not zero-offset, we can specify the index bases. Below we construct the range [1..3]x[1..4].
typedef Range::IndexList IndexList; Range range(SizeList(3, 4), IndexList(1, 1));
We can specify a non-continuous index range by specifying steps to take in each dimension. Below we construct the index range [1, 3, 5]x[1, 4, 7, 10].
Range range(SizeList(3, 4), IndexList(1, 1), IndexList(2, 3));
The range() member function returns the index range for each of the the multidimensional array classes. Below we store the range for an array. We check that the extents, bases, and steps are correct using the MultiIndexRange accessor member functions.
typedef MultiArray<double, 2> MultiArray; typedef MultiArray::SizeList SizeList; typedef MultiArray::IndexList IndexList; typedef MultiArray::Range Range; MultiArray a(SizeList(10, 20)); Range range = a.range(); assert(range.extents() == a.extents()); assert(range.bases() == a.bases()); assert(range.steps() == IndexList(1, 1));
MultiIndexRange is not useful by itself. Rather, we use it to construct index range iterators. MultiIndexRangeIterator is a random access iterator over an index range. Below is a function that sets each element of an array to the sum of its indices.
template<std::size_t _N>
void
setToSum(array::MultiArray<int, _N>* a) {
typedef array::MultiIndexRangeIterator<_N> Iterator;
const Iterator end = Iterator::end(a->range());
for (Iterator i = Iterator::begin(a->range()); i != end; ++i) {
(*a)(*i) = sum(*i);
}
} We construct iterators to the beginning and end of the index range with the static member functions begin() and end(). Dereferencing the iterator yields an index.
Next we consider a more sophisticated example. Given a N-D arrays a and b, we set the boundary elements of b to the same value as in a, and set the interior elements to the average of the adjacent neighbors in a. In 1-D, bi = (ai-1 + ai+1) / 2 for the interior elements. In 2-D we have bi,j = (ai-1,j + ai+1,j + ai,j-1 + ai,j+1) / 4.
template<typename _T, std::size_t _N>
void
laplacianAverage(const array::MultiArray<T, _N>& a, array::MultiArray<_T, _N>* b) {
assert(a.range() == b->range());
typedef array::MultiIndexRangeIterator<_N> Iterator;
typedef typename Iterator::IndexList IndexList;
// Get the boundary values by copying the entire array.
*b = a;
// Skip if there are no interior points.
if (min(a.extents()) <= 2) {
return;
}
// The range for the interior elements.
const array::MultiIndexRange<_N> range(a.extents() - 2, a.bases() + 1);
// Compute the interior values.
_T s;
IndexList index;
const Iterator end = Iterator::end(range);
for (Iterator i = Iterator::begin(range); i != end; ++i) {
s = 0;
for (std::size_t n = 0; n != _N; ++n) {
index = *i;
index[n] -= 1;
s += a(index);
index[n] += 2;
s += a(index);
}
(*b)(*i) = s / _T(2 * _N);
}
} The MultiArrayRef and MultiArrayConstRef classes are multidimensional arrays that reference externally allocated data. Their constructors differ from MultiArray in that the first argument is a pointer to the data. Below is a function that receives C arrays for the data, extents and bases and constructs a 3-D MultiArrayRef.
void
foo(double* data, const std::size_t[3] extents, const int[3] bases) {
typedef MultiArrayRef<double, 3> MultiArrayRef;
typedef MultiArrayRef::SizeList SizeList;
typedef MultiArrayRef::IndexList IndexList;
MultiArrayRef a(data, SizeList(extents), IndexList(bases));
...
} For MultiArrayConstRef the constructors take const pointers.
MultiArrayView and MultiArrayConstView are array classes that can be used to create views of existing arrays. The easiest way to construct them is to use the view() member function with a specified index range. Below we create a view of the interior elements of an array.
typedef array::MultiArray<double, 2> MultiArray; typedef MultiArray::SizeList SizeList; typedef MultiArray::IndexList IndexList; typedef MultiArray::Range Range; typedef MultiArray::View View; MultiArray a(SizeList(12, 12), IndexList(-1, -1)); View interior(Range(SizeList(10, 10)));
When using the view() function, the index bases of the array view are the same as the index bases of the specified range. Below is a function that returns a view of the interior elements of an array.
template<typename _T, std::size_t _N>
array::MultiArrayView<_T, _N>
interior(array::MultiArray<_T, _N>& a) {
typedef typename array::MultiArray<_T, _N>::Range Range;
assert(min(a.extents()) > 2);
return a.view(Range(a.extents() - 2, a.bases() + 1));
} The array view classes are fully fledged multidimensional arrays, just like MultiArray, MultiArrayRef, and MultiArrayConstRef. In particular, they have iterators. Below we use the interior() function be defined above and set the interior elements of an array to zero.
MultiArray a(SizeList(12, 12), IndexList(-1, -1));
MultiArrayView<double, 2> = interior(a);
typedef MultiArrayView<double, 2>::iterator iterator;
const iterator end = interior.end();
for (iterator i = interior.begin(); i != end; ++i) {
*i = 0;
} These iterators are standard random access iterators. In addition, you can get the index list that corresponds to the iterator position with the indexList member function. Below we verify that this works correctly for the interior array.
for (iterator i = interior.begin(); i != end; ++i) {
assert(*i == interior(i->indexList));
} Note that the iterators for MultiArrayView and MultiArrayConstView are less efficient than those for the other multidimensional array types. This is because the data for the view classes is not contiguous in memory. The view array classes use MultiViewIterator, while the rest of the classes use raw pointers for iterators.
MultiArrayBase is a base class for the multidimensional arrays. It stores the index extents and bases as well as the number of elements.
MultiArrayConstView derives from MultiArrayBase. It adds accessors (like begin() and end()) that make it a constant, random access container. It also adds constant array indexing with operator().
MultiArrayView derives from MultiArrayConstView using virtual inheritance. It adds the mutable interface.
MultiArrayConstRef derives from MultiArrayConstView using virtual inheritance. It provides more efficient constant iterators and adds constant container indexing.
MultiArrayRef derives from both MultiArrayConstRef and MultiArrayView. It provides more efficient mutable iterators and adds mutable container indexing.
MultiArray derives from MultiArrayRef. Other than the constructors (and related functions), it hase the same interface as its base. It adds memory allocation.
When writing functions which take these arrays as arguments, it is best to use the least derived class which supplies the needed interface. Consider the following function.
template<typename _T, std::size_t _N> void foo(array::MultiArray<_T, _N>& a);
Since the function only uses the constant interface, we can make the function more general by accepting either a MultiArrayConstRef or a MultiArrayConstView. If the function uses iterators on the array, then use the former because arrays with contiguous data have more efficient iterators.
template<typename _T, std::size_t _N> void foo(array::MultiArrayConstRef<_T, _N>& a);
Otherwise use the latter.
Consider a function that modifies an array.
template<typename _T, std::size_t _N> void bar(array::MultiArray<_T, _N>* a);
We can make the function more general by accepting either a MultiArrayRef or a MultiArrayView. Again the choice is usually based on whether we use the array iterators. One needs to use a MultiArray as the argument type only when changing the size of the array.
Why, oh why, did I write a multidimensional array package? Hasn't someone already done that? Yes, there are a couple of good libraries. I wrote this package because I wasn't entirely satisfied with other libraries. Also, I was worried about depending on external libraries. In writing this package I have stolen many ideas from Blitz++ and the Boost multidimensional array library.
Blitz++ is a full-featured library and has good documentation. It is a powerful library that uses template meta-programming to generate efficient code. One thing that I don't like about Blitz++ is that it is not entirely STL-compliant and breaks some usual C++ conventions. For instance, the copy constructor does not copy an array, but references the data. (This is a feature to keep naive programmers from shooting themselves in the foot.) Also, I've had some compilation problems with the template expressions.
Multidimensional arrays are provided in the multi_array package of the Boost C++ libraries. It is a nice design and has some cool tricks. However, in using this library I found it hard to write dimension-independent code, i.e. write functions and classes where the dimension is a template parameter.