矮小

井の中の蛙

AOJ CGL_1_A 射影

大学の図書館に入っていたのでこの本をざっと読んでいました．AOJにある幾何の問題を全然解いていないことに気付いたので着手しました．簡単な問題から．

問題

平面上にある点からある直線に下ろした垂線の足を求める．

解釈

幾何的にも解けるけど線形代数を使って解いてみる．下図の \vec{p_1x} が分かれば \vec{x} = \vec{p_1} + \vec{p_1x} より求められる．

\vec{p_1x} は点 p1, p2 を通る直線上にあるため，

\vec{p_1x} = k\ \vec{p_1p_2} (k \in \bf{R})

と表せる．そのため

|\vec{p_1x}| = k\ |\vec{p_1p_2}| \cdots (1)

である（ k の正負は \cos \theta の値によって変わる）．また， \vec{p_1p} と \vec{p_1p_2} の内積を考えると

\vec{p_1p} \cdot \vec{p_1p_2} = |\vec{p_1p}||\vec{p_1p_2}| \cos \theta

だが，図中の直角三角形に着目すると

|\vec{p_1p}| \cos \theta = |\vec{p_1x}|

であることが分かる．したがって

\vec{p_1p} \cdot \vec{p_1p_2} = |\vec{p_1p_2}||\vec{p_1x}| より \displaystyle |\vec{p_1x}| = \frac{\vec{p_1p} \cdot \vec{p_1p_2}}{|\vec{p_1p_2}|} \cdots (2) である．式(1), (2)より，

\displaystyle k = \frac{\vec{p_1p} \cdot \vec{p_1p_2}}{|\vec{p_1p_2}|^2}

よって

\displaystyle \vec{p_1x} = \frac{\vec{p_1p} \cdot \vec{p_1p_2}}{|\vec{p_1p_2}|^2} \vec{p_1p_2} だから， \displaystyle \vec{x} = \vec{p_1} + \frac{\vec{p_1p} \cdot \vec{p_1p_2}}{|\vec{p_1p_2}|^2} \vec{p_1p_2}

\vec{p_1p_2} の係数がスカラーであることに気付かないと一見「あれ？」ってなる．とにかくここまで来れば最後の式を実装するだけ．

import std.stdio;
import std.string;
import std.algorithm;
import std.conv;
import std.range;
import std.numeric;

void main() {
    auto point = readln.chomp.split.map!(to!double).array;
    auto p1 = point[0..2];
    auto p2 = point[2..4];

    int q = readln.chomp.to!int;
    for (int i = 0; i < q; i++) {
        auto p = readln.chomp.split.map!(to!double).array;
        double[2] p1p = p[] - p1[];
        double[2] p1p2 = p2[] - p1[];
        double[2] ph = p1p2[] * dotProduct(p1p, p1p2) / dotProduct(p1p2, p1p2)
                       + p1[];
        writefln("%.10f %.10f", ph[0], ph[1]);
    }
}

\vec{p_1p},\ \vec{p_1p_2} を求めるには，各点同士の座標を引き算するだけ．D言語には配列のスライスを使ったベクトル演算が可能なため，p[] - p1[]のようなベクトルの定義に沿った書き方が出来る（p[] += 3などのようなベクトルに対するスカラー演算も可能です）．dotProductはstd.numericに定義されている．D言語は神．

参考