hypotの近似

数学は強くないので、変なところがあるかも。

hypot(x, y)とは、sqrt(x*x+y*y)のことだけども、xが正なら、以下のように変形できる。

  sqrt(x*x+y*y)
= sqrt((y/x)^2+1)*x

q=sqrt(k^2+1)は、k (つまりy/x) の値が大きくなるほど+1の要素が小さくなるため、kの値が大きくなるほどsqrt(k^2)=kに近づく。そこで、sqrt(k^2)-kとする。すると、曲線が残るので、定数/kの曲線で近似するため、Octaveで差を0に近づけるようにソルバーをかける (kの範囲はてきとう)。

octave:1> k = (0.01:0.01:100000);
octave:2> q = @(n)max(sqrt(k.^2 + 1) - k - n./k)./k;
octave:3> fsolve(q, 0.1)
ans =  0.50000

すると、定数が求まるため、それを代入して、更に変形すると以下になる。

  (y/x + 0.5*x/y)*x
= y + 0.5 * x * x / y

と、ここまできて、既にstackoverflow.comに載ってることに気付いた。

For the record, here are a few more approximations, listed in roughly increasing order of complexity and accuracy. All these assume 0 ≤ a ≤ b.

h = b + 0.337 * a // max error ≈ 5.5 %
h = max(b, 0.918 * (b + (a>>1))) // max error ≈ 2.6 %
h = b + 0.428 * a * a / b // max error ≈ 1.04 %

係数の違いを見て、限界を感じる。近似アルゴリズム本が欲しいなぁ。
# にしても、近所の本屋はどこも品揃えが悪くて、中身見て買う派の自分は困る。図書館の方がまだ揃ってるとか何かがおかしい気がする…。