boost::multi_array で快適多次元配列生活

これは C++ Advent Calendar 2017 の17日目の投稿です。

なぜ boost::mutli_array

普段は数値計算をぼちぼちとやっている人です。

数値計算では多次元配列（特に2次元）を使う事が多く、普段は Eigen を使ってコードを書いていました。しかし、コードの中で3次元配列が必要になり、さてこれはどうしようかということで boost::mutli_array を使ってみる事にしました。

実際に使ってみると結構便利な感じで、趣味のコーディングでも良い感じに使えそうです。というわけで、布教がてらアドベントカレンダーの記事にしてみました。日本語化されたドキュメントも見れます。翻訳感謝！
https://boostjp.github.io/archive/boost_docs/libs/multi_array/user.html

今までのC++での多次元配列

C++で多次元配列を実現する方法は多岐に渡りますが、

const int N1 = 100;
const int N2 = 200;
double **my_array = new double*[N1];
for(int i = 0; i < N1; ++i) {
  my_array[i] = new double[N2];
}
# ...
for(int i = 0; i < N1; ++i) {
  delete[] my_array[i];
}
delete[] my_array;

と言った書き方は new, delete を意識しなければならないですし、流石に前時代的な書き方です。最近のC++ではほとんどnewをしないでプログラムを書くのが当たり前ですからね。多くの人は、以下のようなコードを書いているのではないでしょうか。

const int N1 = 100;
const int N2 = 200;
std::vector<std::vector<double>> my_array(N1);
for(int i = 0; i < N1; ++i) {
  my_array[i].resize(N2);
}

これで delete をコードから排除出来ました。しかし multi_array[i].size() が i によらない場合は少し冗長に感じますし、型名がちょっと長すぎます。何よりもこの方法では multi_array 全体に渡ってメモリ空間が連続である保証がありません。なので、外部ライブラリで行列計算をしたい場合や、マシン性能をフルに活用したい場合はこのパターンを用いるのは厳しいでしょう。二次元配列ならば Eigen などを用いれば便利ですが、3次元配列以上になったら困りますよね。

boost::multi_array の基本的な使い方

というわけで、こういう場合は boost::multi_array を使いましょう。

const int N1 = 100;
const int N2 = 200;
boost::multi_array<double, 2> my_array(boost::extents[N1][N2])

とっても簡単ですね。テンプレート引数には、型と多次元配列の次元を指定します。配列の大きさは、boost::extents で与えます。
次に、配列をループしてみましょう。

for(int i = 0; i < my_array.shape()[0]; ++i) {
  for(int j = 0; j < my_array.shape()[1]; ++j) {
    std::cout << my_array[i][j] << std::endl;
  }
}

このように shape を用いてサイズを知ることができるので、便利です。ところで、インデックスの型に何も考えずに int を使っていますが、実は long long でなければならないほど、配列サイズが大きいかもしれません。そういう場合にも対応できるように、インデックスの型は boost::multi_array<double, 2>::index を用いた方が良いでしょう。

using ma_type = boost::multi_array<double, 2>;
using ma_index = ma_type::index;
ma_type my_array(boost::extents[3'000'000'000ll][2]);
for(ma_index i = 0; i < my_array.shape()[0]; ++i) {
  for(ma_index j = 0; j < my_array.shape()[1]; ++j) {
    my_array[i][j] = i + j;
  }
}

配列サイズを変える

関数に与えた boost::multi_array の配列サイズを変更したい場合があるかもしれません。そのような場合は resize を利用しましょう。

void resize_and_hoge(int n, boost::multi_array<double, 2> my_array) {
  int s = n * n + 1;
  decltype(my_array)::extent_gen extents;
  my_array.resize(extents[s][2]);
  for(decltype(my_array)::index i = 0; i < my_array.shape()[0]; ++i) {
    my_array[i][0] = i * 2 + 0;
    my_array[i][1] = i * 2 + 1;
  }
}

このサンプルでは使っていませんが、地味に凄いことに、以前に代入していた配列要素にもアクセスすることが出来ます（勿論配列サイズが縮小されて、消えてしまった領域はアクセスできなくなります）。そのために、内部ではコピーを行っていたりするので、多用しないほうがいいでしょうが、使いどころによってはかなり便利だと思います。

もし、単に配列の解釈の仕方を変えたい場合は reshape が便利です。

boost::multi_array<int, 2> A(boost::extents[2][3]);
boost::array<int, 2> dims = {{3, 2}};       
A.reshape(dims);

ただし、全体の要素数が変化しないように注意しましょう。

Fortran like なメモリ配置

C++ではCと同様に一番左側の配列アクセスが連続になるようにメモリに配置されますが、これを Fortran のように、右側の配列アクセスが連続になるような形に変えることが出来ます。これは例を示した方が分かりやすいでしょう。

boost::multi_array<int, 2> A(boost::extents[2][2]);
A[0][0] = 1;
A[0][1] = 2;
A[1][0] = 3;
A[1][1] = 4;
for(int i = 0; i < 4; ++i) {
  std::cout << *(A.data() + i) << ", ";
}
std::cout << std::endl;
# 1, 2, 3, 4, と表示されます。

boost::multi_array<int, 2> B(boost::extents[2][2], boost::fortran_storage_order());
B[0][0] = 1;
B[0][1] = 2;
B[1][0] = 3;
B[1][1] = 4;
for(int i = 0; i < 4; ++i) {
  std::cout << *(B.data() + i) << ", ";
}
std::cout << std::endl;
# 1, 3, 2, 4, と表示されます。

一部を切り出す

いわゆるスライスとかビューとか言われる機能が使えます。これも例を示したほうが早いでしょう。

using range = boost::multi_array_types::index_range;
using array_type = boost::multi_array<int, 3>;
array_type my_array(boost::extents[2][3][4]);
for(int i = 0; i < 2 * 3 * 4; ++i) {
  *(my_array.data() + i) = i;
}
array_type::array_view<3>::type my_view =
    my_array[ boost::indices[range(0,2)][range(1,3)][range(0,4,2)] ];

for(int i = 0; i != 2; ++i)
  for(int j = 0; j != 2; ++j)
    for(int k = 0; k != 2; ++k) {
      assert(my_view[i][j][k] == my_array[i][j+1][k*2]);
      my_view[i][j][k] = 0;
    }
for(int i = 0; i < 2 * 3 * 4; ++i) {
  std::cout << *(my_array.data() + i) << ", ";
}
std::cout << std::endl;

上記の例では、1つ目の添字についてはそのまま [0, 2) の範囲、2つ目の添字については [1,3) の範囲に、3つめの添字については [0, 4) の範囲で2個飛ばしでアクセスするビューを作成しています。また、ビューとか言っているのですが、実際のところ上記の例のように書き込みが出来ます。ビューで書き換えている場所が、0になっていることを確認できるはずです。

これからの多次元配列の動向

なんでこんな基本的な機能が標準ライブラリにないんだと思っていたのですが、どうやら将来的には標準ライブラリに追加されるかもしれません。
現在、 boost::multi_array とほぼ同等の機能が提供されると思われる mdspan が策定中のようです。もし、自分がターゲットにしている環境で問題なく動作するようならばこちらの使用も検討してみてはいかがでしょうか。

まとめ

多次元配列を使いたいなら boost::multi_array を使えば直感的にコードが書けます。以下のようなメリットがあります。

• 普通の多次元配列と同じ使用感
• 連続的なメモリ領域の保証
• 関数に多次元配列を渡すのが簡単（配列サイズなどを別途渡さずに済む）
• サイズの変更可能
• ビューの作成可能
• Fortran を意識した関数に渡す配列を扱うのも簡単

以上のようなメリットがあるので、是非 boost::multi_array を使ってみてください。

所感

Boost 関係の記事なので、C++ アドベントカレンダーの内容として適切かわかりませんが、何でもありらしいので、きっと大丈夫でしょう。
今回の記事を書いてる途中で、冒頭で紹介したドキュメントの Sample にバグが見つかったりしました（修正済み）。それだけでも記事を書いた意義があったかなという気分です。

最後に、ここまで読んでくれたあなたに感謝を。