1

I'm new to C and having a problem with array types when embedded in structures. The following is an example of my issue:

typedef struct {
    double J[151][151];
} *UserData;

static int PSolve(void *user_data, N_Vector solution)
{
UserData data;
data = (UserData) user_data;

double J[151][151];
J = data->J;

/* Solve a matrix equation that uses J, stored in 'solution' */

return(0);
}

When I try to compile this I get error: incompatible types when assigning to type ‘double[151][151]’ from type ‘double (*)[151]’

My current workaround for this has been to replace 'J[x][y]' by 'data->J[x][y]' in the code to solve the matrix equation, but profiling has shown this to be less efficient.

Changing the arguments to PSolve is not an option because I'm using the sundials-cvode solver which prescribes the type and order of the arguments.

3
  • possible duplicate of C array declaration and assignment?
    – outis
    Feb 21, 2011 at 21:25
  • Couldn't you just read off of data->J instead of copying it to a local array? Or will you be modifying it within that function? Feb 21, 2011 at 21:32
  • I mentioned that as a workaround I was reading off data->J[x][y] in the matrix equation. But this function gets called a few billion times during an execution. The data->J is appearing to be a performance bottleneck. I'm trying to get a local copy as a workaround.
    – Sevenless
    Feb 21, 2011 at 21:39

4 Answers 4

5
typedef struct {
    double J[151][151];
} UserData; // This is a new data structure and should not a pointer!

static int PSolve(void *user_data, N_Vector solution)
{
UserData* data; // This should be a pointer instead!
data = (UserData*) user_data;

double J[151][151];
memcpy(J, data->J, sizeof(double) * 151 * 151); // use memcpy() to copy the contents from one array to another

/* Solve a matrix equation that uses J, stored in 'solution' */

return(0);
}
5
  • You can simply use sizeof J in the memcpy().
    – caf
    Feb 21, 2011 at 21:47
  • I updated my code to reflect the new structure definition. But memcpy() is not a viable option as it added about 20% to the total computational time relative to not having a local copy and referencing J[x][y] as data->J[x][y]. The primary issue is data->J is appearing to be a bottleneck, so i'm going to try @caf 's suggestion of using they keyword 'restrict'.
    – Sevenless
    Feb 21, 2011 at 22:37
  • +1 for the suggestion to not make the structure a pointer. I'm not totally sure why, but it shaved some time off of operations elsewhere in the program.
    – Sevenless
    Feb 21, 2011 at 22:38
  • @Sevenless Dude, don't hijack your own thread. You completely changed the question and invalidated our answers!!! Rollback and leave the old question intact. It has been already answered. If you have another question related to this code (and clearly you have), then ask another question, please! Feb 22, 2011 at 1:29
  • Apologies, I was trying to add new information to get at the root of my question. It's now rolled back. I'll repost the revised question.
    – Sevenless
    Feb 22, 2011 at 4:34
2

The root cause of your issue is that array types cannot be directly assigned in C. You must explicitly use memcpy(), as karlphillip's answer shows.

Note however that performing the copy may wipe out the optimisation gains you make in the remainder of the function. Presumably, the solution argument is a pointer, and the optimiser is worried about potential aliasing between user_data / data and solution. As an alternative, if you have a compiler that supports C99's restrict keyword, is to use that qualifier on the arguments:

static int PSolve(void * restrict user_data, N_Vector restrict solution)

This promises the compiler that those pointers do not alias, and should allow you to directly use data->J without sacrificing the compiler optimisations.

Some compilers make the restrict keyword available in C89 mode under a spelling such as __restrict - consult your compiler documentation for more details.

3
  • I've added the restrict keyword, but I get: cc1: error: unrecognized command line option "-std=C99" when I compile with ` gcc44 -c -std=C99 -I/home/matteson/sundials/include/ -I/misc/linux/64/opt/pkg/matlab/R2010b/extern/include -I/misc/linux/64/opt/pkg/matlab/R2010b/simulink/include -DMATLAB_MEX_FILE -D_GNU_SOURCE -fexceptions -fPIC -fno-omit-frame-pointer -pthread -fopenmp -DMX_COMPAT_32 -O3 -DNDEBUG "het_kry.c"` I checked the manpage for gcc and it says that's the correct flag.
    – Sevenless
    Feb 21, 2011 at 22:42
  • 1
    @Sevenless: It's a lower-case "c" - -std=c99 (but you might actually want -std=gnu99, which is a C99 version of the normal default gnu89).
    – caf
    Feb 21, 2011 at 23:10
  • @Sevenless: Also, it may be beneficial to add some appropriate const qualifiers. Can you update your question to show the definition of N_Vector?
    – caf
    Feb 21, 2011 at 23:14
0

You can't assign to an array. You should look up about the difference between arrays and pointers in C.

What you need in your code is just something like this:

double (*J)[151] = data->J;

that is a pointer to arrays of length 151. Or if you'd like with typedefs

typedef double line[151];
line *J = data->J;

that's all, you should not copy the data but just the pointer to the data.

Edit: But seeing the whole thread of answers, I think all of that is pure speculation where your bottleneck might be. It could just be anywhere, e.g that you are accessing the matrix the "wrong way round", column wise, or whatever. Or that pumping the data from memory is just dominating your computation.

Perhaps consider to look into the assembler that your compiler produces (option -S) to see if there is something fishy.

5
  • This is unlikely to behave any differently to direct use of data->J later in the function.
    – caf
    Feb 21, 2011 at 21:56
  • @caf, it is one indirection less. So it depends on how data is really declared. Feb 21, 2011 at 23:50
  • @caf @Jens-Gustedt: I edited the question to include the typedef for UserData and to clarify the performance element of the question.
    – Sevenless
    Feb 21, 2011 at 23:59
  • Gustedt: It is exactly the same number of indirections, because data->J is at a constant offset (of zero) from data itself.
    – caf
    Feb 22, 2011 at 0:38
  • @caf, in that case you are probably right. But generally, this works only if the compiler is able to do constant propagation for this offset. Otherwise, this constant offset might hinder direct addressing in the assembler. Feb 22, 2011 at 8:27
-1

Try defining local variable double J[151][151] as double **J.

1
  • no dice, it will compile with a warning about incompatible pointer types, but it segfaults on execution. Thanks for the try. (I wasn't the -1)
    – Sevenless
    Feb 21, 2011 at 21:35

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.