Cast Array to Matrix in C/C++?

Question

If I have a function which internally uses a 3x3 matrix, but due to API limitations must pass it a contiguous array of nine elements, is there a way to essentially cast an array to a 2d array?

The following approach doesn't work, but hopefully it conveys what I'm trying to do.

I'm aware that packages like Eigen make this trivial, but unfortunately it's not an option for me in this case.

void transpose_3x3( double x[3][3] ) {
    ...
}

void foo( double x[9], double y[3][3] ) {

    transpose_3x3(y)  // OK    

    double *my_x[3][3] = &x[0];

    transpose_3x3(&my_x);  // How?

}

Edit

Thanks to Daniel for showing me the correct target for a cast, I was able to use the C++ reinterpret_cast function to accomplish this. While the union and memcpy solutions work, a single explicit cast seems to be the most straightforward solution, IMHO.

#include <iostream>

void print_matrix( double x[3][3] ) {

    for(unsigned int i=0; i< 3; i++) {
        for(unsigned int j=0; j<3; j++) {
            std::cout << x[i][j] << "  ";
        }
        std::cout << std::endl;
    }
}

int main() {

    double a[9] = { 1, 2, 3, 4, 5, 6, 7, 8, 9};    
    double (*b)[3] = reinterpret_cast<double (*)[3]>(a);  

    print_matrix(b);
    return 0;
}

Answer 1

I can't think of a way to do an explicit cast.

However you could just memcpy it. Calls to memcpy are not as stupid as you may think. The compiler will often see that the 2 bits of data actually represent the same thing and that the memcpy is a fixed size and optimise out the actual copy.

I must say I've never tried it with something quite so large as the matrix you are doing but I can't see why it wouldn't work.

Edit: In fact I thought I'd give it a try. I wrote the following code:

void transpose_3x3( double (*x)[3][3] ) 
{
    const double t01    = (*x)[0][1];
    const double t02    = (*x)[0][2];
    const double t12    = (*x)[1][2];
    (*x)[0][1] = (*x)[1][0];
    (*x)[0][2] = (*x)[2][0];
    (*x)[1][0] = t01;
    (*x)[1][2] = (*x)[2][1];
    (*x)[2][0] = t02;
    (*x)[2][1] = t12;
}

void foo() 
{
    double x[9] = { 1.0f, 2.0f, 3.0f,
                    4.0f, 5.0f, 6.0f,
                    7.0f, 8.0f, 9.0f };

    double y[3][3];
    memcpy( y, x, sizeof( double ) * 9 );

    transpose_3x3( &y );

    printf( "%f, %f, %f\n", y[0][0], y[0][1], y[0][2] );
    printf( "%f, %f, %f\n", y[1][0], y[1][1], y[1][2] );
    printf( "%f, %f, %f\n", y[2][0], y[2][1], y[2][2] );
}

And built it in release mode with VS2010.

The resulting assembly is as follows:

void foo() 
{
00E11000  push        ebp  
00E11001  mov         ebp,esp  
00E11003  and         esp,0FFFFFFC0h  
00E11006  sub         esp,0B8h  
    double x[9] = { 1.0f, 2.0f, 3.0f,
00E1100C  fld1  
00E1100E  push        esi  
00E1100F  fstp        qword ptr [esp+2Ch]  
00E11013  push        edi  
00E11014  fld         qword ptr [__real@4000000000000000 (0E12138h)]  
                    4.0f, 5.0f, 6.0f,
                    7.0f, 8.0f, 9.0f };

    double y[3][3];
    memcpy( y, x, sizeof( double ) * 9 );

    transpose_3x3( &y );

    printf( "%f, %f, %f\n", y[0][0], y[0][1], y[0][2] );
00E1101A  sub         esp,18h  
00E1101D  fstp        qword ptr [esp+50h]  
00E11021  mov         ecx,12h  
00E11026  fld         qword ptr [__real@4008000000000000 (0E12130h)]  
00E1102C  lea         esi,[esp+48h]  
00E11030  fstp        qword ptr [esp+58h]  
00E11034  lea         edi,[esp+90h]  
00E1103B  fld         qword ptr [__real@4010000000000000 (0E12128h)]  
00E11041  fst         qword ptr [esp+60h]  
00E11045  fld         qword ptr [__real@4014000000000000 (0E12120h)]  
00E1104B  fstp        qword ptr [esp+68h]  
00E1104F  fld         qword ptr [__real@4018000000000000 (0E12118h)]  
00E11055  fstp        qword ptr [esp+70h]  
00E11059  fld         qword ptr [__real@401c000000000000 (0E12110h)]  
00E1105F  fst         qword ptr [esp+78h]  
00E11063  fld         qword ptr [__real@4020000000000000 (0E12108h)]  
00E11069  fstp        qword ptr [esp+80h]  
00E11070  fld         qword ptr [__real@4022000000000000 (0E12100h)]  
00E11076  fstp        qword ptr [esp+88h]  
00E1107D  rep movs    dword ptr es:[edi],dword ptr [esi]  
00E1107F  fstp        qword ptr [esp+10h]  
00E11083  fstp        qword ptr [esp+8]  
00E11087  fld         qword ptr [esp+90h]  
00E1108E  fstp        qword ptr [esp]  
00E11091  mov         esi,dword ptr [__imp__printf (0E120A0h)]  
00E11097  push        offset string "%f, %f, %f\n" (0E120F4h)  
00E1109C  call        esi  
    printf( "%f, %f, %f\n", y[1][0], y[1][1], y[1][2] );
00E1109E  add         esp,4  
00E110A1  fld         qword ptr [esp+0C8h]  
00E110A8  fstp        qword ptr [esp+10h]  
00E110AC  fld         qword ptr [esp+0B0h]  
00E110B3  fstp        qword ptr [esp+8]  
00E110B7  fld         qword ptr [__real@4000000000000000 (0E12138h)]  
00E110BD  fstp        qword ptr [esp]  
00E110C0  push        offset string "%f, %f, %f\n" (0E120F4h)  
00E110C5  call        esi  
    printf( "%f, %f, %f\n", y[2][0], y[2][1], y[2][2] );
00E110C7  fld         qword ptr [esp+0D4h]  
00E110CE  add         esp,4  
00E110D1  fstp        qword ptr [esp+10h]  
00E110D5  fld         qword ptr [__real@4018000000000000 (0E12118h)]  
00E110DB  fstp        qword ptr [esp+8]  
00E110DF  fld         qword ptr [__real@4008000000000000 (0E12130h)]  
00E110E5  fstp        qword ptr [esp]  
00E110E8  push        offset string "%f, %f, %f\n" (0E120F4h)  
00E110ED  call        esi  
00E110EF  add         esp,1Ch  
}

You'll note that there is no memcpy. Actually all its doing is manually copying the matrix from x into y and then printing it in a transposed fashion. Basically its interesting seeing what things the compiler will do to optimise things out ...

Edit 2: Of course thinking a bit further after seeing paddy's excellent response it occurs to me that you could just case the thing directly

transpose_3x3( (double (*)[3][3])&x );

Which works without memcpy or the union :D

Answer 2

Well, this is a little dirty, but you could use a union with casting. Someone's going to tell me off for this but I'll post it anyway =)

#include <iostream>

using namespace std;

typedef union {
    double linear[9];
    double square[3][3];
} Matrix;

void transpose_3x3( double x[3][3] )
{
    for( int i = 0; i < 3; i++ ) {
        for( int j = 0; j < 3; j++ ) {
            cout << i << "," << j << ": " << x[i][j] << endl;
        }
    }
}

void foo( double x[9], double y[3][3] )
{
    transpose_3x3( y );
    transpose_3x3( ((Matrix*)x)->square );
}

int main()
{
    double x[9] = { 1, 2, 3, 4, 5, 6, 7, 8, 9 };
    double y[3][3] = { {1, 2, 3}, {4, 5, 6}, {7, 8, 9} };
    foo( x, y );
    return 0;
}