Hi,

I am really interested in using Eigen for software development for the geostatistical software in Stanford.

I would like a class that has all those implementations but for a 3D matrix instead of just 2D. For example I would like to have a 3D model (i.e. a cube) and be able to store the values in a matrix and be able to perform all those operations available in Eigen.

For example: I would like something like this:

std::vector<std::vector<std::vector>>> matrix3D;

so I can access its 3-dimensional elements like this matrix3D(i,j,k)

How can I use Eigen (or other libraries) to do this ? What should I change ? How can I get the best performance too ?

Best,
Mehrdad

I've been using boost::multi_array so far for anything with > 2 dimensions. You really don't want to be using nested std::vector's or std::arrays! Pain lies that way. I also need to map 1D of multi_array's into Eigen::VectorXds but I haven't worked out how to do this yet.

I'd be interested to see if Eigen's developers might add in support for this, or if (as I suspect) their intention is to keep Eigen purely for linear algebra related concepts (which is fair enough!).

Frequently I need to be able to store tensor type quantities, such as second order derivatives of functions in 3 dimensions with symmetry so support would be nice.

Eigen in principle doesn't go beyond 2D (usual matrices), however in the development branch we have a Array type which is like vectors/matrix except that all operations work coefficient-by-coefficient. Array can be used as the Scalar type for Matrix. This is perfectly efficient for fixed sizes, not so much for dynamic sizes (where it will cause 1 malloc per coefficient). So you can do, for example, Matrix<Array4f,4,4> can be seen as a "cube" 4x4x4 3D matrix.

Tensor images are possible, I am using them myself. Afaik, the only thing you need is to define Eigen::NumTraits (see also here: http://eigen.tuxfamily.org/api/Customiz ... ScalarType)

Code: Select all

  template <>
  struct NumTraits< Vector2f >
  {
    typedef NumTraits<Vector2f::RealScalar>::Real Real;
    typedef Vector2f FloatingPoint;
    enum {
      IsComplex = 0,
      HasFloatingPoint = NumTraits<Real>::HasFloatingPoint,
      ReadCost = 2,
      AddCost = 2 * NumTraits<Real>::AddCost,
      MulCost = 2 * NumTraits<Real>::MulCost
    };
    inline static float dummy_precision() { return NumTraits<Real>::dummy_precision(); }
  };

  typedef Matrix<Vector2f, Dynamic, Dynamic> TensorImage;


The more important question is what you gain by doing this? You still don't have arithmetic - you can't view slices as matrices, etc. So if you do it, you are just storing data in a 2D grid...

- Hauke

Thanks for your helps.
In my case, I only need to read-write to the coefficients of this 3D cube. More precisely, I need to have access to different slices (3D sub-matrices) of the cube, and I can store them as row-vectors and use Eigen afterwards to do my calculations.

Since it is Dynamic size, I prefer not to use Matrix of Arrays. After long search all around, I think my best options are either: (1) boost-multiArray (2) Armadillo (3) macstl (although for mac, but it can be used for every environment). I don't know which was has the best performance for multi-dimensional arrays ?

PS: I will use them to read/write to my cube and reshape the 3D submatrices into a VectorXf, and then, use Eigen for all subsequent calculations.

Mehrdad

mehrdadh wrote:PS: I will use them to read/write to my cube and reshape the 3D submatrices into a VectorXf, and then, use Eigen for all subsequent calculations.

If you anyways plan to do that you may as well use what I proposed before; you just need to specialize the NumTraits for VectorXf.

You access calls might look a little bit strange (e.g. m(row,col)(tensor_index)) but you can overcome this by adding a single method through the plugin system.

- Hauke