逆グーデルマン関数の計算式のJavaScriptにおける精度

メルカトル図法の投影座標の計算に逆グーデルマン関数を用いる(メルカトル図法 - Wikipedia)が、等価な計算式が複数ある。

JavaScriptで計算する際に、どの計算式を用いるのが誤差が少なくなるのか調査した。

比較条件

等価な数式は、Wikipedia(グーデルマン関数 - Wikipedia)にある以下のものとした。

  • (1): x => Math.log(Math.tan(Math.PI*0.25+x*0.5))
  • (2): x => Math.atanh(Math.sin(x))
  • (3): x => Math.asinh(Math.tan(x))
  • (4): x => Math.log(Math.tan(x)+1.0/Math.cos(x))
  • (5): x => Math.log((1.0+Math.sin(x))/Math.cos(x))
  • (6): x => 0.5*Math.log((1.0+Math.sin(x))/(1.0-Math.sin(x)))

このそれぞれの関数に対して、0 ≤ x ≤ π/2の区間を100,000分割したxを与えたときの計算結果を、Mathematicaで30桁精度で計算した結果と比較した。

また1,000,000回計算したときの計算時間も合わせて計測した(ループにかかる時間は引いたもの)。

環境はGoogle Chrome 93.0.4577.63とFirefox 92.0(いずれもUbuntu 20.04 LTS)で行った。

結果

各計算式の誤差の絶対値の比率、および計算時間をまとめたものは下表(太字は最も良い結果)。

関数式 誤差のRMS(Chrome) 誤差の最大値(Chrome) 計算時間(Chrome) (ms) 誤差のRMS(Firefox) 誤差の最大値(Firefox) 計算時間(Firefox) (ms)
(1) 2.87e-16 6.17e-12 54.2 2.87e-16 6.17e-12 65
(2) 4.18e-11 1.27e-8 59.4 4.18e-11 1.27e-8 53
(3) 8.28e-17 3.06e-16 71.9 8.25e-17 2.68e-16 80
(4) 3.33e-14 7.96e-12 68.5 3.33e-14 7.96e-12 79
(5) 2.37e-14 6.17e-12 52.8 2.37e-14 6.17e-12 69
(6) 4.18e-11 1.27e-8 36.4 4.18e-11 1.27e-8 43

また、xの値による誤差の出方は下図(横軸がx、縦軸が誤差)。

関数式 Chrome Firefox
(1)
(2)
(3)
(4)
(5)
(6)

これらから、(3)の計算式が最も誤差が少なく、またxの値によらず安定している。

Math.asinhはInternet Explorerには実装されていない(Math.asinh() - JavaScript | MDN)ので、 IEでも動作する(1)、(4)、(5)が次点として選択肢として入るが、

  • (1)はx = 0のとき-1.11e-16と0でない値が返ってくる
  • 計算速度が(5)が最も速い

の2点から、(5)が良いと考えられる。

Google ChromeとFirefoxの違いは(3)のみで、他の5つはほとんど同じ計算結果となった。

精度より計算速度を優先する場合は、(6)が1番速い結果となった。

結論

  • 計算式は(3)のx => Math.asinh(Math.tan(x))が最も精度がよく、Chrome、Firefox、Edgeなどのモダンブラウザではこれを使うのが良い。
  • Math.asinhが実装されていないIEをサポートする場合は(5)のx => Math.log((1.0+Math.sin(x))/Math.cos(x))
  • 計算速度を優先する場合は(6)のx => 0.5*Math.log((1.0+Math.sin(x))/(1.0-Math.sin(x)))