title: sum of overlapping arrays, auto-vectorization, and restrict
body: <p>Arstechnia recently had an article <a href="http://arstechnica.com/information-technology/2014/05/ask-ars-why-are-some-programming-languages-faster-than-others/" rel="nofollow noreferrer">Why are some programming languages faster than others</a>.  It compares Fortran and C and mentions summing arrays.  In Fortran it's assumed that arrays don't overlap so that allows further optimization.  In C/C++ pointers to the same type may overlap so this optimization can't be used in general.  However, in C/C++ one can use the <code>restrict</code> or <code>__restrict</code> keyword to tell the compiler not to assume the pointers overlap.  So I started looking into this in regards to auto-vectorization.</p>
  
  <p>The following code vectorizes in GCC and MSVC</p>
  
  <pre><code>void dot_int(int *a, int *b, int *c, int n) {
      for(int i=0; i&lt;n; i++) {
          c[i] = a[i] + b[i];
      }
  }
  </code></pre>
  
  <p>I tested this with and without overlapping arrays and it gets the correct result. However, the way I would vectorize this loop manually with SSE does not handle overlapping arrays.</p>
  
  <pre><code>int i=0;    
  for(; i&lt;n-3; i+=4) {
      __m128i a4 = _mm_loadu_si128((__m128i*)&amp;a[i]);
      __m128i b4 = _mm_loadu_si128((__m128i*)&amp;b[i]);
      __m128i c4 = _mm_add_epi32(a4,b4);
      _mm_storeu_si128((__m128i*)c, c4);
  }
  for(; i&lt;n; i++) {
      c[i] = a[i] + b[i];
  }
  </code></pre>
  
  <p>Next I tried using <code>__restrict</code>.  I assumed that since the compiler could assume the arrays don't overlap it would not handle overlapping arrays but both GCC and MSVC still get the correct result for overlapping arrays even with <code>__restrict</code>.</p>
  
  <pre><code>void dot_int_restrict(int * __restrict a, int * __restrict b, int * __restrict c, int n) {
      for(int i=0; i&lt;n; i++) {
          c[i] = a[i] + b[i];
      }
  }
  </code></pre>
  
  <p><strong>Why does the auto-vectorized code with and without <code>__restrict</code> get the correct result for overlapping arrays?</strong>.</p>
  
  <p>Here is the full code I used to test this:</p>
  
  <pre><code>#include &lt;stdio.h&gt;
  #include &lt;immintrin.h&gt;
  void dot_int(int *a, int *b, int *c, int n) {
      for(int i=0; i&lt;n; i++) {
          c[i] = a[i] + b[i];
      }
      for(int i=0; i&lt;8; i++) printf("%d ", c[i]); printf("\n"); 
  }
  
  void dot_int_restrict(int * __restrict a, int * __restrict b, int * __restrict c, int n) {
      for(int i=0; i&lt;n; i++) {
          c[i] = a[i] + b[i];
      }
      for(int i=0; i&lt;8; i++) printf("%d ", c[i]); printf("\n"); 
  }
  
  void dot_int_SSE(int *a, int *b, int *c, int n) {
      int i=0;    
      for(; i&lt;n-3; i+=4) {
          __m128i a4 = _mm_loadu_si128((__m128i*)&amp;a[i]);
          __m128i b4 = _mm_loadu_si128((__m128i*)&amp;b[i]);
          __m128i c4 = _mm_add_epi32(a4,b4);
          _mm_storeu_si128((__m128i*)c, c4);
      }
      for(; i&lt;n; i++) {
          c[i] = a[i] + b[i];
      }
      for(int i=0; i&lt;8; i++) printf("%d ", c[i]); printf("\n"); 
  }
  
  int main() {
      const int n = 100;
      int a[] = {1,1,1,1,1,1,1,1};
      int b1[] = {1,1,1,1,1,1,1,1,1};
      int b2[] = {1,1,1,1,1,1,1,1,1};
      int b3[] = {1,1,1,1,1,1,1,1,1};
  
      int c[8];
      int *c1 = &amp;b1[1];
      int *c2 = &amp;b2[1];
      int *c3 = &amp;b3[1];
  
      dot_int(a,b1,c, 8);
      dot_int_SSE(a,b1,c,8);
  
      dot_int(a,b1,c1, 8);
      dot_int_restrict(a,b2,c2,8);
      dot_int_SSE(a,b3,c3,8);
  
  }
  </code></pre>
  
  <p>The output (from MSVC)</p>
  
  <pre><code>2 2 2 2 2 2 2 2 //no overlap default
  2 2 2 2 2 2 2 2 //no overlap with manual SSE vector code
  2 3 4 5 6 7 8 9 //overlap default
  2 3 4 5 6 7 8 9 //overlap with restrict
  3 2 2 2 1 1 1 1 //manual SSE vector code
  </code></pre>
  
  <p>Edit:</p>
  
  <p>Here is another inserting version which produces much simpler code</p>
  
  <pre><code>void dot_int(int * __restrict a, int * __restrict b, int * __restrict c, int n) {
      a = (int*)__builtin_assume_aligned (a, 16);
      b = (int*)__builtin_assume_aligned (b, 16);
      c = (int*)__builtin_assume_aligned (c, 16);
      for(int i=0; i&lt;n; i++) {
          c[i] = a[i] + b[i];
      }
  }
  </code></pre>