M.Hiroi's Home Page

お気楽 Ruby プログラミング入門

番外編 : 複素関数

[ PrevPage | R u b y | NextPage ]

複素関数

Ruby の場合、複素数に対応した数学関数 (複素関数) はモジュール CMath に用意されています。

CMath は標準ライブラリではないのでインストールする必要があります。公式ページには gem による方法が記載されていますが、Ubuntu 系の OS であれば以下のコマンドでインストールすることができます。

sudo apt install ruby-cmath

使用するときは cmath を require してください。

require "cmath"

●指数関数と対数関数

複素数を引数にとる指数関数はオイラー (Euler) の公式から導くことができます。

e = cos θ + i sin θ -- オイラーの公式
ex+iy = ex * eiy = ex * (cos y + i sin y)

複素数の対数関数は複素数 z を絶対値 |z| と偏角 θ を使って導くことができます。

x + iy = |z| * e
loge (x + iy) = loge (|z| * e) 
              = loge |z| + loge e
              = loge |z| + iθ, (-pi <= θ <= pi)

複素数 x, y のべき乗 xy は ey*log x で求めることができます。

モジュール CMath の関数 exp と log は複素数に対応しています。簡単な例を示しましょう。

irb> require "cmath"
=> true
irb> Math::PI
=> 3.141592653589793
irb> CMath.exp(0.0)
=> 1.0
irb> CMath.exp(Complex(0.0, 0.0))
=> (1.0+0.0i)
irb> CMath.exp(Complex(0.0, Math::PI / 2))
=> (6.123233995736766e-17+1.0i)
irb> CMath.exp(Complex(0.0, Math::PI))
=> (-1.0+1.2246467991473532e-16i)
irb> CMath.exp(Complex(0.0, -Math::PI))
=> (-1.0-1.2246467991473532e-16i)
irb> CMath.exp(Complex(1.0, 1.0))
=> (1.4686939399158851+2.2873552871788423i)

irb> CMath.log(1.0)
=> 0.0
irb> CMath.log(Complex(1.0, 0.0))
=> (0.0+0.0i)
irb> CMath.log(Complex(1.0, 1.0))
=> (0.3465735902799727+0.7853981633974483i)
irb> CMath.log(Complex(0.0, 1.0))
=> (0.0+1.5707963267948966i)
irb> CMath.log(Complex(1.0, -1.0))
=> (0.3465735902799727-0.7853981633974483i)
irb> CMath.log(Complex(1e300, 1e300))
=> (691.1221014884936+0.7853981633974483i)

irb> Complex(1.0, 1.0) ** 0
=> (1+0i)
irb> Complex(1.0, 1.0) ** 1
=> (1.0+1.0i)
irb> Complex(1.0, 1.0) ** 2
=> (0.0+2.0i)
irb> Complex(1.0, 1.0) ** 3
=> (-2.0+2.0i)
irb> Complex(1.0, 2.0) ** Complex(3.0, 4.0)
=> (0.129009594074467+0.03392409290517014i)

関数 log z (z = x + iy) は負の実軸 (-∞ < x < 0) において、x + 0.0j と x - 0.0j では値が異なります。

irb> CMath.log(Complex(-1.0, 0.0))
=> (0.0+3.141592653589793i)
irb> CMath.log(Complex(-1.0, -0.0))
=> (0.0-3.141592653589793i)
irb> CMath.log(Complex(-2.0, 0.0))
=> (0.6931471805599453+3.141592653589793i)
irb> CMath.log(Complex(-2.0, -0.0))
=> (0.6931471805599453-3.141592653589793i)
irb> CMath.log(Complex(-1e300, 0.0))
=> (690.7755278982137+3.141592653589793i)
irb> CMath.log(Complex(-1e300, -0.0))
=> (690.7755278982137-3.141592653589793i)

このように、関数 log z は負の実軸上で 2 つの値を持ちます。数学では値を一つ返す関数を「一価関数」、複数の値を返す関数を「多価関数」といいます。ここで、定義域を制限することで多価関数を一価関数にみなすことを考えます。関数 log z の場合、負の実軸を定義域から取り除けば、log z を一価関数とみなすことができるわけです。

参考 URL 7 によると、この取り除いた領域を branch cut と呼ぶそうです。プログラミングでは branch cut を定義域から取り除くのではなく、その領域では不連続な関数とするそうです。Python のドキュメントには、log は「分枝切断を一つもち、0 から負の実数軸に沿って -∞ へと延びており」と記述されています。Python のマニュアルに合わせて、本ページでも branch cut を「分枝切断」と記述することにします。

プログラミング言語の場合、0.0 と -0.0 を区別する処理系であれば、Python や Ruby のように 2 つの値を区別することができます。0.0 と -0.0 を区別しない処理系では、偏角 θ の範囲を -pi < θ <= pi に制限することで、log z の返り値を (-pi を取り除いて) 一つにすることができます。

●三角関数

複素数の三角関数の定義は、オイラーの公式から導かれる式の θ を複素数 z に変えたものになります。

e = cos θ + i sin θ  -- (1)
ei(-θ) = cos θ + i sin -θ  (sin -θ = -sin θ, cos -θ = cos θ)
      = cos θ - i sin θ -- (2)
(1) + (2) 
e + e-iθ = 2 cos θ
cos θ = (e + e-iθ) / 2
(1) - (2) 
e - e-iθ = 2i sin θ
sin θ = (e - e-iθ) / 2i

θ を複素数 z に置き換えた式が三角関数の定義になる
sin z = (eiz - e-iz) / 2i
cos z = (eiz + e-iz) / 2

sin z, cos z に純虚数 ix を与えると双曲線関数 (sinh x, cosh x) になります。

双曲線関数の定義
sinh x = (ex - e-x) / 2 
cosh x = (ex + e-x) / 2

sin ix = (eiix - e-iix) / 2i
       = (e-x - ex) / 2i  (分子と分母に -i を掛ける)
       = i (ex - e-x) / 2
       = i sinh x
cos ix = (eiix + e-iix) / 2
       = (e-x + ex) / 2
       = cosh x

これに三角関数の加法定理 [*1] を使うと次の式が導かれます。

sin(x + iy) = sin x cos iy + cos x sin iy
            = sin x cosh y + i cos x sinh y
cos(x + iy) = cos x cos iy - sin x sin iy
            = cos x cosh y - i sin x sinh y
tan(x + iy) = (sin 2x + sin 2iy) / (cos 2x + cos 2iy)
            = (sin 2x + i sinh 2y) / (cos 2x + cosh 2y) 

モジュール CMath にある三角関数 (sin, cos, tan) は複素数に対応しています。簡単な実行例を示します。

irb> a = [Complex(0.0, 1.0),Complex(0.0,-1.0),Complex(1.0,1.0),Complex(1.0,-1.0)]
=> [(0.0+1.0i), (0.0-1.0i), (1.0+1.0i), (1.0-1.0i)]

irb> for x in a
irb>   print CMath.sin(x), "\n", CMath.sin(x).abs, "\n"
irb> end
0.0+1.1752011936438014i
1.1752011936438014
0.0-1.1752011936438014i
1.1752011936438014
1.2984575814159773+0.6349639147847361i
1.4453965766582495
1.2984575814159773-0.6349639147847361i
1.4453965766582495
=> [(0.0+1.0i), (0.0-1.0i), (1.0+1.0i), (1.0-1.0i)]

irb> for x in a
irb>   print CMath.cos(x), "\n", CMath.cos(x).abs, "\n"
irb> end
1.5430806348152437-0.0i
1.5430806348152437
1.5430806348152437+0.0i
1.5430806348152437
0.8337300251311491-0.9888977057628651i
1.2934544550420957
0.8337300251311491+0.9888977057628651i
1.2934544550420957
=> [(0.0+1.0i), (0.0-1.0i), (1.0+1.0i), (1.0-1.0i)]

irb> for x in a
irb>   print CMath.tan(x), "\n", CMath.tan(x).abs, "\n"
irb> end
0.0+0.7615941559557649i
0.7615941559557649
0.0-0.7615941559557649i
0.7615941559557649
0.2717525853195118+1.0839233273386946i
1.1174700207060704
0.2717525853195118-1.0839233273386946i
1.1174700207060704
=> [(0.0+1.0i), (0.0-1.0i), (1.0+1.0i), (1.0-1.0i)]
-- note --------
[*1] 三角関数の公式は引数が複素数でも成り立ちます。ただし、|sin x| <= 1, |cos x| <= 1 という関係式は、x が実数だと成立しますが複素数では成立しません。

●双曲線関数

複素数の双曲線関数の定義は、実数の定義で引数 x を複素数 z に変えたものになります。

双曲線関数の定義 (z は複素数)
sinh z = (ez - e-z) / 2 
cosh z = (ez + e-z) / 2

sinh z, cosh z に純虚数 ix を与えると三角関数 (sin x, cos x) になります。

sinh ix = (eix - e-ix) / 2
        = i (eix - e-ix) / 2i  (分子と分母に i を掛ける)
        = i sin x
cosh ix = (eix + eix) / 2
        = cos x

これに双曲線関数の加法定理を使うと、次の式が導かれます。

双曲線関数の加法定理
sinh(x + y) = sinh x cosh y + cosh x sinh y
cosh(x + y) = cosh x cosh y + sinh x sinh y
tanh(x + y) = sinh(x + y) / cosh(x + y)
            = (sinh 2x + sinh 2y) / (cosh 2x + cosh 2y)

sinh(x + iy) = sinh x cos y + i cosh x sin y
cosh(x + iy) = cosh x cos y + i sinh x sin y
tanh(x + iy) = sinh(x + iy) / cosh(x + iy)
             = (sinh 2x + i sin 2y) / (cosh 2x + cos 2y)

モジュール CMath では双曲線関数 (sinh, cosh, tanh) を使用することができます。簡単な使用例を示します。

irb> a
=> [(0.0+1.0i), (0.0-1.0i), (1.0+1.0i), (1.0-1.0i)]

irb> for x in a
irb>   print CMath.sinh(x), "\n"
irb> end
0.0+0.8414709848078965i
0.0-0.8414709848078965i
0.6349639147847361+1.2984575814159773i
0.6349639147847361-1.2984575814159773i
=> [(0.0+1.0i), (0.0-1.0i), (1.0+1.0i), (1.0-1.0i)]

irb> for x in a
irb>   print CMath.cosh(x), "\n"
irb> end
0.5403023058681398+0.0i
0.5403023058681398-0.0i
0.8337300251311491+0.9888977057628651i
0.8337300251311491-0.9888977057628651i
=> [(0.0+1.0i), (0.0-1.0i), (1.0+1.0i), (1.0-1.0i)]

irb> for x in a
irb>   print CMath.tanh(x), "\n"
irb> end
0.0+1.557407724654902i
0.0-1.557407724654902i
1.0839233273386946+0.2717525853195118i
1.0839233273386946-0.2717525853195118i
=> [(0.0+1.0i), (0.0-1.0i), (1.0+1.0i), (1.0-1.0i)]

●複素数の平方根

複素数 z の平方根は次の式で求めることができます。

z = x + iy, |z| = √(x2 + y2), θ = z の偏角 (-pi <= θ <= pi)
とすると
√(x + iy) =
√(|z| * e) = √|z| * eiθ/2  -- (1)
√(|z| * eiθ+ 2pi) = √|z| * eiθ/2 + pi  -- (2)

式 (1) を平方根の主値といいます。角度は 2pi を足すと同じ角度になるので、式 (2) がもう一つの解になります。三角関数の半角の公式を使うと、式 (1) から次の式が導かれます。

三角関数の半角の公式
sin2(θ/2) = (1 - cos θ) / 2
cos2(θ/2) = (1 + cos θ) / 2

y >= 0 の場合
√|z| * eiθ/2 = √|z| * (cos(θ/2) + i sin(θ/2))
            = √|z| * (√((1 + cos θ)/2) + i √((1 - cos θ)/2))
            = √((|z| + |z|cos θ)/2 + i √((|z| - |z|cos θ)/2) 
|z|cos θ = x だから
            = √((|z| + x) / 2) + i √((|z| - x) / 2)
y < 0 の場合、虚部の符号が - になる
            = √((|z| + x) / 2) - i √((|z| - x) / 2)

モジュール CMath にある sqrt は複素数に対応しています。簡単な実行例を示します。

irb> CMath.sqrt(1.0)
=> 1.0
irb> CMath.sqrt(2.0)
=> 1.4142135623730951
irb> CMath.sqrt(-1.0)
=> (0+1.0i)
irb> CMath.sqrt(-2.0)
=> (0+1.4142135623730951i)

irb> b = CMath.sqrt(Complex(1.0, 1.0))
=> (1.09868411346781+0.4550898605622274i)
irb> b * b
=> (1.0+1.0000000000000002i)
irb> b = CMath.sqrt(Complex(1.0, -1.0))
=> (1.09868411346781-0.4550898605622274i)
irb> b * b
=> (1.0-1.0000000000000002i)

sqrt は log と同じ分枝切断を持っています。x を負の整数とすると sqrt(x) の解は i √x になりますが、もうひとつ -i √x という解があります。

irb> CMath.sqrt(Complex(-1.0, 0.0))
=> (0.0+1.0i)
irb> CMath.sqrt(Complex(-1.0, -0.0))
=> (0.0-1.0i)
irb> CMath.sqrt(Complex(-2.0, 0.0))
=> (0.0+1.4142135623730951i)
irb> CMath.sqrt(Complex(-2.0, -0.0))
=> (0.0-1.4142135623730951i)
irb> CMath.sqrt(Complex(-1e300, 0.0))
=> (0.0+1.0e+150i)
irb> CMath.sqrt(Complex(-1e300, -0.0))
=> (0.0-1.0e+150i)

●逆三角関数

三角関数の逆関数を「逆三角関数 (inverse trigonometric function)」といいます。以下にモジュール Math にある逆三角関数を示します。

ここでは引数 x を実数とします。asin x は引数 x が与えられたとき sin w = x となる角度 w を求めます。同様に acos x は cos w = x となる角度 w を、atan x は tan x = w となる角度 w を求めます。三角関数には周期性があるので、上式を満たす角度 w は無数に存在します。つまり、逆三角関数の返り値は無数にあることになりますが、通常は一つの値を返すように範囲を制限します。これを「主値」といいます。

逆三角関数の主値を以下に示します。

簡単な実行例を示します。

irb> for x in [-1.0, -0.5, 0.0, 0.5, 1.0]
irb>   print Math.asin(x), "\n"
irb> end
-1.5707963267948966
-0.5235987755982989
0.0
0.5235987755982989
1.5707963267948966
=> [-1.0, -0.5, 0.0, 0.5, 1.0]

irb> for x in [-1.0, -0.5, 0.0, 0.5, 1.0]
irb>   print Math.acos(x), "\n"
irb> end
3.141592653589793
2.0943951023931957
1.5707963267948966
1.0471975511965979
0.0
=> [-1.0, -0.5, 0.0, 0.5, 1.0]

irb> for x in [-1e300, -1.0, 0.0, 1.0, 1e300]
irb>   print Math.atan(x), "\n"
irb> end
-1.5707963267948966
-0.7853981633974483
0.0
0.7853981633974483
1.5707963267948966
=> [-1.0e+300, -1.0, 0.0, 1.0, 1.0e+300]

モジュール Math には 2 引数の関数 atan2 も用意されています。これは他のプログラミング言語、たとえばC言語の数学関数 atan2 と同じです。

atan2(y, x) => 角度 (ラジアン)

引数 x, y は実数です。atan2 は直交座標系においてベクトル (x, y) と x 軸との角度を求める関数です。複素平面で考えると、複素数 x + iy の偏角 θ を求めることと同じです。返り値 (角度 θ) の範囲は -pi <= θ <= pi になります。

簡単な実行例を示します。

irb> Math.atan2(0.0, 1.0)
=> 0.0
irb> Math.atan2(1.0, 1.0)
=> 0.7853981633974483
irb> Math.atan2(1.0, 0.0)
=> 1.5707963267948966
irb> Math.atan2(1.0, -1.0)
=> 2.356194490192345
irb> Math.atan2(0.0, -1.0)
=> 3.141592653589793
irb> Math.atan2(-1.0, 1.0)
=> -0.7853981633974483
irb> Math.atan2(-1.0, 0.0)
=> -1.5707963267948966
irb> Math.atan2(-1.0, -1.0)
=> -2.356194490192345
irb> Math.atan2(-0.0, -1.0)
=> -3.141592653589793

●複素数の逆三角関数

複素数の逆三角関数の定義は、複素数の三角関数の定義から導くことができます。asin z の定義は次のようになります。

asin z = w (z, w は複素数)
z = sin w = (eiw - e-iw) / 2i
両辺に 2ieiw を掛けると
2iz(eiw) = (eiw)2 - 1
(eiw)2 - 2iz(eiw) - 1 = 0
eiw の二次方程式と考えて解くと
eiw = iz ± √(1 - z2)
両辺の対数をとって -i を掛け算し、平方根の主値 (+) を選ぶ
w = asin z = -i log(iz + √(1 - z2))

acos z と atan z は定義だけを示します。

acos z = pi/2 - asin z
atan z = i(log(1 - iz) - log(1 + iz)) / 2

asin, acos, atan は次に示す分枝切断を持っています。

モジュール CMath の逆三角関数 (asin, acos, atan) は複素数に対応しています。簡単な実行例を示します。

irb> for x in [-1.0, -0.75, -0.5, -0.25, 0.0, 0.25, 0.5, 0.75, 1.0]
irb>   print CMath.asin(x), "\n"
irb> end
-1.5707963267948966
-0.848062078981481
-0.5235987755982989
-0.25268025514207865
0.0
0.25268025514207865
0.5235987755982989
0.848062078981481
1.5707963267948966
=> [-1.0, -0.75, -0.5, -0.25, 0.0, 0.25, 0.5, 0.75, 1.0]
irb> CMath.asin(Complex(0.0, 1.0))
=> (0.0+0.8813735870195428i)
irb> CMath.sin(CMath.asin(Complex(0.0, 1.0)))
=> (0.0+0.9999999999999997i)
irb> CMath.asin(Complex(1.0, 1.0))
=> (0.6662394324925153+1.0612750619050355i)
irb> CMath.sin(CMath.asin(Complex(1.0, 1.0)))
=> (1.0+0.9999999999999998i)

irb> CMath.asin(Complex(2.0, 0.0))
=> (1.5707963267948966-1.3169578969248166i)      # 虚部の符号が逆
irb> CMath.asin(Complex(2.0, -0.0))
=> (1.5707963267948966-1.3169578969248166i)
irb> CMath.asin(Complex(-4.0, 0.0))
=> (-1.5707963267948966+2.0634370688955617i)
irb> CMath.asin(Complex(-4.0, -0.0))
=> (-1.5707963267948966+2.0634370688955617i)     # 虚部の符号が逆

irb> for x in [-1.0, -0.75, -0.5, -0.25, 0.0, 0.25, 0.5, 0.75, 1.0]
irb>   print CMath.acos(x), "\n"
irb> end
3.141592653589793
2.4188584057763776
2.0943951023931957
1.8234765819369754
1.5707963267948966
1.318116071652818
1.0471975511965979
0.7227342478134157
0.0
=> [-1.0, -0.75, -0.5, -0.25, 0.0, 0.25, 0.5, 0.75, 1.0]
irb> CMath.acos(Complex(0.0, 1.0))
=> (1.5707963267948966-0.8813735870195429i)
irb> CMath.cos(CMath.acos(Complex(0.0, 1.0)))
=> (8.659560562354932e-17+1.0i)
irb> CMath.acos(Complex(1.0, 1.0))
=> (0.9045568943023813-1.0612750619050357i)
irb> CMath.cos(CMath.acos(Complex(1.0, 1.0)))
=> (1.0000000000000002+1.0i)

irb> CMath.acos(Complex(2.0, 0.0))
=> (0.0+1.3169578969248164i)                     # 虚部の符号が逆
irb> CMath.acos(Complex(2.0, -0.0))
=> (0.0+1.3169578969248164i)
irb> CMath.acos(Complex(-4.0, 0.0))
=> (3.141592653589793-2.0634370688955608i)
irb> CMath.acos(Complex(-4.0, -0.0))
=> (3.141592653589793-2.0634370688955608i)       # 虚部の符号が逆

irb> for x in [2.0, 1.5, 1.0, 0.5, 0.0, -0.5, -1.0, -1.5, -2.0]
irb>   print CMath.atan(x), "\n"
irb> end
1.1071487177940904
0.982793723247329
0.7853981633974483
0.4636476090008061
0.0
-0.4636476090008061
-0.7853981633974483
-0.982793723247329
-1.1071487177940904
=> [2.0, 1.5, 1.0, 0.5, 0.0, -0.5, -1.0, -1.5, -2.0]

irb> CMath.atan(Complex(1.0, 1.0))
=> (1.0172219678978514+0.4023594781085251i)
irb> CMath.tan(CMath.atan(Complex(1.0, 1.0)))
=> (1.0+0.9999999999999999i)
irb> CMath.atan(Complex(1.0, -1.0))
=> (1.0172219678978514-0.40235947810852507i)
irb> CMath.tan(CMath.atan(Complex(1.0, -1.0)))
=> (1.0000000000000002-0.9999999999999999i)

irb> CMath.atan(Complex(0.0, 2.0))
=> (-1.5707963267948966+0.5493061443340549i)      # 実部の符号が逆
irb> CMath.atan(Complex(-0.0, 2.0))
=> (1.5707963267948966+0.5493061443340549i)       # 実部の符号が逆
irb> CMath.atan(Complex(0.0, -4.0))
=> (1.5707963267948966-0.25541281188299536i)
irb> CMath.atan(Complex(-0.0, -4.0))
=> (-1.5707963267948966-0.25541281188299536i)

M.Hiroi が使っている CMath (ver 1.0.0) の逆三角関数には分枝切断に不具合があるようです。

●逆双曲線関数

双曲線関数の逆関数を「逆双曲線関数 (inverse hyperbolic function)」といいます。モジュール Math で用意されている逆双曲線関数には asinh, acosh, atanh があります。双曲線関数と逆双曲線関数の定義域と値域を示します。

x と y は実数です。x = cosh y の逆関数 y = acosh x を満たす y の値は 2 つありますが、ここでは y >= 0 を主値として選ぶことにします。

逆双曲線関数の定義は双曲線関数の定義から導くことができます。

asinh x = y
y = sinh x
  = (ex - e-x) / 2
両辺に 2ex を掛けると
2yex = e2x - 1
e2x - 2yex - 1 = 0, (ex の 2 次方程式として解く)
ex = (2y ± √(4y2 + 4))/2 = y ± √(y2 + 1)
ex > 0 だから平方根の符号に + を選んで対数をとる
x = log (y + √(y2 + 1))
acosh x = y
y = cosh x
  = (ex + e-x) / 2
両辺に 2ex を掛けると
2yex = e2x + 1
e2x - 2yex + 1 = 0, (ex の 2 次方程式として解く)
ex = (2y ± √(4y2 - 4))/2 = y ± √(y2 - 1)
x = log (y ± √(y2 - 1))
平方根の符号 + を主値としてとると
x = log (y + √(y2 - 1)), (y >= 1)
atanh x = y
y = tanh x 
  = sinh x / cosh x
  = (ex - e-x) / (ex + e-x)
  = (e2x - 1) / (e2x + 1)
y(e2x + 1) = (e2x - 1)
(1 - y)e2x = 1 + y
e2x = (1 + y) / (1 - y)
x = (1 / 2) * log((1 + y) / (1 - y))

また、次の関係式を使って atanh から asinh と acosh を求めることができます。

asinh x = atanh(x / √(x2 + 1))
acosh x = atanh(√(1 - x2) / x)

簡単な実行例を示します。

irb> for x in [2.0, 1.5, 1.0, 0.5, 0.0, -0.5, -1.0, -1.5, -2.0]
irb>   print Math.asinh(x), "\n"
irb> end
1.4436354751788103
1.1947632172871094
0.881373587019543
0.48121182505960347
0.0
-0.48121182505960347
-0.881373587019543
-1.1947632172871094
-1.4436354751788103
=> [2.0, 1.5, 1.0, 0.5, 0.0, -0.5, -1.0, -1.5, -2.0]
irb> for x in [2.0, 1.5, 1.0, 0.5, 0.0, -0.5, -1.0, -1.5, -2.0]
irb>   print Math.sinh(Math.asinh(x)), "\n"
irb> end
1.9999999999999998
1.5000000000000004
1.0
0.5
0.0
-0.5
-1.0
-1.5000000000000004
-1.9999999999999998
=> [2.0, 1.5, 1.0, 0.5, 0.0, -0.5, -1.0, -1.5, -2.0]

irb> for x in [1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0]
irb>   print Math.acosh(x), "\n"
irb> end
0.0
0.9624236501192069
1.3169578969248166
1.566799236972411
1.762747174039086
1.9248473002384139
2.0634370688955608
=> [1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0]
irb> for x in [1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0]
irb>   print Math.cosh(Math.acosh(x)), "\n"
irb> end
1.0
1.5
1.9999999999999998
2.5
3.0000000000000004
3.5000000000000004
4.000000000000001
=> [1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0]

irb> for x in [-0.75, -0.5, -0.25, 0.0, 0.25, 0.5, 0.75]
irb>   print Math.atanh(x), "\n"
irb> end
-0.9729550745276566
-0.5493061443340548
-0.25541281188299536
0.0
0.25541281188299536
0.5493061443340548
0.9729550745276566
=> [-0.75, -0.5, -0.25, 0.0, 0.25, 0.5, 0.75]
irb> for x in [-0.75, -0.5, -0.25, 0.0, 0.25, 0.5, 0.75]
irb>   print Math.tanh(Math.atanh(x)), "\n"
irb> end
-0.75
-0.49999999999999994
-0.25
0.0
0.25
0.49999999999999994
0.75
=> [-0.75, -0.5, -0.25, 0.0, 0.25, 0.5, 0.75]

●複素数の逆双曲線関数

複素数の逆双曲線関数の定義は、実数の定義で引数 x を複素数 z に変えたものになります。

asinh z = log (z + √(z2 + 1))
acosh z = log (z + √(z2 - 1)) -- (1)
        = log (z + (√(z + 1)) * (√(z - 1))) -- (2)
atanh z = (1 / 2) * log((1 + z) / (1 - z))
        = (1 / 2) * (log(1 + z) - log(1 - z))

asinh, acosh, atanh は次に示す分枝切断を持っています。

モジュール CMath の逆双曲線関数 (asinh, acosh, atanh) は複素数に対応していますが、CMaht.acosh は式 1 で定義されているため、虚軸上の値が他の処理系と異なることがあります。参考文献, URL 8 によると、acosh(z) を式 (1) でプログラムすると「分枝切断線」が複雑になるため、他の式 (たとえば式 2 など) でプログラムする処理系が多いようです。

簡単な実行例を示します。

irb> for x in [2.0, 1.5, 1.0, 0.5, 0.0, -0.5, -1.0, -1.5, -2.0]
irb>   print CMath.asinh(Complex(x, 0.0)), "\n"
irb> end
1.4436354751788103+0.0i
1.1947632172871094+0.0i
0.8813735870195429+0.0i
0.48121182505960347+0.0i
0.0+0.0i
-0.48121182505960336+0.0i
-0.8813735870195428+0.0i
-1.1947632172871097+0.0i
-1.4436354751788099+0.0i
=> [2.0, 1.5, 1.0, 0.5, 0.0, -0.5, -1.0, -1.5, -2.0]
irb> for x in [2.0, 1.5, 1.0, 0.5, 0.0, -0.5, -1.0, -1.5, -2.0]
irb>   print CMath.sinh(CMath.asinh(Complex(x, 0.0))), "\n"
irb> end
1.9999999999999998+0.0i
1.5000000000000004+0.0i
1.0+0.0i
0.5+0.0i
0.0+0.0i
-0.4999999999999999+0.0i
-0.9999999999999997+0.0i
-1.5000000000000004+0.0i
-1.999999999999999+0.0i
=> [2.0, 1.5, 1.0, 0.5, 0.0, -0.5, -1.0, -1.5, -2.0]
irb> CMath.asinh(Complex(1.0, 1.0))
=> (1.0612750619050357+0.6662394324925153i)
irb> CMath.sinh(CMath.asinh(Complex(1.0, 1.0)))
=> (1.0+1.0000000000000002i)

irb> CMath.asinh(Complex(0.0, 2.0))
=> (1.3169578969248166+1.5707963267948966i)
irb> CMath.asinh(Complex(-0.0, 2.0))
=> (-1.3169578969248164+1.5707963267948966i)
irb> CMath.asinh(Complex(0.0, -4.0))
=> (2.0634370688955608-1.5707963267948966i)
irb> CMath.asinh(Complex(-0.0, -4.0))
=> (-2.0634370688955617-1.5707963267948966i)

irb> for x in [1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0]
irb>   print CMath.acosh(Complex(x, 0.0)), "\n"
irb> end
0.0+0.0i
0.9624236501192069+0.0i
1.3169578969248166+0.0i
1.5667992369724109+0.0i
1.762747174039086+0.0i
1.9248473002384139+0.0i
2.0634370688955608+0.0i
=> [1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0]
irb> for x in [1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0]
irb>   print CMath.cosh(CMath.acosh(Complex(x, 0.0))), "\n"
irb> end
1.0+0.0i
1.5+0.0i
1.9999999999999998+0.0i
2.499999999999999+0.0i
3.0000000000000004+0.0i
3.5000000000000004+0.0i
4.000000000000001+0.0i
=> [1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0]
irb> CMath.acosh(Complex(1.0, 1.0))
=> (1.0612750619050357+0.9045568943023813i)
irb> CMath.cosh(CMath.acosh(Complex(1.0, 1.0)))
=> (1.0000000000000002+1.0i)

irb> CMath.acosh(Complex(0.0, 0.0))
=> (0.0+1.5707963267948966i)
irb> CMath.acosh(Complex(0.0, -0.0))
=> (0.0-1.5707963267948966i)
irb> CMath.acosh(Complex(-4.0, 0.0))
=> (-2.0634370688955617+3.141592653589793i)       # 実部の符号が逆
irb> CMath.acosh(Complex(-4.0, -0.0))
=> (-2.0634370688955617+3.141592653589793i)       # 実部と虚部の符号が逆

irb> for x in [-0.75, -0.5, -0.25, 0.0, 0.25, 0.5, 0.75]
irb>   print CMath.atanh(Complex(x, 0.0)), "\n"
irb> end
-0.9729550745276567+0.0i
-0.5493061443340549+0.0i
-0.25541281188299536+0.0i
0.0+0.0i
0.25541281188299536+0.0i
0.5493061443340549+0.0i
0.9729550745276566+0.0i
=> [-0.75, -0.5, -0.25, 0.0, 0.25, 0.5, 0.75]
irb> for x in [-0.75, -0.5, -0.25, 0.0, 0.25, 0.5, 0.75]
irb>   print CMath.tanh(CMath.atanh(Complex(x, 0.0))), "\n"
irb> end
-0.7500000000000001+0.0i
-0.5000000000000001+0.0i
-0.25+0.0i
0.0+0.0i
0.25+0.0i
0.5000000000000001+0.0i
0.7500000000000001+0.0i
=> [-0.75, -0.5, -0.25, 0.0, 0.25, 0.5, 0.75]
irb> CMath.atanh(Complex(1.0, 1.0))
=> (0.4023594781085251+1.0172219678978514i)
irb> CMath.tanh(CMath.atanh(Complex(1.0, 1.0)))
=> (0.9999999999999999+1.0i)

irb> CMath.atanh(Complex(2.0, 0.0))
=> (0.5493061443340549-1.5707963267948966i)      # 虚部の符号が逆
irb> CMath.atanh(Complex(2.0, -0.0))
=> (0.5493061443340549-1.5707963267948966i)
irb> CMath.atanh(Complex(-4.0, 0.0))
=> (-0.25541281188299536+1.5707963267948966i)
irb> CMath.atanh(Complex(-4.0, -0.0))
=> (-0.25541281188299536+1.5707963267948966i)    # 虚部の符号が逆

M.Hiroi が使っている CMath (ver 1.0.0) の逆双曲線関数には分枝切断に不具合があるようです。

●補足: 分枝切断の不具合

Ruby のモジュール CMath に限らず、複素関数の分枝切断に何らかの不具合がある処理系 (ライブラリ) はけっこう多いです。M.Hiroi が使ったことがある処理系の中で、分枝切断に不具合がなかったのは Python と Julia だけです。他の処理系、たとえばC言語や C++ にも不具合がありました。主な理由は「負のゼロ」の処理 (判定や演算) です。

一般に、プログラミング言語では実数を浮動小数点数 (floating point number) として表します。浮動小数点数には IEEE 754 という標準仕様があり、近代的なプログラミング言語のほとんどは、IEEE 754 に準拠した浮動小数点数をサポートしています。IEEE 754 には通常の数値以外にも、負のゼロ (-0.0)、正負の無限大 (∞, -∞)、NaN (Not a Number, 非数) といった値が定義されています。

負のゼロ (-0.0) は、計算結果が負の極めて小さな値でアンダーフローになったとき発生します。また、正の値を負の無限大で除算する、負の値を正の無限大で除算する、負の値と 0.0 を乗算しても -0.0 が得られます。標準 (IEEE 754) では、演算子 (== など) による 0.0 と -0.0 の比較は等しいと判定されます。Ruby には演算子のほかにも正、負、零を判定するメソッド posiive?, negative?, zero? がありますが、-0.0 はゼロと判定されます。

irb> 0.0 == -0.0
=> true
irb> -0.0 != 0.0
=> false
irb> 0.0 > -0.0
=> false
irb> 0.0 < -0.0
=> false
irb> -0.0.positive?
=> false
irb> -0.0.negative?
=> false
irb> -0.0.zero?
=> true

0.0 と -0.0 を区別するため、負のゼロを判定する関数が用意されています。Ruby ではメソッド phase を使います。

irb> 0.0.phase
=> 0
irb> -0.0.phase
=> 3.141592653589793
irb> 1.0.phase
=> 0
irb> -1.0.phase
=> 3.141592653589793

実数にメソッド phase を適用すると、正であれば 0 を、負であれば PI を返します。Ruby の場合、phase で負のゼロを判定することができます。

このほかに、負のゼロの演算には特別規則があります。加減算の場合、-0.0 + -0.0 は -0.0 になりますが、-0.0 + 0.0 や 0.0 + -0.0 は 0.0 になります。-0.0 - 0.0 は -0.0 になりますが、-0.0 - -0.0 や 0.0 - -0.0 は 0.0 になります。このような規則があるので、複素数を扱うときには注意してください。

ご参考までに、分枝切断の不具合を修正した複素関数 (簡易版) を プログラムリスト に示します。興味のある方は試してみてください。

●参考文献, URL

  1. 奥村晴彦,『C言語による最新アルゴリズム事典』, 技術評論社, 1991
  2. IEEE 754 -- Wikipedia
  3. IEEE 754における負のゼロ - Wikipedia
  4. NaN - Wikipedia
  5. 逆三角関数 - Wikipedia
  6. 逆双曲線関数 - Wikipedia
  7. 逆双曲線関数と逆三角関数の branch cut | 雑記帳

●プログラムリスト

# coding: utf-8
#
# cmath.rb : 複素関数 (簡易版)
#
#            Copyright (C) 2023 Makoto Hiroi
#

module Cmath
  # -0.0 を含んだ負の判定
  def negative?(x) = x.phase > 0.0

  # 実数と複素数の加減算
  # Ruby の x + z は add_fc() と同じ
  # Python は (x + z.real) + (0.0 + z.imag)i と同じになる
  def add_fc(x, z) = Complex(x + z.real, z.imag)
  def sub_fc(x, z) = Complex(x - z.real, -z.imag)

  # 指数関数 exp(z)
  def exp(z)
    a = Math.exp(z.real)
    Complex(a * Math.cos(z.imag), a * Math.sin(z.imag))
  end

  # 対数関数 log(z)
  def log(z)
    Complex(Math.log(z.abs), z.phase)
  end

  # 三角関数
  def sin(z)
    x = z.real
    y = z.imag
    Complex(Math.sin(x) * Math.cosh(y), Math.cos(x) * Math.sinh(y))
  end

  def cos(z)
    x = z.real
    y = z.imag
    Complex(Math.cos(x) * Math.cosh(y), - Math.sin(x) * Math.sinh(y))
  end

  def tan(z)
    x = z.real
    y = z.imag
    d = Math.cos(2.0 * x) + Math.cosh(2.0 * y)
    Complex(Math.sin(2.0 * x) / d, Math.sinh(2.0 * y) / d)
  end

  # 双曲線関数
  def sinh(z)
    x = z.real
    y = z.imag
    Complex(Math.sinh(x) * Math.cos(y), Math.cosh(x) * Math.sin(y))
  end

  def cosh(z)
    x = z.real
    y = z.imag
    Complex(Math.cosh(x) * Math.cos(y), Math.sinh(x) * Math.sin(y))
  end

  def tanh(z)
    x = z.real
    y = z.imag
    d = Math.cosh(2.0 * x) + Math.cos(2.0 * y)
    Complex(Math.sinh(2.0 * x) / d, Math.sin(2.0 * y) / d)
  end

  # 平方根
  def sqrt(z)
    x = z.real
    a = z.abs
    b = Math.sqrt((a - x) / 2.0)
    Complex(Math.sqrt((a + x) / 2.0), negative?(z.imag) ? -b : b)
  end

  # 複素数の逆三角関数
  def asin(z)
    z1 = sqrt(add_fc(1.0, z))
    z2 = sqrt(sub_fc(1.0, z))
    i  = Complex::I
    -1.0 * i * log(i * z + z1 * z2)
  end

  def acos(z) = Math::PI / 2.0 - asin(z)

  def atan(z)
    i = Complex::I
    y = i * z
    x = i * (log(sub_fc(1.0, y)) - log(add_fc(1.0, y))) / 2.0
    z.real == 0.0 && negative?(z.real) && z.imag > 1.0 ? Complex(-x.real, x.imag) : x
  end

  # 複素数の逆双曲線関数
  def asinh(z) = log(z + sqrt(z * z + 1.0))
  def acosh(z) = log(z + sqrt(z + 1.0) * sqrt(z - 1.0))
  def atanh(z) = (log(add_fc(1.0, z)) - log(sub_fc(1.0, z))) / 2.0

  module_function :negative?
  module_function :add_fc
  module_function :sub_fc
  module_function :exp
  module_function :log
  module_function :sin
  module_function :cos
  module_function :tan
  module_function :sinh
  module_function :cosh
  module_function :tanh
  module_function :sqrt
  module_function :asin
  module_function :acos
  module_function :atan
  module_function :asinh
  module_function :acosh
  module_function :atanh
end
リスト : 簡単なテスト (test_cmath.rb)

require "./cmath"

a = [Complex(1.0, 0.0), Complex(0.0, 1.0), Complex(1.0, 1.0), Complex(1.0, -1.0)]

for x in a
  printf "exp(%s) = %s\n", x, Cmath.exp(x)
end
for x in a
  printf "log(%s) = %s\n", x, Cmath.log(x)
end
for x in [Complex(-1.0, 0.0), Complex(-1e300, 0.0)]
  printf "log(%s) = %s\n", x, Cmath.log(x)
  printf "log(%s) = %s\n", x.conj, Cmath.log(x.conj)
end
for x in a
  printf "sin(%s) = %s\n", x, Cmath.sin(x)
end
for x in a
  printf "cos(%s) = %s\n", x, Cmath.cos(x)
end
for x in a
  printf "tan(%s) = %s\n", x, Cmath.tan(x)
end
for x in a
  printf "sinh(%s) = %s\n", x, Cmath.sinh(x)
end
for x in a
  printf "cosh(%s) = %s\n", x, Cmath.cosh(x)
end
for x in a
  printf "tanh(%s) = %s\n", x, Cmath.tanh(x)
end
for x in a
  printf "sqrt(%s) = %s\n", x, Cmath.sqrt(x)
end
for x in [Complex(-1.0, 0.0), Complex(-1e300, 0.0)]
  printf "sqrt(%s) = %s\n", x, Cmath.sqrt(x)
  printf "sqrt(%s) = %s\n", x.conj, Cmath.sqrt(x.conj)
end
for x in a
  printf "asin(%s) = %s\n", x, Cmath.asin(x)
end
for x in a
  printf "acos(%s) = %s\n", x, Cmath.acos(x)
end
for x in a
  printf "atan(%s) = %s\n", x, Cmath.atan(x)
end
for x in [Complex(2.0, 0.0), Complex(-4.0, 0.0)]
  printf "asin(%s) = %s\n", x, Cmath.asin(x)
  printf "asin(%s) = %s\n", x.conj, Cmath.asin(x.conj)
  printf "acos(%s) = %s\n", x, Cmath.acos(x)
  printf "acos(%s) = %s\n", x.conj, Cmath.acos(x.conj)
end
for x in [Complex(0.0, 2.0), Complex(-0.0, 2.0), Complex(0.0, -4.0), Complex(-0.0, -4.0)]
  printf "atan(%s) = %s\n", x, Cmath.atan(x)
end
for x in a
  printf "asinh(%s) = %s\n", x, Cmath.asinh(x)
end
for x in a
  printf "acosh(%s) = %s\n", x, Cmath.acosh(x)
end
for x in a
  printf "atanh(%s) = %s\n", x, Cmath.atanh(x)
end
for x in [Complex(0.0, 2.0), Complex(-0.0, 2.0), Complex(0.0, -4.0), Complex(-0.0, -4.0)]
  printf "asinh(%s) = %s\n", x, Cmath.asinh(x)
end
for x in [Complex(0.0, 0.0), Complex(0.0, -0.0), Complex(-4.0, 0.0), Complex(-4.0, -0.0)]
  printf "acosh(%s) = %s\n", x, Cmath.acosh(x)
end
for x in [Complex(2.0, 0.0), Complex(2.0, -0.0), Complex(-4.0, 0.0), Complex(-4.0, -0.0)]
  printf "atanh(%s) = %s\n", x, Cmath.atanh(x)
end

●実行結果

$ ruby test_cmath.rb
exp(1.0+0.0i) = 2.718281828459045+0.0i
exp(0.0+1.0i) = 0.5403023058681398+0.8414709848078965i
exp(1.0+1.0i) = 1.4686939399158851+2.2873552871788423i
exp(1.0-1.0i) = 1.4686939399158851-2.2873552871788423i
log(1.0+0.0i) = 0.0+0.0i
log(0.0+1.0i) = 0.0+1.5707963267948966i
log(1.0+1.0i) = 0.3465735902799727+0.7853981633974483i
log(1.0-1.0i) = 0.3465735902799727-0.7853981633974483i
log(-1.0+0.0i) = 0.0+3.141592653589793i
log(-1.0-0.0i) = 0.0-3.141592653589793i
log(-1.0e+300+0.0i) = 690.7755278982137+3.141592653589793i
log(-1.0e+300-0.0i) = 690.7755278982137-3.141592653589793i
sin(1.0+0.0i) = 0.8414709848078965+0.0i
sin(0.0+1.0i) = 0.0+1.1752011936438014i
sin(1.0+1.0i) = 1.2984575814159773+0.6349639147847361i
sin(1.0-1.0i) = 1.2984575814159773-0.6349639147847361i
cos(1.0+0.0i) = 0.5403023058681398-0.0i
cos(0.0+1.0i) = 1.5430806348152437-0.0i
cos(1.0+1.0i) = 0.8337300251311491-0.9888977057628651i
cos(1.0-1.0i) = 0.8337300251311491+0.9888977057628651i
tan(1.0+0.0i) = 1.557407724654902+0.0i
tan(0.0+1.0i) = 0.0+0.761594155955765i
tan(1.0+1.0i) = 0.27175258531951174+1.0839233273386948i
tan(1.0-1.0i) = 0.27175258531951174-1.0839233273386948i
sinh(1.0+0.0i) = 1.1752011936438014+0.0i
sinh(0.0+1.0i) = 0.0+0.8414709848078965i
sinh(1.0+1.0i) = 0.6349639147847361+1.2984575814159773i
sinh(1.0-1.0i) = 0.6349639147847361-1.2984575814159773i
cosh(1.0+0.0i) = 1.5430806348152437+0.0i
cosh(0.0+1.0i) = 0.5403023058681398+0.0i
cosh(1.0+1.0i) = 0.8337300251311491+0.9888977057628651i
cosh(1.0-1.0i) = 0.8337300251311491-0.9888977057628651i
tanh(1.0+0.0i) = 0.761594155955765+0.0i
tanh(0.0+1.0i) = 0.0+1.557407724654902i
tanh(1.0+1.0i) = 1.0839233273386948+0.27175258531951174i
tanh(1.0-1.0i) = 1.0839233273386948-0.27175258531951174i
sqrt(1.0+0.0i) = 1.0+0.0i
sqrt(0.0+1.0i) = 0.7071067811865476+0.7071067811865476i
sqrt(1.0+1.0i) = 1.09868411346781+0.4550898605622274i
sqrt(1.0-1.0i) = 1.09868411346781-0.4550898605622274i
sqrt(-1.0+0.0i) = 0.0+1.0i
sqrt(-1.0-0.0i) = 0.0-1.0i
sqrt(-1.0e+300+0.0i) = 0.0+1.0e+150i
sqrt(-1.0e+300-0.0i) = 0.0-1.0e+150i
asin(1.0+0.0i) = 1.5707963267948966-0.0i
asin(0.0+1.0i) = 0.0+0.8813735870195423i
asin(1.0+1.0i) = 0.6662394324925153+1.0612750619050355i
asin(1.0-1.0i) = 0.6662394324925153-1.0612750619050357i
acos(1.0+0.0i) = 0.0+0.0i
acos(0.0+1.0i) = 1.5707963267948966-0.8813735870195423i
acos(1.0+1.0i) = 0.9045568943023813-1.0612750619050355i
acos(1.0-1.0i) = 0.9045568943023813+1.0612750619050357i
atan(1.0+0.0i) = 0.7853981633974483+0.0i
atan(0.0+1.0i) = 0.0+Infinity*i
atan(1.0+1.0i) = 1.0172219678978514+0.4023594781085251i
atan(1.0-1.0i) = 1.0172219678978514-0.4023594781085251i
asin(2.0+0.0i) = 1.5707963267948966+1.3169578969248164i
asin(2.0-0.0i) = 1.5707963267948966-1.3169578969248166i
acos(2.0+0.0i) = 0.0-1.3169578969248164i
acos(2.0-0.0i) = 0.0+1.3169578969248166i
asin(-4.0+0.0i) = -1.5707963267948966+2.0634370688955617i
asin(-4.0-0.0i) = -1.5707963267948966-2.0634370688955608i
acos(-4.0+0.0i) = 3.141592653589793-2.0634370688955617i
acos(-4.0-0.0i) = 3.141592653589793+2.0634370688955608i
atan(0.0+2.0i) = 1.5707963267948966+0.5493061443340549i
atan(-0.0+2.0i) = -1.5707963267948966+0.5493061443340549i
atan(0.0-4.0i) = 1.5707963267948966-0.25541281188299525i
atan(-0.0-4.0i) = -1.5707963267948966-0.25541281188299525i
asinh(1.0+0.0i) = 0.8813735870195429+0.0i
asinh(0.0+1.0i) = 0.0+1.5707963267948966i
asinh(1.0+1.0i) = 1.0612750619050357+0.6662394324925153i
asinh(1.0-1.0i) = 1.0612750619050357-0.6662394324925153i
acosh(1.0+0.0i) = 0.0+0.0i
acosh(0.0+1.0i) = 0.8813735870195432+1.5707963267948966i
acosh(1.0+1.0i) = 1.0612750619050357+0.9045568943023813i
acosh(1.0-1.0i) = 1.0612750619050357-0.9045568943023813i
atanh(1.0+0.0i) = Infinity+0.0i
atanh(0.0+1.0i) = 0.0+0.7853981633974483i
atanh(1.0+1.0i) = 0.4023594781085251+1.0172219678978514i
atanh(1.0-1.0i) = 0.4023594781085251-1.0172219678978514i
asinh(0.0+2.0i) = 1.3169578969248166+1.5707963267948966i
asinh(-0.0+2.0i) = -1.3169578969248164+1.5707963267948966i
asinh(0.0-4.0i) = 2.0634370688955608-1.5707963267948966i
asinh(-0.0-4.0i) = -2.0634370688955617-1.5707963267948966i
acosh(0.0+0.0i) = 0.0+1.5707963267948966i
acosh(0.0-0.0i) = 0.0-1.5707963267948966i
acosh(-4.0+0.0i) = 2.0634370688955608+3.141592653589793i
acosh(-4.0-0.0i) = 2.0634370688955608-3.141592653589793i
atanh(2.0+0.0i) = 0.5493061443340549+1.5707963267948966i
atanh(2.0-0.0i) = 0.5493061443340549-1.5707963267948966i
atanh(-4.0+0.0i) = -0.25541281188299525+1.5707963267948966i
atanh(-4.0-0.0i) = -0.25541281188299525-1.5707963267948966i

Copyright (C) 2023 Makoto Hiroi
All rights reserved.

[ PrevPage | R u b y | NextPage ]