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

これは [C++ Advent Calendar 2017](https://qiita.com/advent-calendar/2017/cpp) の17日目の投稿です。

# なぜ boost::mutli_array

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

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

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

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


```cpp
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をしないでプログラムを書くのが当たり前ですからね。多くの人は、以下のようなコードを書いているのではないでしょうか。

```cpp
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 を使いましょう。

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

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

```cpp
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``` を用いた方が良いでしょう。

```cpp
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``` を利用しましょう。

```cpp
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``` が便利です。

```cpp
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 のように、右側の配列アクセスが連続になるような形に変えることが出来ます。これは例を示した方が分かりやすいでしょう。

```cpp
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, と表示されます。
```

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

```cpp
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 にバグが見つかったりしました（修正済み）。それだけでも記事を書いた意義があったかなという気分です。

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