Boost MultiDimensional Array Library : 101 Example

[Cuchulainn]

// MyFirstBoost.cpp

//

// Initial code to test the feasibility of Boost Multidimensional Array Library

// in finance applications.

//

// Interested in the types 1) Matrices 2) Tensors

//

//

// (C) Datasim Education BV 2007

//



#include "boost/multi_array.hpp"

#include <iostream>

using namespace std;





int main()

{

using namespace boost;



// Define structure of tensor

const int dim = 3; // 3d matrix

multi_array<double, dim> tensor;



// Define the extents of each separate dimension

const int NT = 3;

const int NSIM = 4;

const int NDIM = 3;

tensor.resize(extents[NT][NSIM][NDIM]);



// Define a 3d loop to initialise the data

typedef boost::multi_array<double, dim>::index index;



for (index i = 0; i != NT; ++i)

{

for (index j = 0; j != NSIM; ++j)

{

for (index k = 0; k != NDIM; ++k)

{

tensor[i][j][k] = i + j + k;

}

}

}



// Useful shorthand

typedef multi_array<double, dim> Tensor;

Tensor tensor2(extents [NT][NSIM][NDIM]);



// Create a tensor using a Collections, in this case in Boost

array<Tensor::index, 3> extentsList = { { NT, NSIM, NDIM } };

Tensor tensor3 (extentsList);



tensor3 = tensor2 = tensor;



for (index i = 0; i != NT; ++i)

{

//cout << "\n\nFirst " << i << endl << endl;

for (index j = 0; j != NSIM; ++j)

{

//cout << "\nSecond " << j << endl;

for (index k = 0; k != NDIM; ++k)

{

tensor[i][j][k] = 0.0;

// cout << i << ", " << j << ", " << k << ": " << tensor[i][j][k] << ", ";

}

}

}



// Accessing elements using a Collection of indices

array<Tensor::index, dim> myIndex = { { 0, 0, 0} };

tensor(myIndex) = 3.1415;



for (index i = 0; i != NT; ++i)

{

cout << "\n\nFirst " << i << endl << endl;

for (index j = 0; j != NSIM; ++j)

{

cout << "\nSecond " << j << endl;

for (index k = 0; k != NDIM; ++k)

{

cout << i << ", " << j << ", " << k << ": " << tensor[i][j][k] << ", ";

}

}

}



//cout << "Value at origin " << tensor(myIndex) << endl;



// Changing the Array base (instead of default base index == 0)

typedef Tensor::extent_range range;

Tensor::extent_gen extents;



// Choice #1: defining bases in the constructor

Tensor tensor4(extents[2] [range(1,4)] [range(-1,3)]);



// Choice #2: set all bases to the same index value

int startIndex = 1;

tensor.reindex (startIndex);



// Choice #3: set each base separately using a Collection

Tensor tensor5(extents[NT][NSIM][NDIM]);

array<Tensor::index, dim> bases = { { 0, 1, -1} };

tensor5.reindex(bases);



// Storage regimes

Tensor tensor6(boost::extents[NT][NSIM][NDIM], boost::fortran_storage_order());

Tensor tensor7(extents[NT][NSIM][NDIM], boost::c_storage_order());



// The ability to customise the order in which dimensions are stored and whether

// the dimensions are stored in ascending or descending order



Tensor tensor8(extents[NT][NSIM][NDIM]); // Dimensions are 0, 1 and 2

Tensor::size_type ordering[] = {2, 0, 1}; // Store last dim, then first, then second



// Store first (dimension 0) in descending order

bool orderRegime[] = {false, true, true};

boost::general_storage_order<dim> storage(ordering, orderRegime);



tensor8 = Tensor(extents [NT][NSIM][NDIM], storage);



// Changing the Array's shape; total number of dimensions must remain the same

Tensor tensor9(extents[NT][NSIM][NDIM]);

array<Tensor::index, dim> newDimensions = { { NSIM, NT, NDIM} };

tensor9.reshape(newDimensions);



// Resizing an array; we can increase the extent of each dimension

Tensor tensor10(extents[NT][NSIM][NDIM]);

tensor10.resize(extents[NT-1] [NSIM] [NDIM]);



// Creating Views, 2 basic types



// View type 1: Selection from original structure based on some

// ranges and possibly strides

const int L1 = 1; const int L2 = 1; const int L3 = 1;

const int U1 = NT; const int U2 = NSIM; const int U3 = NDIM;

const int Stride = 1;



// Now create the view (dim == 3)

typedef Tensor::index_range IndexRange;

Tensor::array_view<dim>::type tensorView =

tensor[boost::indices[IndexRange(L1, U1)] [IndexRange(L2, U2, Stride)] [IndexRange(L3, U3)] ];



for (index i = L1; i != U1; ++i)

{

for (index j = L2; j != U2; j += Stride)

{

for (index k = L3; k != U3; ++k)

{

tensorView[i][j][k] = tensor[i][j][k];

}

}

}



for (index i = L1; i != U1; ++i)

{

cout << "\n First: " << i << endl;

for (index j = L2; j != U2; j += Stride)

{

cout << "\n Second: " << j << endl;

for (index k = L3; k != U3; ++k)

{

cout << tensorView[i][j][k] << ", ";

}

}

}



// A view of one dimension less; this is called SLICING

Tensor::index_gen Indices;

const int dim2 = dim - 1;

Tensor::array_view<dim2>::type reducedView =

tensor[Indices[IndexRange(L1, U1)] [1] [IndexRange(L3, U3)]];



return 0;

}

[Cuchulainn]

Here's the same code in file format.

[thijs]

Code: Select all


// MyFirstBoost.cpp

//

// Initial code to test the feasibility of Boost Multidimensional Array Library

// in finance applications.

//

// Interested in the types 1) Matrices 2) Tensors

//

//

// (C) Datasim Education BV 2007

//



#include "boost/multi_array.hpp"

#include <iostream>

using namespace std;





int main()

{

   using namespace boost;



   // Define structure of tensor

   const int dim = 3; // 3d matrix

   multi_array<double, dim> tensor;



   // Define the extents of each separate dimension

   const int NT = 3;

   const int NSIM = 4;

   const int NDIM = 3;

   tensor.resize(extents[NT][NSIM][NDIM]);



   // Define a 3d loop to initialise the data

   typedef boost::multi_array<double, dim>::index index;



   for (index i = 0; i != NT; ++i)

   {

      for (index j = 0; j != NSIM; ++j)

      {

         for (index k = 0; k != NDIM; ++k)

         {

            tensor[i][j][k] = i + j + k;

         }

      }

   }



   // Useful shorthand

   typedef multi_array<double, dim> Tensor;

   Tensor tensor2(extents [NT][NSIM][NDIM]);



   // Create a tensor using a Collections, in this case in Boost

   array<Tensor::index, 3> extentsList = { { NT, NSIM, NDIM } };

   Tensor tensor3 (extentsList);



   tensor3 = tensor2 = tensor;



   for (index i = 0; i != NT; ++i)

   {

      //cout << "\n\nFirst " << i << endl << endl;

      for (index j = 0; j != NSIM; ++j)

      {

         //cout << "\nSecond " << j << endl;

         for (index k = 0; k != NDIM; ++k)

         {

            tensor[i][j][k] = 0.0;

            // cout << i << ", " << j << ", " << k << ": " << tensor[i][j][k] << ", ";

         }

      }

   }



   // Accessing elements using a Collection of indices

   array<Tensor::index, dim> myIndex = { { 0, 0, 0} };

   tensor(myIndex) = 3.1415;



   for (index i = 0; i != NT; ++i)

   {

      cout << "\n\nFirst " << i << endl << endl;

      for (index j = 0; j != NSIM; ++j)

      {

         cout << "\nSecond " << j << endl;

         for (index k = 0; k != NDIM; ++k)

         {

            cout << i << ", " << j << ", " << k << ": " << tensor[i][j][k] << ", ";

         }

      }

   }



   //cout << "Value at origin " << tensor(myIndex) << endl;



   // Changing the Array base (instead of default base index == 0)

   typedef Tensor::extent_range range;

   Tensor::extent_gen extents;



   // Choice #1: defining bases in the constructor

   Tensor tensor4(extents[2] [range(1,4)] [range(-1,3)]);



   // Choice #2: set all bases to the same index value

   int startIndex = 1;

   tensor.reindex (startIndex);



   // Choice #3: set each base separately using a Collection

   Tensor tensor5(extents[NT][NSIM][NDIM]);

   array<Tensor::index, dim> bases = { { 0, 1, -1} };

   tensor5.reindex(bases);



   // Storage regimes

   Tensor tensor6(boost::extents[NT][NSIM][NDIM], boost::fortran_storage_order());

   Tensor tensor7(extents[NT][NSIM][NDIM], boost::c_storage_order());



   // The ability to customise the order in which dimensions are stored and whether

   // the dimensions are stored in ascending or descending order



   Tensor tensor8(extents[NT][NSIM][NDIM]); // Dimensions are 0, 1 and 2

   Tensor::size_type ordering[] = {2, 0, 1}; // Store last dim, then first, then second



   // Store first (dimension 0) in descending order

   bool orderRegime[] = {false, true, true};

   boost::general_storage_order<dim> storage(ordering, orderRegime);



   tensor8 = Tensor(extents [NT][NSIM][NDIM], storage);



   // Changing the Array's shape; total number of dimensions must remain the same

   Tensor tensor9(extents[NT][NSIM][NDIM]);

   array<Tensor::index, dim> newDimensions = { { NSIM, NT, NDIM} };

   tensor9.reshape(newDimensions);



   // Resizing an array; we can increase the extent of each dimension

   Tensor tensor10(extents[NT][NSIM][NDIM]);

   tensor10.resize(extents[NT-1] [NSIM] [NDIM]);



   // Creating Views, 2 basic types



   // View type 1: Selection from original structure based on some

   // ranges and possibly strides

   const int L1 = 1; const int L2 = 1; const int L3 = 1;

   const int U1 = NT; const int U2 = NSIM; const int U3 = NDIM;

   const int Stride = 1;



   // Now create the view (dim == 3)

   typedef Tensor::index_range IndexRange;

   Tensor::array_view<dim>::type tensorView =

   tensor[boost::indices[IndexRange(L1, U1)] [IndexRange(L2, U2, Stride)] [IndexRange(L3, U3)] ];



   for (index i = L1; i != U1; ++i)

   {

      for (index j = L2; j != U2; j += Stride)

      {

         for (index k = L3; k != U3; ++k)

         {

            tensorView[i][j][k] = tensor[i][j][k];

         }

      }

   }



   for (index i = L1; i != U1; ++i)

   {

      cout << "\n First: " << i << endl;

      for (index j = L2; j != U2; j += Stride)

      {

         cout << "\n Second: " << j << endl;

         for (index k = L3; k != U3; ++k)

         {

            cout << tensorView[i][j][k] << ", ";

         }

      }

   }



   // A view of one dimension less; this is called SLICING

   Tensor::index_gen Indices;

   const int dim2 = dim - 1;

   Tensor::array_view<dim2>::type reducedView =

   tensor[Indices[IndexRange(L1, U1)] [1] [IndexRange(L3, U3)]];



   return 0;

}