title: Using nested vectors vs a flatten vector wrapper, strange behaviour
body: <p><strong>The problem</strong></p>
  
  <p>For a long time I had the impression that using a nested <code>std::vector&lt;std::vector...&gt;</code> for simulating an N-dimensional array is in general bad, since the memory is not guarantee to be contiguous, and one may have cache misses. I thought it's better to use a flat vector and map from multiple dimensions to 1D and vice versa. So, I decided to test it (code listed at the end). It is pretty straightforward, I timed reading/writing to a nested 3D vector vs my own 3D wrapper of an 1D vector. I compiled the code with both <code>g++</code> and <code>clang++</code>, with <code>-O3</code> optimization turned on. For each run I changed the dimensions, so I can get a pretty good idea about the behaviour. To my surprise, these are the results I obtained on my machine MacBook Pro (Retina, 13-inch, Late 2012), 2.5GHz i5, 8GB RAM, OS X 10.10.5:</p>
  
  <p><strong>g++ 5.2</strong></p>
  
  <pre><code>dimensions       nested   flat
  X   Y   Z        (ms)     (ms) 
  
  100 100 100  -&gt;  16       24
  150 150 150  -&gt;  58       98
  200 200 200  -&gt;  136     308
  250 250 250  -&gt;  264     746
  300 300 300  -&gt;  440    1537
  </code></pre>
  
  <p><strong>clang++ (LLVM 7.0.0)</strong></p>
  
  <pre><code>dimensions       nested   flat
  X   Y   Z        (ms)     (ms) 
  
  100 100 100  -&gt;  16       18
  150 150 150  -&gt;  53       61
  200 200 200  -&gt;  135     137
  250 250 250  -&gt;  255     271
  300 300 300  -&gt;  423     477
  </code></pre>
  
  <p><br></p>
  
  <p>As you can see, the "flatten" wrapper is never beating the nested version. Moreover, g++'s libstdc++ implementation performs quite badly compared to libc++ implementation, for example for <code>300 x 300 x 300</code> the flatten version is almost 4 times slower than the nested version. libc++ seems to have equal performance. </p>
  
  <p><strong>My questions:</strong></p>
  
  <ol>
  <li>Why isn't the flatten version faster? Shouldn't it be? Am I missing something in the testing code?</li>
  <li>Moreover, why does g++'s libstdc++ performs so badly when using flatten vectors? Again, shouldn't it perform better?</li>
  </ol>
  
  <p><strong>The code I used:</strong></p>
  
  <pre><code>#include &lt;chrono&gt;
  #include &lt;cstddef&gt;
  #include &lt;iostream&gt;
  #include &lt;memory&gt;
  #include &lt;random&gt;
  #include &lt;vector&gt;
  
  // Thin wrapper around flatten vector
  template&lt;typename T&gt;
  class Array3D
  {
      std::size_t _X, _Y, _Z;
      std::vector&lt;T&gt; _vec;
  public:
      Array3D(std::size_t X, std::size_t Y, std::size_t Z):
          _X(X), _Y(Y), _Z(Z), _vec(_X * _Y * _Z) {}
      T&amp; operator()(std::size_t x, std::size_t y, std::size_t z)
      {
          return _vec[z * (_X * _Y) + y * _X + x];
      }
      const T&amp; operator()(std::size_t x, std::size_t y, std::size_t z) const
      {
          return _vec[z * (_X * _Y) + y * _X + x];
      }
  };
  
  int main(int argc, char** argv)
  {
      std::random_device rd{};
      std::mt19937 rng{rd()};
      std::uniform_real_distribution&lt;double&gt; urd(-1, 1);
  
      const std::size_t X = std::stol(argv[1]);
      const std::size_t Y = std::stol(argv[2]);
      const std::size_t Z = std::stol(argv[3]);
  
  
      // Standard library nested vector
      std::vector&lt;std::vector&lt;std::vector&lt;double&gt;&gt;&gt;
          vec3D(X, std::vector&lt;std::vector&lt;double&gt;&gt;(Y, std::vector&lt;double&gt;(Z)));
  
      // 3D wrapper around a 1D flat vector
      Array3D&lt;double&gt; vec1D(X, Y, Z);
  
      // TIMING nested vectors
      std::cout &lt;&lt; "Timing nested vectors...\n";
      auto start = std::chrono::steady_clock::now();
      volatile double tmp1 = 0;
      for (std::size_t x = 0 ; x &lt; X; ++x)
      {
          for (std::size_t y = 0 ; y &lt; Y; ++y)
          {
              for (std::size_t z = 0 ; z &lt; Z; ++z)
              {
                  vec3D[x][y][z] = urd(rng);
                  tmp1 += vec3D[x][y][z];
              }
          }
      }
      std::cout &lt;&lt; "\tSum: " &lt;&lt; tmp1 &lt;&lt; std::endl; // we make sure the loops are not optimized out
      auto end = std::chrono::steady_clock::now();
      std::cout &lt;&lt; "Took: ";
      auto ms = std::chrono::duration_cast&lt;std::chrono::milliseconds&gt;(end - start).count();
      std::cout &lt;&lt; ms &lt;&lt; " milliseconds\n";
  
      // TIMING flatten vector
      std::cout &lt;&lt; "Timing flatten vector...\n";
      start = std::chrono::steady_clock::now();
      volatile double tmp2 = 0;
      for (std::size_t x = 0 ; x &lt; X; ++x)
      {
          for (std::size_t y = 0 ; y &lt; Y; ++y)
          {
              for (std::size_t z = 0 ; z &lt; Z; ++z)
              {
                  vec1D(x, y, z) = urd(rng);
                  tmp2 += vec1D(x, y, z);
              }
          }
      }
      std::cout &lt;&lt; "\tSum: " &lt;&lt; tmp2 &lt;&lt; std::endl; // we make sure the loops are not optimized out
      end = std::chrono::steady_clock::now();
      std::cout &lt;&lt; "Took: ";
      ms = std::chrono::duration_cast&lt;std::chrono::milliseconds&gt;(end - start).count();
      std::cout &lt;&lt; ms &lt;&lt; " milliseconds\n";
  }
  </code></pre>
  
  <p><strong>EDIT</strong></p>
  
  <p>Changing the <code>Array3D&lt;T&gt;::operator()</code> return to</p>
  
  <pre><code>return _vec[(x * _Y + y) * _Z + z];
  </code></pre>
  
  <p>as per <a href="https://stackoverflow.com/a/33093906/3093378">@1201ProgramAlarm's suggestion</a> does indeed get rid of the "weird" behaviour of g++, in the sense that the flat and nested versions take now roughly the same time. However it's still intriguing. I thought the nested one will be much worse due to cache issues. <strong>May I just be lucky and have all the memory contiguously allocated?</strong> </p>