座標ベクトルと汎化について

私の書くプログラムでは、趣味でも仕事でも、低次元(大抵は2か3次元)でのベクトルを扱うことが多い。そうなると普通は座標ベクトルを意味する構造体を作って管理するだろう。

基本的に、ベクトルとしての基本的な演算、ベクトル同士の足し引き、スカラー倍や長さの計算などはスカラーの型が何であるか、また次元がどの程度かということとは関係がない。 なのでtemplateにした方がよい。この程度ならそんなに手間が増えることもないし、N次元ベクトルやfloat型等について変な接頭辞や接尾辞をつけながら逐一コピペするよりずっと良い。

template<typename T, std::size_t N>
struct coordinate
{
    typedef T value_type;

    coordinate& operator+=(const coordinate<T, N>& rhs) noexcept;
    coordinate& operator-=(const coordinate<T, N>& rhs) noexcept;
    coordinate& operator*=(const value_type& rhs) noexcept;
    coordinate& operator/=(const value_type& rhs) noexcept;

    value_type& operator[](std::size_t i)       noexcept;
    value_type  operator[](std::size_t i) const noexcept;
    value_type& at(std::size_t i);
    value_type  at(std::size_t i) const;

    std::array<T, N> values_;
};

template<typename T, std::size_t N>
coordinate<T, N>
operator+(const coodinate<T, N>& lhs, const coordinate<T, N>& rhs) noexcept;
//...

だが、最初から3次元しか必要ないと決まっている、或はこのプロジェクトでは永遠に2次元と3次元しか使わないので丸ままでコピペしても余裕で管理できるなら、例えば以下のように実装しても良いはずだ。

template<typename T>
struct xyz
{
    typedef T value_type;
    value_type x, y, z;
};

こうすることのデメリットは上で述べた通り(3次元と決め打ちしている)だが、メリットがないわけではない。 読むときにpos[0]でアクセスするよりも、pos.xでアクセスする方がわかりやすいという人は少なからずいるはずだ。 というか、それによってよりわかり難くなるということはないだろう。 それに、これだとコンストラクタを書くのが非常に自明になる。

template<typename T>
struct xyz
{
    xyz(T x_, T y_, T z_) noexcept : x(x_), y(y_), z(z_){}
};

わかりやすい。もし同様のものを先に述べたクラスについて作ろうと思えば、先のクラスは次元の数に対して一般化されているので、N次元のどれででも同様の機能性を持たせないと筋が悪くなってしまう。一般化の弊害だ。 解答はいくつかある。一つは、3次元の場合について特殊化すること。だがこれはコピペと同じだ。必要なものを再定義することになるのだから、それをするくらいなら最初から決め打ちすればいい。 N次元のベクトルという文脈に沿うような、より適切な解答は、しかし少々技巧的になってしまう。一つ目の、比較的わかりやすい解答は、std::initializer_listだろう。

template<typename T, std::size_t N>
struct coordinate
{
    typedef T value_type;

    coordinate(std::initializer_list<T> il)
    {
        if(il.size() != N) {throw std::invalid_argument("invalid size");}
        std::copy(il.begin(), il.end(), values_.begin());
    }
};

しかしこれだと、サイズがおかしかったときには実行時にエラーを吐くことになる。引数の個数はコンパイル時に確定しているのに、実行するまでその問題がわからないとは! あまり喜ばしいことではない。折角コンパイラを使うのだから、問題は全て見つけてほしいではないか。次元数はコンパイル時にわかっているというのに。

引数の個数を間違えるとコンパイルエラーにするなら、少々面倒なことになる。

namespace meta
{
template<typename ... Ts>             struct and_;
template<typename T>                  struct and_<T>: T {};
template<typename T, typename ... Ts> struct and_<T, Ts...>
: std::integral_constant<bool,
    (T::value ? and_<Ts...>::value : false)> {};
}//meta

template<typename T, std::size_t N>
struct coordinate
{
    typedef T value_type;

    template<typename ... Ts, std::enable_if<sizeof...(Ts) == N &&
        meta::and_<std::is_convertible<Ts, value_type>...>::value,
        std::nullptr_t>::type = nullptr>
    coordinate(Ts&& ... xs): values_{static_cast<value_type>(xs)...} {}
};

こうすることで、任意のN次元ベクトルについて、その要素型value_typeに変換可能なあらゆる型(それぞれ異なっていても良い)を引数としてちょうどN個取るようなコンストラクタが作れる。

基本的な発想は、まずvariadic templateによって任意個の引数を取れるようにした後、enable_ifによって適切な数の引数(sizeof...(Ts) == N)と適切な型(value_typeis_convertibleな型)の実体化以外を弾いているわけだ。ここで渡された任意個のstd::integral_constant<bool, /*~~*/>風の値を再帰的に評価するand_を使っている。

だが、これは「わかりやすい」だろうか? 私はそうは思わない。これは突っ込んだ技術を学ぶ気がある人なら少し時間をかければ誰でもいずれ思いつくようになる小技だが、逆にプログラミングの世界は表層だけ知れば十分で、奥深くに時間をかけて分け入っていくことを時間の無駄だと切り捨てる人には永遠に辿りつけないだろう。そして潜在的なユーザーには後者も少なからずいる。

そうでなくとも、特定の言語について長い時間をかけて学ぶことのバランスは難しい。ある程度手に馴染んで、今まで知っていた言語にも影響がでるくらいに、考え方まで感化されるようにならないと学ぶ意味はないに等しいが(使わないといけないライブラリがある場合はプログラミングの世界の外にそれを学ぶ意味がある)、かといってコードの書き方や考え方に影響を及ぼすほど手に馴染ませるには結構な時間をかけてその言語の突っ込んだことを学ばないといけない。他の言語でのスタイルを持ち込んでいてはだめで、その言語のスタイルと特徴的な問題解決法を身につけなければならないからだ。だが、その中にもやはり本質的に重要でない部分もあるわけで、そこを分けてやっていくのはやはり難しい。というか無理だ。

コードが単純になるが筋が悪くなるか、コードが技巧的になるが意味が一貫するかなら後者を取りたい。だが最近は「3次元しか使わないなら、拡張するときはプロジェクト全体を非推奨にすると決めて汎用性なんか捨てて3次元のものだけを単純に書けばいいのでは」と思うようにもなってきた。趣味のプログラムでは好きにやっていくけれども。C++を学び始めて2年半経って、ようやく1周回ってきた気がする。理由は主に現実を知ったからなのだが。ただ、現実を知ることと現実をそのままで肯定することは違う。現実を知ってある程度妥協できるようになっても、良くないコードとそれを生む土壌を憎む気持ちは保たなくてはならない。そう簡単なことではないが……。

私としては、技巧的なことが必要なくなれば皆が良質なコードを書けると思っているので(善良な人間が多数派であるという前提がある)、N4072が採択されてこのようなコードを書くのが非常に簡単になってくれるのが一番よいと思う。まあそれも新しいルールを追加することにはなるので、結局言語全体は複雑になっていくのだが。C++が大多数のプログラマの手に負えないほど複雑になった理由がまさにこの蓄積なら、この世に神はないことがわかる。そしてそれは悲しいことに事実だ。後は、知見が溜まってルールが増えすぎたときに、あまり良くなかったものをもっともっと切り捨てる勇気が、残ったものをうまくまとまるように再設計する発想が、それを広めていく地道な活動が、必要なのだろう。そしてそれに最も近いものがある。Rustだ。

話が逸れた。個人的な主張が爆発してしまうのはよくない。これは技術記事で、愚痴や宗教勧誘ではなかったはずだ。

問題としたかったのは、多次元のベクトルを使うようなライブラリを書くときにどうなるかという話だ。 そこで使う多くのコードがベクトルの内部表現に関係しないとする。というか、アーキテクチャ依存の特殊な命令を使うときなどを除いて、通常は処理がデータの表現方法に関係あることなどない。 つまり、上で挙げたような表現のうちどんなものが来ても、やることは変わらないはずで、汎用的なライブラリを書く場合はできるだけどんな場合でも対応できてほしい。

しかし値へのアクセス方法が完璧に異なっているようなクラスを同様に扱うのは難しい。このような場合、その違いを緩衝するような関数を一度噛ませることで解決できる。

template<typename T, std::size_t I>
struct access;

template<typename T, std::size_t N, std::size_t I>
struct access<coordinate<T, N>>
{
    static_assert(N > I, "out_of_range");
    static T&       invoke(coordinate<T, N>& p)       noexcept {return p[I]};
    static T const& invoke(coordinate<T, N> const& p) noexcept {return p[I]};
};

template<typename T>
struct access<xyz<T>, 0>
{
    static T&       invoke(xyz<T>& p)       noexcept {return p.x};
    static T const& invoke(xyz<T> const& p) noexcept {return p.x};
};
template<typename T>
struct access<xyz<T>, 1>
{
    static T&       invoke(xyz<T>& p)       noexcept {return p.y};
    static T const& invoke(xyz<T> const& p) noexcept {return p.y};
};
template<typename T>
struct access<xyz<T>, 2>
{
    static T&       invoke(xyz<T>& p)       noexcept {return p.z};
    static T const& invoke(xyz<T> const& p) noexcept {return p.z};
};

そして、内部で値にアクセスする必要が生じた際には、.xでも[0]でもなく、これを使うようにするのだ。 例えば、簡単な例として、座標間の距離を計算するdistance関数を書いてみる。

namespace detail
{

template<typename T, std::size_t N>
struct distance_sq_impl
{
    static typename coordinate_traits<T>::scalar_type
    invoke(const T& lhs, const T& rhs)
    {
        const auto dr = access<T, N>::invoke(lhs) -
                        access<T, N>::invoke(rhs);
        return dr * dr + distance_sq_impl<T, N-1>(lhs, rhs);
    }
};

template<typename T>
struct distance_sq_impl<T, 0>
{
    static typename coordinate_traits<T>::scalar_type
    invoke(const T& lhs, const T& rhs)
    {
        const auto dr = access<T, 0>::invoke(lhs) -
                        access<T, 0>::invoke(rhs);
        return dr * dr;
    }
};

}//detail

template<typename T>
inline typename coordinate_traits<T>::scalar_type
distance(const T& lhs, const T& rhs)
{
    return std::sqrt(detail::distance_sq_impl<
        T, coordinate_traits<T>::size-1>::invoke(lhs, rhs));
};

非常に面倒なことになったが、汎用的で質の良いライブラリの開発者は見えないところでこういったことを適切にやっている。 どの関数を使うかはコンパイル時に確定する上に、アクセス関数は全て非常に単純で、またメンバ関数をその場に書いているので暗黙のinline指定が入る。 まともなコンパイラならインライン化によってこれらの緩衝関数を実際のコードから消し去ってくれるだろう。 template再帰コンパイル時に決まるので、適切なインライン化が入れば、上のdistanceは以下のような実装と同じになる。

template<typename T>
inline typename coordinate_traits<T>::scalar_type
distance(const T& lhs, const T& rhs)
{
    const auto dx = lhs[0] - rhs[0];
    const auto dy = lhs[1] - rhs[1];
    const auto dz = lhs[2] - rhs[2];
    return std::sqrt(dx * dx + dy * dy + dz * dz);
};

だがこれではまだ問題が残っている。わかるだろうか。上記のaccessは、値を変更するときに参照を返している。 これは、受け取る座標ベクトル型がどんな形であれ全ての座標の値をどこかに格納していて、それを書き換えた時の影響がそこに留まるという前提を置いている。 今まで上げてきた座標ベクトルの例ではそれは正しい。しかし、そうでないような例を思いつくことは可能だ。例えば以下のような型はどうだろう。

template<typename T>
struct strange
{
// x, y座標はそのまま持っているが、z座標は持たず、ベクトルの長さを持っている
    T x, y, r;

    T x() const noexcept {return x;}
    T y() const noexcept {return y;}

// なので、z座標の値は都度計算して返す
    T z() const {return std::sqrt(r*r - x*x - y*y);}

// もちろんどの座標を変更するときもrに影響が出る
    void set_x(T x_) const
    {
        const auto z_ = this->z();
        x = x_;
        r = std::sqrt(x_ * x_ + y * y + z_ * z_);
        return;
    }
};

誰がこんな奇妙な表現をするか、という反論は理にかなっている。これは極端な例だ。 ただし、ここまで奇妙な表現でなくとも、ベクトルが自身の長さを単にキャッシュしている場合はあり得る。 非常に頻繁に長さの計算をするが、ベクトルの実際の値はあまり変更されない、またはされてもスカラー倍が多い場合、キャッシュしておくことに意味があるケースは存在しうる。と思う。

そのようなベクトルは、まともな設計なら、set_xのような関数を持っていて適切にキャッシュしておいた長さの値を更新するだろう。 しかしそのような場合、単にxにアクセスしてその値を書き換えてしまうと、事前に計算しておいた長さの値と齟齬が発生する。 だから、そういうベクトルの値を書き換えるときは、ライブラリの側も特別な関数を使わなければならない。なので適切なaccessは以下のようになる。

template<typename T, std::size_t I>
struct access
{
    typedef typename coordinate_traits<T>::scalar_type scalar_type;
    static void        set(T& p, scalar_type s);
    static scalar_type get(T const& p);
};

このような処理をベースにその後の処理を積み上げていくことで、あらゆるユーザーコードに(ほんの少しの面倒はあるが)引っかかることなくするりと入り込んでいくようなライブラリが書ける。 あとは、templateの特殊化などの現実的で面倒なことをユーザーから隠すため、何種類かのマクロを提供するだけだ。 マクロだと全ての状況には対応できないかもしれないが、もちろん、心得のあるユーザーはこのようなものが裏にあることを見抜くので、必要あらば手作業で特殊化するだろう。

何故こんな記事を書いているかというと、今まさにこのようなことをするかどうか迷っているからだ。 求める処理を提供していると思ったライブラリがどれもことごとく肝心の部分を外しており、とりあえずその中で一番まともだったものに少しトリックを入れて動かしているが、パフォーマンスやリーダビリティを考えると一から書いた方がよいのは間違いない。 あまり自明でない上に速度を落とすだろうトリックが、奥に隠蔽したとはいえ複数個入っているからだ。

だが自分で一から汎用的になるように書くとなると上のようなことをした方が良いのだろうか、と思い悩んでいるのだ。 実際には、ここまでするのは非常に激烈な場合だけだ。普通はベクトル型もセットで提供してしまうか、ユーザー定義ベクトル型が満たすべき要求を提示(し、可能ならコンパイル時にチェック)するだけで構わないだろう。 大抵のライブラリはこの理想的だが無限に面倒な実装をするどころか、固定サイズのベクトルをstd::vectorで管理していたりする。 流石にそれはどうかと思うが、自分の生産性がまず第一だ。完成しなければ話にならない。この辺りが本当に難しい。

ところで、Boost.QVMというライブラリはその目的にこのような面倒を代行することを含んでいる。 また、Boost.geometryは上記の「理想的だが無限に面倒な実装」を位置ベクトル以外にも非常に多数の対象について実際に行っている。 この記事はかなりの分Boost.geometryの実装に影響を受けており、ほぼその設計の基本的なアイデアの部分的な解説と言って過言ではない(中盤で爆発していた個人的な感情とBoost.Geometryは無関係だ)。

さて、しかしどうしようか。