CUDAの挙動を再現するC++のプログラム

#pragma once

struct COORD {
  int x, y;
  COORD(int _x, int _y) {
    x = _x;
    y = _y;
  }
  COORD(int _x) {
    x = _x;
    y = 1;
  }
  COORD() {
    x = 0;
    y = 0;
  }
};

using MyBlockIdx = COORD;
using MyThreadIdx = COORD;

using DimBlock = COORD;
using DimGrid = COORD;

struct thData {
  MyBlockIdx blockIdx;
  DimBlock blockDim;
  MyThreadIdx threadIdx;
};

//! @brief grid x block回関数fを呼び出す
//! @param grid グリッド
//! @param block ブロック
//! @param f 呼び出す関数
//! @param args N個の引数
template<class F, class ...Args>
void fakeCudaKernel(DimBlock grid, DimGrid block, F&& f, Args&&... args) {

  thData _thread;
  _thread.blockDim = block;

  for (int gx = 0; gx < grid.x; gx++) {
    for (int gy = 0; gy < grid.y; gy++) {
      for (int bx = 0; bx < block.x; bx++) {
        for (int by = 0; by < block.y; by++) {
          _thread.blockIdx.x = gx;
          _thread.blockIdx.y = gy;
          _thread.threadIdx.x = bx;
          _thread.threadIdx.y = by;
          f(_thread, args...);

        }
      }

    }
  }
}

#include <iostream>
#include <algorithm>

#include "fakeCUDA.hpp"

#pragma warning(disable:4996)

//カーネルに該当する関数
void line(
  thData t,
  unsigned int* pimage,
  const int width,
  const int height,
  const int sx,
  const int sy,
  const int ex,
  const int ey) {

  //アクセス法
  //このスレッドが担当する画素のx位置を二次元座標で求める
  int xpos = t.blockIdx.x * t.blockDim.x + t.threadIdx.x;
  //int ypos = blockIdx.y * blockDim.y + threadIdx.y;

  double a = (ey - sy) / double(ex - sx);//傾き
  double b = sy - a * sx;//切片

  //xposは0～になっているので、offsetを足す
  int x = xpos +(std::min)(sx,ex);

  if (x < 0 || x >= width)return;
  if (x < (std::min)(sx, ex) || x >= (std::max)(sx, ex))return;

  int y = a * x + b;

  if (y < 0 || y >= height)return;
  if (y < (std::min)(sy, ey) || y >= (std::max)(sy, ey))return;


  int pos = y * width + x;

  pimage[pos] = 1;

}

void pbmP1_Write(const char* const fname, const int width, const int height, const unsigned int* const p);

int main()
{

  int width = 100;
  int height = 50;

  int sx = 40; // 線分の始点・終点を設定
  int sy = 5;
  int ex = 10;
  int ey = 20;

  //ループ回数を決定
  int xlen = (std::max)(sx, ex) - (std::min)(sx, ex);

  // xlenを grids×blocksで表現
  int blocks = 512;
  int grids = (xlen + 511) / blocks;

  DimBlock block(blocks, 1);
  DimGrid grid(grids, 1);

  unsigned int* p_cpu = new unsigned int[width*height];
  for (int i = 0; i < width*height; i++) {
    p_cpu[i] = 0;
  }

  fakeCudaKernel(grid, block, line, p_cpu, width, height,sx,sy,ex,ey);


  ///////////////////////////////////////////////
  // 結果出力

  pbmP1_Write("cudaline.pbm", width, height, p_cpu);


  delete[] p_cpu;

}

//! @brief PBM(1byte,テキスト)を書き込む
//! @param [in] fname ファイル名
//! @param [in] width 画像の幅
//! @param [in] height 画像の高さ
//! @param [in] p 画像のメモリへのアドレス
//! @details 1画素1Byteのメモリを渡すと、0,1テキストでファイル名fnameで書き込む
void pbmP1_Write(const char* const fname, const int width, const int height, const unsigned int* const p) { // PPM ASCII

  FILE* fp = fopen(fname, "wb");
  fprintf(fp, "P1\n%d\n%d\n", width, height);

  size_t k = 0;
  for (size_t i = 0; i < (size_t)height; i++) {
    for (size_t j = 0; j < (size_t)width; j++) {
      fprintf(fp, "%d ", p[k] /*? 1 : 0*/);
      k++;
    }
    fprintf(fp, "\n");
  }

  fclose(fp);
}