[ander]

Hi
I am going to develop a 2D structured grid explicit Euler solver in C++. I am planning to use 1D dynamic arrays to store the conserved variables. My question is if it is preferable to use an array of structures (aos) or structure of arrays (soa) in this case.

The aos approach could look like the following

struct Conserved {
    double density, mom_x, mom_y, energy;
};

Conserved* field = new Conserved[N*N];

The soa approach could look something like this:

struct ConservedField {
    double* density;
    double* mom_x;
    double* mom_y;
    double* energy;
};

ConservedField field;
//allocate memory for field

Which approach is likely to yield the best performance?

[sbaffini]

Depends on how you actually access your variables the most and how the compiler (or language) promises to lay down variables in the aos case.

If you mostly access them only one at the time, the soa is definitely what you need.

But what if you use a solver that mostly needs to access all of them at once (i.e., a coupled solver)? In this case you would probably prefer the aos version, PROVIDED that consecutive array elements are stored consecutively. But I am no expert in C++, so I have no idea if that ever happens or under what conditions it does.

As usual, testing is needed for the specific case.

[aerosayan]

Aight. You're where I started one year ago.

Finally a topic in which I consider myself somewhat experienced

Store all of values together using AOS.

Most of your solver's time will be spent in calculating the fluxes and residuals. So you need the data to be together in cache.

Also, just focus on parallelizing your code.

Don't even bother with improving the performance of each thread, because it doesn't matter that much.

Your code will be limited by the memory bandwidth, and if you use SOA, you're going to be doing lots of reads and writes to very far away positions in memory.

If you choose SOA, your fluxes might also be in SOA form, which will be bad because you will be doing lots of write operations during the flux calculation stage.

Go through my post history in this forum.

You'll see that I had the same questions and after some research I found out that memory bandwidth is the limiting factor.

My code was fast, too fast, so much that the memory bandwidth couldn't keep up and I was getting lots of cache misses.

That's what you want to avoid.

Look up roofline model of performance analysis, you'll see that your performance is limited by how many reads you do.

sbaffini had previously said to calculate all things on the fly, which only made sense after I did my own analysis too. So, additionally, calculate the volume, cell area, etc. on the fly. That is, don't even bother storing them. Those reads will slow you down. Just calculate them directly from the coordinates and then use them.

[ander]

Luckily I choose AOS randomly before I read this, so I won't have to bother with rewriting everything
I'll check out your post history and the other stuff, when i have time!
I also have a question about operator overloading. Is this acceptable to use? I'm thinking about using it in the cases where the same operations are done on all 4 entries in the vectors, like for instance when timestepping. It would make the code easier to read and write, but I have heard they are slow.

[aerosayan]

Modern compilers can optimize operator overloading. Although check the disassembly to make sure you're getting optimized code. If you don't know what optimized code looks like in assembly, you don't need to worry about it right now.

If you're operator overloading, don't go too deep. Nested operator overloads might not get optimized correctly if you have 3 or 4 levels of indirection.