std::ratioについて

C++にはコンパイル時の有理数計算ライブラリがある。個人的な理解では、これは<chrono>のために入ったライブラリで、<chrono>durationを綺麗で速い(使いやすいとは言っていない)やり方で実装するためのものだ。コンパイル時に比率を計算できるので、実行時にかかるオーバーヘッドはせいぜい整数か浮動小数点数の掛け算くらいしか残らないことになる。

ご存知の通り、浮動小数点数は我々の知る実数とは少しだけ異なる。有限のメモリ領域を使ってできる限り広い範囲の数を表現しようとしてはいるが、もちろん完璧ではない。特に、我々は指が10本なので単位系を10の倍数ベースに作っており、そして悲しいことに0.1は2進数では循環少数になる。なので計算途中に丸めが入ることになる。これを避けるためには、整数同士の比で計算するしかない。整数は(実行時に伸びる多倍長整数を使えば)事実上無限の長さを扱えるので、有理数で計算している限り丸め誤差は入る余地がない。まあ、コンパイル時に計算する場合は固定長でなければならないのだが。コンパイル時多倍長演算を実装するという手は置いておいて。

なのでstd::ratioは2つのstd::intmax_tを取る。それぞれが分子と分母に相当する。そして演算もサポートされていて、ratio_addratio_subtractratio_multiplyratio_divideが用意されている。これらを使えば、約分まで済んだstd::ratioが得られる。また、比較も用意されていて、ratio_equalその他がそれぞれ定義されている。

また、主目的が単位換算なので、SI接頭辞(std::kilostd::milliなど)も定義されている。

cpprefjp.github.io

さて、単位には非常に大きな数値が登場するものがあったりする。例えば、分子量は6.02e+23の大きさがあり、原子1個あたりの重さは例えば12 / 6.02e+23グラムになったりする。ところで、log10(264)の値は20に届かない。なのでおそらくこの単位はほとんどの環境でオーバーフローを引き起こしてしまい、上記の枠組みに乗らなくなってしまうだろう。

となると、この際丸め誤差には目をつぶって、コンパイル時に浮動小数点演算を行うしかない。幸いなことにC++11以降ではconstexprがあるので、浮動小数点演算を行うことが可能だ。すると、既存のstd::ratioと相互運用可能なfratioとでも言うべきものを作る必要がある。

面倒なことにC++のnon-type template argumentは浮動小数点数を持てないので、static constexpr double value = /**/;のようにする必要がある。

struct mole
{
    static constexpr double value = 6.02e22;
};

template<typename Numer, typename Denom>
struct fratio
{
    static constexpr double num = Numer::value;
    static constexpr double den = Denom::value;
    static constexpr double value = num / den;
};

しかしながらstd::ratiostd::ratio::valueを持っていないので、共通で使うには::valueを使ってはいけない。とはいえfratioの中で::valueを使わざるを得ない以上(再帰的に定義できた方が便利だ)、共通のインターフェースを用意しておく必要がある。

template<typename T>
struct value_of
{
    static constexpr double value = T::value;
};
// std::ratioへの部分特殊化
template<std::intmax_t N, std::intmax_t D>
struct value_of<std::intmax_t<N, D>>
{
    // 有理数をdoubleに変換する
    static constexpr double value = static_cast<double>(N) / D;
};

これを噛ませれば、std::ratioを共通で使っていくことが可能になる。

正確な数値という利便性を捨ててしまうことになるが、多少の丸め誤差が問題ない場合はこれでことが足りる。もしどうしても正確な数値が必要なら、std::intmax_tのペアもしくは可変長引数を取ることによって多倍長を実現しても構わないのではあるが、使うのがとても難しくなってしまうので、難しいところだ。

Permission denied (publickey) on Travis.CI

タイトルで察した方は帰って結構です。

git submoduleというのがある。これは、別レポジトリをレポジトリの一部として管理できる機能で、要するに依存しているレポジトリのURLと対応するコミットをコードの一部として管理するものだ(大雑把すぎる)。すると依存レポジトリのバージョンをgitで一括で管理できるようになる。

$ git submodule add <repo URL>
# その後
$ git submodule init
$ git submodule update
# または、上記を一括で
$ git submodule update --init
# submodule の submodule まで再帰的に
$ git submodule update --init --recursive

submoduleも一つの独立したgitディレクトリとして振る舞うので、そこでpullしてから上でaddすることによって「依存関係にあるレポジトリのバージョンを上げる」という変更をgitで記録できる。

さて、このようなことを実践したレポジトリをpushしたところ、Travisが落ちた。Permission denied (publickey)と言われている。gitのsshプロトコルでアクセスした場合、sshの鍵がないとダウンロードできない(sshでログイン出来ないので)。もちろん、自分の秘密鍵を公開レポジトリに置くような馬鹿はいない。よって、submoduleをダウンロードできず、CIは失敗する。

さてどうしたものか。 とりあえずググると、Gistに解決策が転がっていた。まず、Travisが自動で行うgit clone --recursiveをやめさせる(git: submodules: false)。続いて、sedgit@github.comhttps://github.comに変換する。力技だ。

Travis-CI submodules · GitHub

確かにこれはLinuxでなら動くが、OS Xだと動かない。OS XUnixでありLinuxではないので、コマンドの挙動がほんの少しずつ異なる。以下に同じ問題に苦しめられた人間の怨嗟の声がある。

Sed: 'sed: 1: invalid command code R' on Mac OS X

何が起きているかというと、sed -iはファイルをその場で(in-place)書き換える。これは怖いので、sedはバックアップファイルを作りたいと思っている。そのバックアップファイルの拡張子を-iオプションの引数として要求されているのだが、我々はそれを渡さず、代わりにsedスクリプトと変更して欲しいファイルを渡した。するとsedスクリプトをバックアップファイルの拡張子と思い込み、続くファイル名をsedスクリプトだと思って実行しようとして、「そんなコマンドはない」と返しているというわけだ。

回避するには、適当な拡張子をくれてやればいい。あるいは、空の文字列を渡せばそもそもバックアップファイルは作られない。

# この""が必要
$ sed -i "" 's/hoge/foo/g' foobar.dat

あとはこれで、Travisで環境を見て分岐してから実行して終わりだ。

アライメント解説

しばらく放置してしまっていた。中々書くべきことが見つからなかったので。Nintendo Switchを購入してしまい、遊んでいたというのもあるが。

今回は、ちょっと普段気にしないアライメントのことについて話してみようと思う。

以前、C++17でのアライメント指定付き動的メモリ確保の話をした時少し話した気がするが、CPUはメモリのどんな位置にでも好き勝手にアクセスできるわけではない。例えば、一回のメモリアクセスでは16バイト分のデータを16の倍数のアドレスからしか読み込めない、というような感じの制限がある(数字は適当)。なので、以下のようになって困ることがある。

メモリアクセス
 v v v v
| | |d|a|t|a| | |

というわけで、メモリ上でのデータの配置には守るべきルールが存在し、型Talignof(T)(例えば8バイト)の倍数にあたるメモリアドレスから始まる領域にしか置くことができない。

これは、以下のようなコードを書くときに問題になってくる。

  1. バイナリファイルを読み書きする場合
  2. 何らかの理由でデータのビット表現を取得する必要があるとき
  3. 普通でないアライメント指定をして動的メモリ確保するとき

他にもあるかもしれない。何にせよ、例えばバイナリファイルを読み書きしている時、以下のようなことをすると一発で地雷を踏みぬいてしまう。

const char* stream;
const double d = *reinterpret_cast<const double*>(stream);

ここで、streamの位置はバイト単位で変わる。例えばデータ型のタグとして1バイト使っていたとすると、その次からdoubleのデータが始まったりするだろう。すると、メモリ上に確かにsizeof(double)分の領域があり、そこにdoubleのビット列が入っているのだが、メモリアドレスがalignof(double)の倍数になっていないという状況になりかねない。その場合、上のコードでデリファレンスしているポインタは、アライメント要求を満たさない不適切なポインタになってしまう。 つまりこうだ。

0x0000
 v
| |d|a|t|a| | |
   ^ 0x0001 (8の倍数でないため、`double`のポインタにできない)

ではどうするか。C規格では(unsigned|signed) charとならどんな型でも変換可能ということになっている。C++17ではstd::byteもここに仲間入りした。これは、char*からT*への変換が常に成立するという意味ではなく(これは常には成立しない)、T*charの配列へのポインタへ再解釈して構わないということである。charのアライメント要求が最小ということなのだろう。 なので、より安全なコードは以下のようになる。

const char* stream;
double d;
std::memcpy(reinterpret_cast<void*>(std::addressof(d)),
            reinterpret_cast<void*>(stream), sizeof(double));

上記コードではdoubleへのポインタを(memcpyは内部でunsigned charへのポインタに変換するので)unsigned charへのポインタに変換し、両方をunsigned char[N]だと思ってデータをコピーしている。 まったく、面倒この上ない。

他に、データのビット表現を取得するとき、というのは、例えばfloatの値をuint32_tへビットをそのままに変換するという意味だ。なぜそんなことをするのかというと、謎のビット演算によっていくつかの高コストな浮動小数点演算が近似できるからだ。このテクニックは以下によくまとまっている。

GitHub - keon/awesome-bits: A curated list of awesome bitwise operations and tricks

この場合、以下のようにしたくなる。

float f = /**/;
std::uint32_t i = *reinterpret_cast<std::uint32_t*>(&f);

実際、以前このブログでも書いてしまった気がする。だがこれは通らない(大抵動くのだが、間違ったコードなので動く保証はない)。これは、floatstd::uint32_tのアライメント要求が一致している前提で書かれたコードだが、そんな保証はない。保証されているのは前述の通りcharstd::byteへの変換のみである。

なので、正しくは以下だ。

float f = /**/
std::uint32_t i;
std::memcpy(reinterpret_cast<void*>(std::addressof(f)),
            reinterpret_cast<void*>(std::addressof(i)),
            sizeof(float));

よく考えると、floatが4バイトなのは定義されているのだろうか? std::uint32_tはその点便利で、これは2の補数表現でピッタリ32ビットであることが保証されている。そのような型をサポートしていないアーキテクチャでは、この型が定義されない。なので知らずにコンパイルしてもエラーになってくれる。

ちょっと見てみたが、IEEE754準拠は必須ではない(__STDC_IEC_559__が1ならIEEE 754に対応しているとわかる)。ということは、floatの中身については何も仮定できない。まあintがそうなのだから当たり前か。

(ため息)

さて、上記の3つの状況の最後のものは、普通でないアライメントを指定して何かする時だ。

まず、なぜそんなことをする必要があるのかというと、これもいくつか理由がある。

  1. SIMDなどの特殊な命令を使う
  2. パフォーマンス最適化

まず、SIMDとはSingle Instruction Multiple Dataの略で、一つの命令を複数のデータに同時に適用することで速度向上を図るものだ。例えば、レジスタに4つの数値をロードしておいて、その全てに加算命令を発行して同時に処理すれば、1命令で4つ分の計算ができ、速度が理論上4倍になる。このとき、ロードするデータが通常より大きなアライメントを満たしていれば、ロード速度が向上する。そうでなければ複数回のメモリアクセスが発生し、少し時間がかかってしまう。

他に、キャッシュの更新を抑える目的でアライメントを調整することがある。根本に立ち戻ると、現代のアーキテクチャではメモリはCPUに比べると遅い。なので、CPUの計算よりもメモリアクセスが律速になってしまうことが多い。それをなんとか隠蔽しようとして、CPUは内部に高速にアクセスできるが容量の少ないメモリを持っておき、メモリ上の使いそうなものをそこに先に持ってきておくようになった。CPUは一回につきある決まった量のデータをキャッシュしておくのだが、これはキャッシュラインと呼ばれている。

賢い戦略ではあるのだが、やはり面倒なことはある。例えば並列化で複数のCPUが同じメモリ領域を見ているとすると、一つのCPUでのキャッシュへの書き込みが他のCPUのキャッシュに同期されなければ、データがおかしくなってしまう。なので書き込みが生じた場合、全CPUとメインメモリでデータを同期する必要がある(キャッシュコヒーレンシ)。さて、もし一つのキャッシュライン上に凄まじく頻繁に更新される値が一つと、全く更新されないとわかっている値が沢山乗っていたらどうだろう。ほとんどの値は変更されないのに、同じキャッシュラインに乗っているたった一つの値のせいでキャッシュライン全体の動機が発生する(false sharing)。これを回避する方法はいくつかある。キャッシュラインサイズに沿ったアライメント指定をして別のラインに乗るようにするものと、構造体にダミーデータでパディングして別のキャッシュラインに入るようにするというものだ。

他に外部デバイスと通信するときに気をつけないといけないという話を聞いたことがあるが、あまり書いたことがないのでよくわかっていない。GPUのテクスチャメモリ(か、別の特殊なメモリだったかもしれない)だと何か要求されていたような気がする。

さて、上記のようなコードを書く場合は同時にアライメントについて調べると思うので、あまり変な事にはならないのではないかと思う。ただ、落とし穴があるのは動的にメモリ確保をする時だ。通常、mallocnewはメモリ領域のサイズしか知らされず、アライメントまでは指定されない。mallocにアライメント値を渡した記憶のある方はいるだろうか? 受け取らないのだからいるはずがない。ではどうしているのかというと、デフォルトのアライメントが決まっており、それに合わせたメモリ領域が返ってくるのだ。このデフォルト値は大抵、基本型の要求するアライメントのうち最大のものになっている。アライメントは何かの倍数のアドレスに置くことを要求するものなので、最大の値の倍数(最大公倍数と言えばいいのだが、大抵全て2のN乗なので最大の値に等しい)は他のものを満たすという便利な性質がある。ではデフォルト値より大きな値を設定するとどうなるか? どうにもならない。mallocnewもそんなことは知らない。なのでアライメントは単に無視されてしまう。

このために、C11ではaligned_allocが、C++17ではアライメント指定されたデータの動的メモリ確保が入った。これらはアライメント要求を受け取り、それに適した領域をアロケートして返す。これ以前だと、大きめにアロケートして先頭ポインタが要求を満たすところまでずらすか、posix_memalign_aligned_mallocなどの処理系定義関数を使う必要があった。

面倒この上ない。だが、これらは計算機科学の奥底に眠っている問題ではなく、バイナリファイルを読もうとしたりちょっとした最適化をしようとしたときにすぐに首をもたげる問題なので、意識しておくのもいいかもしれない。

開発ブログの必要性

何かを開発していると、何かの目的(ランタイムパフォーマンスなど)で少しわかりにくいハックをする必要が出てくることがある。もちろんそういう場所ではコメントを入れるわけだが、長くなりすぎたり、一つのソースファイルから複数のソースファイルへ言及する羽目になったりして、むしろ見た人を面食らわせることになったりもする。

例えば、ある値を先に計算しておくことで、1.計算コストを減らす、2.キャッシュ局所性を(結果的に)上げるなどが達成できるとしよう。だがそれを実装するには、先に計算した値を入れておくコンテナ、それを使って実際に計算するクラス、そのクラスとコンテナを管理するクラスの概ね3つに影響が出る。もちろん設計によっては1つのクラスで済むかもしれないし、もういくつかのレイヤーが登場するかもしれない。となると、似たようなコメントを複数箇所に入れるのは冗長な気もするし、冗長でなくしようとすると説明不足になりそうな気もする。

そういうときに全体の流れを書くことのできるdevlogがあれば、便利なのかも知れない。関連技術の普及にも繋がるかも知れないし。

今日、外に出ると晩夏の匂いがした。湿気と、草木と、熱気の残り香。夏が終わりかけていることと、それでもまだ夏であるということが、アパートの廊下の外を見るまでもなくわかった。外を歩きながら、今年の夏が殆ど夏らしいことができないまま終わろうとしていることにちょっとした焦りを覚えつつ、これまでの夏にした「夏らしいこと」について思いを馳せていた。

数年前、瀬戸内海のある島に旅行したことを思い出す。とても暑い日だった。一周するのに徒歩でも数時間という小さい島だったというのに、てっぺんまで登っただけで汗だくになってしまっていた。非常に理想的な形をした島だったので、展望台からは360度海が見え、対岸の中国地方の沿岸部や、遠くにある他の島、その間を行き来する船などもよく見えた。その時のことはよく覚えている。

当時、簡単な物理シミュレーションが書けるようになった頃で、なので自分の持っているそこそこのスペックのノートPCでどのくらいの規模の系を計算できるかがある程度わかってきていた。眼下に広がる海の規模が、乗ってきた船の規模を足がかりに把握できる。その表面に立つさざなみが、それに反射してきらめく光が、どれほどの規模なのかがわかってくる。そのスケール感についての直感が突如湧き上がってきて、しばし圧倒されていた。たまにしかないことだが、本気で、心の底から自然界に対しての畏敬の念を覚えていた。

思えば、その光景を私が認識するまでのパスウェイの方が小さいながらもよっぽど複雑で、そちらの方がある意味では驚嘆すべきことだったのかもしれない。しかし、そういった諸々を吹き飛ばすだけの物量がその景色には含意されていた。単純な距離や広がりを超えた、世界の本当の意味での広さのようなものが感じられた気がした。もちろん人類はまだミクロスケールの全てを解明したわけではないし(もちろん私は現時点での最先端は理解していない)、そもそもそこに単純な答えがあるかどうかすらわからないのだが、人類が、いや私がある程度自信を持っているレイヤーだけから見ても、世界は凄まじく広大だ。

その後も、そんな気持ちになることがある。山に登ったり海に行ったりしなくても、単に川辺を散歩している時にも、それどころかベンチでふと上を見上げただけでも、この感覚が襲ってくる。光子が大量に飛んできて、膨大な数の水分子にぶつかり、たまに電子を励起させたりしつつ私の目まで飛んでくる。「複雑ではない、ただ数が多いだけ」と言ったのはFeynmanだったか。数が多いだけでも、人を圧倒することはある。

自然を理解することで自然への畏敬の念が減ずる、というのはよく反駁されるべき言説として登場するし、私もこれは間違っていると思う。理解し、作ろうとして始めて、自分が何に相対していたのかがわかってくる。実際、あれほどの恐ろしさを感じたことはそれまでの人生でなかった。当たり前に過ごしていた日常の裏側で何が起きていたか、という認識の転換。人生でそう何回もできないだろう経験の一つが、あの瞬間だったろう。他にも何かの現象の感覚を掴んでものの見方が変わったことはあるが、その中でもあれは特異な例だったように思う。かなり初等的な物理しか想定せずとも、認識はああも変わるのだ。

今日もコップの中の珈琲に、木の葉の影の中に、浮かぶ雲にそれを見る。まだ手の届かない世界にため息をつきながら。

Expression templateとfmaについて

Fused Multiply Addという命令がある。これは、a * b + cという形の演算を1命令で処理するものだ。丸めが一度しか走らないので精度がよくなる。以前、理論値と数値計算を比較してテストするときにこいつが使われるかどうかでテストが通ったり落ちたりして困った思い出がある。

さて、必ずこれを使うようにするにはどうしたらよいか。簡単だ。C99ではmath.hfmafmaffmalが定義されており、この関数を使えばよい。C規格(N1256)の§F.9.10.1 "The fma functions"で、

fma(x,y,z) computes xy + z, correctly rounded once.

とのことなので、CPUがサポートしていなかったらソフトウェア側でエミュレーションするのだろう。実際FP_FAST_FMAというマクロがあり(§7.12 Mathematics <math.h>)、これは積と和を計算するよりもFMAの方が早い場合のみ定義される。FMAFとFMAL版もある。これはほぼ専用命令があれば定義されることと同義なので、これが定義されているかどうかで使い分けもできるだろう。もしCPUがFPUを持っておらず、浮動小数点数演算を全てソフトウェアエミュレーションしていれば(そしてエミュレーションの実装でFMAの方が速ければ)、FMA専用命令がなくとも定義されるかもしれない。いまどきそんなCPUあるか?

さて、C++C++11からstd::fmaがサポートされた。これはCのfmaと同等である。また、SIMDFMAでは以下のようなイントリンシックがある。

__m128 _mm_fmadd_ps(__m128 a, __m128 b, __m128 c);

FMA関係のイントリンシックはもちろんこの一つではなく、幅(__m128, __m256) かける 精度(_ps_pd) 足す 非SIMD_ss_sd) で10通りある。AVX512ではマスクがあるのでもっと色々ある。

さて、これを制御したいとする。ライブラリ側で、ユーザーには気づかれないうちに。

目的からして、ユーザーがどんな演算をするか何も仮定できない。とすると、どういう演算をするかという情報を取り出せる場所は一つしかない。ユーザーのコードだ。

先に書いておいたライブラリがユーザーのコードから情報を取り出す方法としては、まずExpression Template(ET)が思いつく。これは本来はBlitz++という線形代数ライブラリが必要のない要素の計算を避けるために開発した技法で、一種の遅延評価である。実装をかいつまんで書くと、まず以下のようなクラスを作って、演算子が計算結果の一次オブジェクトでなく、両辺への参照を持った特殊なクラスを返すようにする。

template<typename L, typename R>
struct added
{
    using value_type = decltype(
        std::declval<typename L::value_type>() +
        std::declval<typename R::value_type>()
    );
    value_type operator[](std::size_t i) const noexcept
    {
        return l[i] + r[i]};
    }

    L const& l;
    R const& r;
};

このクラスは2つの型(左辺のオブジェクトと右辺のオブジェクトになる)をとって、アクセスしたときに必要な値を計算して返すクラスだ。このクラスはアクセスされるまでは計算を行わない。よって、例えば巨大な行列の積を計算するが対角要素以外いらないとか、そういう時のために使われてきた。

さて、これを使うと、以下の式のdの型はどうなるだろうか。それぞれの変数はmatrix型だとする。

auto d = a * b + c;

まずa*bが行われるので、以下のようになる。

auto d = multiplied<matrix, matrix> + c;

その後、足し算が行われるので、以下のようになる。

auto d = added<multiplied<matrix, matrix>, matrix>;

というわけで、dの型はadded<multiplied<matrix, matrix>, matrix>だ。これは木構造になっている。というか、抽象構文木(AST)になっているのだ。

本来は実装の手間を省くために、std::plusなどを取る型パラメータを付け加えるためよりASTっぽくなるが、今回は例が複雑になるためやめた。

このように、ETを使うことでユーザーコードで出てくる変数の型をASTにすることができる。ここで、templateは部分特殊化、つまりパターンマッチが使えるということを思い出そう。すると、以下のような操作ができる。

template<typename L, typename M, typename R>
struct added<multiplied<L, M>, R>
{
    using value_type = decltype(
        std::declval<typename L::value_type>() *
        std::declval<typename M::value_type>() +
        std::declval<typename R::value_type>()
    );
    value_type operator[](std::size_t i) const noexcept
    {
        return std::fma(l[i], m[i], r[i]};
    }

    L const& l;
    M const& m;
    R const& r;
}:

これは、addedの左辺が何かしらのmultipliedだった場合への特殊化だ。そういうクラスが現れた時は、通常のaddedmultipliedではなくこちらで実体化される。

――という風にすれば、ユーザーが意識せずに書いてもfmaが勝手に生成されるという風にできるなあ、と思いつつも実装していない。面倒なのである。

だが、ライブラリ側でユーザーコードのASTを変換できる可能性というのは面白いのではないかなあと思う。コンパイラが最適化のためにいじっているのも何らかの中間表現なので、理論上は同等のものが実装できるはずだ。最近はETはコンパイラの自動SIMD化を阻害する(解析が難しくなる)ので下火という話も聞くが、こういう方向の発展もありなのではないか。

……それはコンパイラの最適化を阻害しつつ、自分でコンパイラがやってるのと同じ最適化をするってことだから、めっちゃくちゃ特殊な場合を除いて意味ないんじゃないの。

maveの中身について

前回SIMDベクトルライブラリを作った話をしたが、そこでやっていることについてかいつまんで書いておく。

まず、以下のようなクラスがある。

namespace mave
{
template<typename T, std::size_t R, std::size_t C>
class matrix;

template<typename T, std::size_t N>
using vector = matrix<T, N, 1>;
}

これについて、ごく普通の演算子の他に、以下のような演算子オーバーロードを行った。

namespace mave
{
template<typename T, std::size_t N>
std::tuple<vector<T, N>, vector<T, N>>
operator+(std::tuple<const vector<T, N>&, const vector<T, N>&> lhs,
          std::tuple<const vector<T, N>&, const vector<T, N>&> rhs)

template<typename T, std::size_t N>
std::tuple<vector<T, N>, vector<T, N>>
operator+=(std::tuple<      vector<T, N>&,       vector<T, N>&> lhs,
           std::tuple<const vector<T, N>&, const vector<T, N>&> rhs)
}

そう、あまり見る機会はないが、実はoperator+=もクラス外で定義できる。通常この演算子はメンバへの書き込みを伴うのでクラス内でメンバとして定義する(ことが多い)。しかし今回は、std::tupleにはメンバを追加できないのでfreeな関数として実装している。まあこの場合はコンストラクタで与えられるものが全てなので、メンバ関数として実装する必要がそもそもないのだが。

念頭にあるのは以下のようなユーザーコードである。

const mave::vector<float, 3> v1(/*...*/), v2(/*...*/);
const mave::vector<float, 3> w1(/*...*/), w2(/*...*/);
auto [u1, u2] = std::tie(v1, v2) + std::tie(w1, w2);

まず、std::tieは参照を格納したstd::tupleを作って返す。これは本来、タプルのパターンマッチをライブラリレベルでサポートするための機能だ。普通は以下のように使う。

std::tuple<int, double, std::string> f();

int main()
{
    int i;
    double d;
    std::string s;
    std::tie(i, d, s) = f();
}

これは、多値を返す関数の戻り値を分解しつつ受け取るコードだ。何が起きるかというと、まずstd::tieによってstd::tuple<int&, double&, std::string&>が作られ、そこにstd::tuple<int, double, std::string>が代入される。参照への代入なのでもともとのオブジェクトへの書き込みが行われ、idsのそれぞれが書き換わるという寸法だ。その後std::tieによって作成されたオブジェクト(のprvalue)の寿命は人知れず切れる。

ついでにこれはstd::ignoreを使うことで一部の値を無視することができる。

int i;
double d;
std::tie(i, d, std::ignore) = f();

このstd::ignoreは処理系定義なクラスのインスタンスと定義されている。代入演算子template化してどんな型でも代入できるようにしておき、実際には代入演算子は何もしないようにすれば実装できるだろう。こんな感じで。

struct ignore_t
{
    template<typename T>
    ignore_t& operator=(const T&) noexcept {return *this;}
};

これらの機能は便利だったのだが、std::tieだと初期化にコストのかかるオブジェクトを受け取りたい時に不必要なコストがかかってしまうことがある(先にデフォルトコンストラクタが走り、その後std::tieで受け取るときムーブコンストラクタが走る)。また、そもそもデフォルトコンストラクタを持っていないとstd::tieは使いにくい。適当な引数を与えてあとで代入することで使えないわけではないが、明らかに意味不明だろう。というわけで、結局コア言語側で構造化束縛が導入されることになった。

閑話休題std::tieは参照のタプルを作るので、参照のタプルを受け取って使うことにすればよい、というのは使い方を決めると自動的に決まる。要素数が2個の場合はstd::pairを使えるようにしようかと一瞬思ったが、微妙に複雑になってしまう(結局tupleもサポートしないといけないので転送するだけの関数が増えてしまう)のでやめた。

さて、これで大丈夫だろう、と思ったのだが、少し困ったことが起きた。先ほどの状態でこれは通ったのだが、

const mave::vector<float, 3> v1(/*...*/), v2(/*...*/);
const mave::vector<float, 3> w1(/*...*/), w2(/*...*/);
auto [u1, u2] = std::tie(v1, v2) + std::tie(w1, w2);

以下は通らなかった。

mave::vector<float, 3> v1(/*...*/), v2(/*...*/);
mave::vector<float, 3> w1(/*...*/), w2(/*...*/);
auto [u1, u2] = std::tie(v1, v2) + std::tie(w1, w2);

例えば、以下のようなtemplate関数があったとする。

template<typename T1, typename T2>
void f(std::tuple<T1 const&, T2 const&> tied)
{
    std::cout << std::get<0>(tied) << ", "
              << std::get<1>(tied) << std::endl;
}

これにstd::tuple<int&, double&>などを渡すと、const T1intがマッチしないのでこれは候補にならない。なので、以下は通らない。

int i = 42;
double pi = 3.14;
f(std::tie(i, pi)); // no matching function!

std::tie(T1, T2)std::tuple<T1&, T2&>を返すからだ。だが、もちろんT&const T&に変換できるので、templateをやめると通る。

void f(std::tuple<int const&, double const&> tied)
{
    std::cout << std::get<0>(tied) << ", "
              << std::get<1>(tied) << std::endl;
}

int main()
{
    int i = 42;
    double pi = 3.14;
    f(std::tie(i, pi)); // conversion from int& to int const&
    return 0;
}

これが通ってtemplateにすると通らないのはハマりどころな気がする。実際、書いてエラーを見るまでtemplateでも通ると思っていた。一度、呼び出す関数の解決順序を見直した方が良さそうだ。

さて、これをどうするかだ。多分templateをやめてmave::vector<float, 3><double,3>matrix<float, 3, 3>……とそれぞれに演算子オーバーロードを用意すれば良いのだが、その場合用意していないサイズのmatrixは、作れるが演算だけできないというよくわからないことになる。全てのサイズ、型についてオーバーロードを用意することはできない。そういう問題を回避する為にtemplateがあるのだ。では非templatetemplate演算子を共存させれば、と考えるわけだが、当然template版と非template版の両方がマッチするので、conflictしてオーバーロード解決が失敗する。

解決策として、こうした。

#define MAVE_GENERATE_FORWARDING_BINARY_FUNCTIONS_MATRIX_MATRIX_MATRIX(FUNC_NAME, L_MODIFICATION, R_MODIFICATION)\
  //...
MAVE_GENERATE_FORWARDING_BINARY_FUNCTIONS_MATRIX_MATRIX_MATRIX(operator+, MAVE_EMPTY_ARGUMENT, MAVE_EMPTY_ARGUMENT)
MAVE_GENERATE_FORWARDING_BINARY_FUNCTIONS_MATRIX_MATRIX_MATRIX(operator+, const&             , MAVE_EMPTY_ARGUMENT)
MAVE_GENERATE_FORWARDING_BINARY_FUNCTIONS_MATRIX_MATRIX_MATRIX(operator+, &                  , MAVE_EMPTY_ARGUMENT)
MAVE_GENERATE_FORWARDING_BINARY_FUNCTIONS_MATRIX_MATRIX_MATRIX(operator+, MAVE_EMPTY_ARGUMENT, &)
MAVE_GENERATE_FORWARDING_BINARY_FUNCTIONS_MATRIX_MATRIX_MATRIX(operator+, const&             , &)
MAVE_GENERATE_FORWARDING_BINARY_FUNCTIONS_MATRIX_MATRIX_MATRIX(operator+, &                  , &)
MAVE_GENERATE_FORWARDING_BINARY_FUNCTIONS_MATRIX_MATRIX_MATRIX(operator+, MAVE_EMPTY_ARGUMENT, const&)
MAVE_GENERATE_FORWARDING_BINARY_FUNCTIONS_MATRIX_MATRIX_MATRIX(operator+, &                  , const&)

グワーッ!

渡されうる引数はconst&&(値渡し)の3通りがありえて、2引数取る場合はその組み合わせで増えるので大変なことになる。 これ、何かいい解決法がないだろうか。考えているのだが、浮かばなくて困っている。