メモリアライメントを意識した動的確保を行う際は posix_memalign が良いぽい

posix_memalign() が良い理由

本ブログでは今まで _mm_malloc() を使用してメモリアライメントを揃えていましたが，posix_memalign() を使用したほうが良いぽいです．理由は _mm_malloc() が Intel CPU 環境でのみ動作するのに対し，posix_memalign() は標準化されているためメインストリームな環境であれば概ね動作するためです．また，posix_memalign() で確保したメモリは free() で解放することができ，_mm_free() が混在しなくなるためコードの保守性の向上が考えられます．

posix_memalign の使用方法

簡単な使用方法のみ掲載します．
・インクルードするファイル

#include <stdlib.h>

・posix_memalign() の使用方法

int posix_memalign(void **memptr, size_t alignment, size_t size)
//memptr: メモリの先頭
//alignment: 揃えたいメモリアライメント
//size: 確保したいメモリのサイズ

・メモリの解放方法
free() でメモリ解放可能

posix_memalign() に直した総和計算 with SIMD

/*--------------------
  posix_memalign_sum.cpp
  --------------------*/

#include <iostream>
#include <cstdlib>
#include <ctime>
#include <numeric>
#include <immintrin.h>

float vec_sum(std::size_t const n, float const* x){
  std::size_t static const para = 8;
  std::size_t const end = n / para;
  std::size_t const beg = end * para;

  auto vsum = _mm256_setzero_ps();
  auto sum = reinterpret_cast<float*>(&vsum);
  auto vx = reinterpret_cast<__m256 const*>( x );

  for(std::size_t i=0; i<end; ++i)
    vsum = _mm256_add_ps(vsum, vx[i]);
  for(std::size_t i=beg; i<n; ++i)
    sum[0] += x[i];
  sum[0] += sum[1] + sum[2] + sum[3] + sum[4] + sum[5] + sum[6] + sum[7];
  
  return sum[0];
}

double vec_sum(std::size_t const n, double const* x){
  std::size_t static const para = 4;
  std::size_t const end = n / para;
  std::size_t const beg = end * para;

  auto vsum = _mm256_setzero_pd();
  auto sum = reinterpret_cast<double*>(&vsum);
  auto vx = reinterpret_cast<__m256d const*>( x );

  for(std::size_t i=0; i<end; ++i)
    vsum = _mm256_add_pd(vsum, vx[i]);
  for(std::size_t i=beg; i<n; ++i)
    sum[0] += x[i];
  sum[0] += sum[1] + sum[2] + sum[3];
  
  return sum[0];
}

int main(void){
  std::size_t const n = 1024;
  float* x;
  double* y;

  //x = static_cast<float*>( _mm_malloc(sizeof(float) * n, 32) );
  //y = static_cast<double*>( _mm_malloc(sizeof(double) * n, 32) );
  posix_memalign(reinterpret_cast<void**>(&x), 32, sizeof(float) * n);
  posix_memalign(reinterpret_cast<void**>(&y), 32, sizeof(double) * n);

  srand( static_cast<unsigned>(time(NULL)) );
  for(std::size_t i=0; i<n; ++i) x[i] = static_cast<float>( rand() ) / RAND_MAX;
  for(std::size_t i=0; i<n; ++i) y[i] = static_cast<double>( rand() ) / RAND_MAX;

  std::cout << "Single Precision\n";
  std::cout << "-Accumulate: " << std::accumulate(x, x+n, 0.0) << std::endl;  
  std::cout << "-SIMD: " << vec_sum(n, x) << std::endl;

  std::cout << "\nDouble Precision\n";
  std::cout << "-Accumulate: " << std::accumulate(y, y+n, 0.0) << std::endl;  
  std::cout << "-SIMD: " << vec_sum(n, y) << std::endl;

  //_mm_free(x);
  //_mm_free(y);
  free(x);
  free(y);

  return 0;
}