title: Calculate average using SSE with STL vectors
body: <p>I'm trying to learn about vectorisation, and rather than reinvet the wheel I'm using <a href="http://www.agner.org/optimize/#vectorclass" rel="nofollow">Agner Fog's vector library</a></p>
  
  <p>Here's my original C++/STL code</p>
  
  <pre><code>#include &lt;vector&gt;
  #include &lt;vectorclass.h&gt;   
  template&lt;typename T&gt;
  double mean_v1(T begin,T end) {
          float mean = 0;
          std::for_each(begin,end,[&amp;mean](const double&amp; d) { mean+=d; });
  
      return mean / std::distance(begin,end);
  }
  
  double mean_v2(T begin,T end) {
      float mean = 0;
      const int distance = std::distance(begin,end); // This is expensive
      const int loop = ( distance &gt;&gt; 2)+1; // divide by 4
      const int partial = distance &amp; 2; // remainder 4
      Vec4d vec;
      for(int i = 0; i &lt; loop;++i) {
          if(i == (loop-1)) {
              vec.load_partial(partial,&amp;*begin);
              mean = horizontal_add(vec);
          }
          else  {
              vec.load(&amp;*begin);
              mean = horizontal_add(vec);
              begin+=4; // This is expensive
          }
      }
      return mean / distance;
  }
  
  int main(int argc,char**argv) {
      using namespace boost::assign;
      std::vector&lt;float&gt; numbers;
      // Note 13 numbers, which won't fit into a sse register perfectly
      numbers+=39.57,39.57,39.604,39.58,39.61,31.669,31.669,31.669,31.65,32.09,33.54,32.46,33.45;
  
      const float mean1 = mean_v1(numbers.begin(),numbers.end());
      const float mean2 = mean_v2(numbers.begin(),numbers.end());
  
  
      return 0;
  }
  </code></pre>
  
  <p>Both v1 and v2 work correctly and they both take about the same time. However profiling it shows the the std::distance() and moving the iterator along takes almost 45% of the total time. The vector adds is just 0.8% which is significantly faster than v1.</p>
  
  <p>Searching the web, all the examples seem to deal with perfect number of values that fit precisely into the SSE registers. How do people deal with odd numbers of values eg for this example where setting up the loop is taking a lot longer than the calculation.</p>
  
  <p>I'm thinking there must be best practices or ideas on how to deal with this scenario.</p>
  
  <p>Assume I can't change the interface of mean() to take float[], but must use iterators</p>