ひとまず何か書いてコンパイルしてみる。プロジェクト名やソースコード名は適当でいい。ここ（http://tadaoyamaoka.hatenablog.com/entry/2017/05/25/234934）のコードを参考にして以下。ひとつのモジュールに複数の関数を定義することも可能。
#define BOOST_PYTHON_STATIC_LIB
#define BOOST_NUMPY_STATIC_LIB
//上記2行でBoostライブラリを静的にリンクすることを明示
#include
#include
#include
#include
namespace p = boost::python;
namespace np = boost::python::numpy;
char const* hello() {
return "hello, world";
}
char const* hello2() {
return "hello, hell";
}
BOOST_PYTHON_MODULE(test_hoge) {
Py_Initialize();
np::initialize();
p::def("hello", hello);
p::def("hello2", hello2);
}
ビルドが正常にすめば、プロジェクト内の/64/Release ディレクトリにプロジェクト名.dllファイルが生成される。生成されたファイルの拡張子をdllからpydに変更し、名前もtest_hogeに変更してtest_hoge.pydとした上で、pythonの作業ディレクトリにコピー（この管理方法は雑）。Pythonでためしにimportして実行。
>import test_hoge
>test_hoge.hello()
'hello'
というわけでひとまず動いた。
Numpyのndarray型を引数にとって、計算結果もndarray型で返したいわけで、どうすればいいのかね。ひとまずBoostのチュートリアルを覗いてみる。
ndarray型でNxNのゼロ行列を生成する関数をtest_hogeモジュールに追加してみる。
np::ndarray new_zero1(unsigned int N) {
p::tuple shape = p::make_tuple(N,N); //正方行列のサイズ指定
np::dtype dtype = np::dtype::get_builtin(); //要素の組み込み型を指定
return np::zeros(shape, dtype); // ndarrayをnumpy.hppの関数で生成
}
BOOST_PYTHON_MODULE(test_hoge) {
~~
~~
p::def("new_zero1", new_zero1);
ndarrayの扱いを始めるまえに、ndarray型について復習。listと違って、ndarray型は連続したメモリアドレスを持つ。すべての要素は同じ型で、テンソルでいうところの各モードの要素数はすべて同じ。異なるサイズのlistをlist化したもののように不ぞろいなサイズの配列にはできない。加えて、ndarray型は一度宣言してメモリ上に確保されると、サイズを変えるためには全部コピーなりで配列を確保しなおすらしい。
さて、これでndarray型多次元配列を宣言できたわけですが、要素にアクセスしたり線形演算したり、引数にndarray型をとるにはどうしたらよいのか？ ここ（https://qiita.com/termoshtt/items/81eeb0467d9087958f7f）になにやら書いていることが役に立ちそうな感じ。
void mult_two(np::ndarray a) {
int nd = a.get_nd();
if (nd != 1)
throw std::runtime_error("a must be 1-dimensional");
size_t N = a.shape(0);
if (a.get_dtype() != np::dtype::get_builtin())
throw std::runtime_error("a must be float64 array");
double *p = reinterpret_cast(a.get_data());
std::transform(p, p + N, p, [](double x) { return 2 * x; });
//transform()の第1,2引数は始点・終点位置、第3引数は計算結果の格納位置、3つ目が処理。
}
この関数は1次元ndarray型の配列を引数とし、要素の値を2倍にする。配列の次元、サイズ、要素の型はそれぞれget_nd(), shape(), get_dtype()で得られる。get_data()はraw pointerを返しているらしい。ndarray.hppを見ると以下のコメントがあった。というわけでreinterpret_castでdouble型ポインタに。最後の行ではpからp+Nまでの要素に対してラムダ式で書かれた2倍する処理を行っている。というわけで生ポインタとして扱ってデータアクセスができている。宣言されたndarray型の配列が連続したメモリアドレスを持つからできるわけね・・・。
/**
* @brief Return the array's raw data pointer.
*
* This returns char so stride math works properly on it. It's pretty much
* expected that the user will have to reinterpret_cast it.
*/
最後に、多次元のndarray型変数を宣言してアクセス・演算することを考える（あとちょっと）。先ほど参考にしたサイト（ここ（https://qiita.com/termoshtt/items/81eeb0467d9087958f7f）の最後の部分を参考にする。2次元ndarray型変数を受け取ってその要素を順番に表示する関数ぽい。インデックスi,jの要素にアクセスしたいときに、i,j方向のメモリアドレス上の移動量をstride[0], stride[1]として取得し、ポインタを動かしている。strideについてはここ（https://docs.scipy.org/doc/numpy/reference/generated/numpy.ndarray.strides.html#numpy.ndarray.strides）にいろいろ書いてある。
void print(np::ndarray a) {
int nd = a.get_nd();
if (nd != 2)
throw std::runtime_error("a must be two-dimensional");
if (a.get_dtype() != np::dtype::get_builtin())
throw std::runtime_error("a must be float64 array");
auto shape = a.get_shape();
auto strides = a.get_strides();
std::cout << "i j val\n";
for (int i = 0; i < shape[0]; ++i) {
for (int j = 0; j < shape[1]; ++j) {
std::cout << i << " " << j << " "
<< *reinterpret_cast(a.get_data() + i * strides[0] + j * strides[1]) << std::endl;
}
}
}
ここでのi,jはそれぞれ行・列のインデックスに相当している。というわけで、最低限の数値計算はできそうな気配が出てきた。
numpyの2次元のndarray型を引数として受け取ったとき、列ベクトルを取り出したり、列ベクトル同士で線形演算するメソッドが何か用意されているとうれしいのだが・・・。Boost.NumPyのヘッダファイルを見てみるか・・・。
help関数を使うと、モジュール内の関数一覧と各外部仕様を表示できるっぽい。C++で作成したモジュールの場合、C++の戻り型、引数の型とかが表示されるのね。でも結局これだけじゃようわからんわな・・・。
$ help('test_hoge')
