That is certainly doable, although the arbitrary number of dimensions is a bit of a pain.
Is there a reason you punish yourself with linear arrays? Or is this just code that is being converted from Fortran to C? Or are you or your colleagues more proficient in Fortran than C? (I'm not teasing or critizising you, I'm just trying to understand why you're going in this direction.)
I would use a data structure instead:
Code:
#define _ISO_C99_SOURCE
typedef double data_t;
typedef struct {
long size;
long refcount;
data_t data[];
} owner_t;
typedef struct {
long dimensions;
long *size;
long *step;
owner_t *owner;
data_t *data;
} array_t;
/* Each array_t variable needs to be initialized via
* array_t a = ARRAY_INIT;
* or
* array_t a;
* array_init(&a);
*/
#define ARRAY_INIT { 0, NULL, NULL, NULL, NULL }
static void array_init(array_t *const arrayptr)
{
if (arrayptr)
*arrayptr = (array_t)ARRAY_INIT;
}
/* array_ptr() returns a pointer to the element,
* or 'outside' if the element is outside the array.
*/
static inline data_t *array_ptr(const array_t *const arrayptr, const long coord[], data_t *const outside)
{
long index = 0;
long i;
/* No array? */
if (!arrayptr || arrayptr->dimensions < 1L)
return outside;
for (i = 0; i < arrayptr->dimensions; i++)
if (coord[i] >= 0 || coord[i] < arrayptr->size[i])
index += coord[i] * arrayptr->step[i];
else
return outside;
return arrayptr->origin + index;
}
/* array_get() and array_set() return array_outside,
* if the referred element is outside bounds.
*/
static data_t array_outside = 0.0;
static inline data_t array_get(const array_t *const arrayptr, const long coord[])
{
return *array_ptr(arrayptr, coord, &array_outside);
}
static inline data_t array_set(const array_t *const arrayptr, const long coord[], const data_t value)
{
data_t *const ptr = array_ptr(arrayptr, coord, NULL);
if (ptr)
return *ptr = data;
else
return array_outside;
}
If you look carefully, you can see that each dimension has their own
step. (The step is in data units, not bytes; so 1 means one
data_t element.)
This means you can have your data in any order (compare row-major and column-major for 2D arrays), without modifying the data. You can even swap dimensions.
The fact that the data is stored separately from the array means you can have more than one array referring to the same data. (As long as you destroy every array you no longer need, the data is discarded when the last array referring to it is discarded.)
In your case, if your consumers were to consume data in
array_t format, then you wouldn't need to copy any data; you would just create another
array_t referring to the same data. (These are usually called
views, by the way.)
As long as the view is regular, you can use any data order. You can pick any plane, or even a diagonal. This is like C or Python slicing on steroids.
What this data structure cannot do, is fold multidimensional data into lower dimensions. It can trivially
select a lesser-dimensional slice, it just cannot fold it; it cannot create a linear sequence from more than one dimension. That's why all consumers would have to consume the
array_t format data on input.