今回は引数に複素数を受け取る数学関数 (複素関数) を取り上げます。複素関数は .NET のライブラリ System.Numerics.Complex に用意されていますが、私達でもプログラムすることができます。ここでは、拙作のページ F# 入門: 「複素数」で作成した構造体 Complex に複素関数を追加することにします。なお、複素数の四則演算は float と Complex でも可能なように多重定義しています。詳細はプログラムリストをお読みくださいませ。
複素数の対数関数は複素数 z を絶対値 \(|z|\) と偏角 \(\theta\) を使って導くことができます。
複素数 x, y のべき乗 \(x^y\) は次式で求めることができます。
これをそのままプログラムすると次のようになります。
リスト : 指数関数、対数関数、べき乗 // 指数関数 Exp(z) static member Exp(z: Complex) = let a = exp(z.Real) Complex(a * cos(z.Imag), a * sin(z.Imag)) // 対数関数 Log(z) static member Log(z: Complex) = Complex(log(z.Magnitude()), z.Phase()) // べき乗 Pow(x, y) static member Pow(x: Complex, y: Complex) = Complex.Exp(y * Complex.Log(x))
簡単な例を示しましょう。
> let pi = System.Math.PI;; val pi: float = 3.141592654 > Complex.Exp(Complex(0.0, pi/4.0));; val it: Complex = complex(0.707107, 0.707107) { ... } > Complex.Exp(Complex(0.0, pi/2.0));; val it: Complex = complex(6.12323e-17, 1) { ... } > Complex.Exp(Complex(0.0, pi));; val it: Complex = complex(-1, 1.22465e-16) { ... } > Complex.Exp(Complex(1.0, 1.0));; val it: Complex = complex(1.46869, 2.28736) { ... } > Complex.Log(Complex(1.0, 1.0));; val it: Complex = complex(0.346574, 0.785398) { ... } > Complex.Log(Complex(1.0, 0.0));; val it: Complex = complex(0, 0) { ... } > Complex.Log(Complex(0.0, 1.0));; val it: Complex = complex(0, 1.5708) { ... } > Complex.Log(Complex(1.0, -1.0));; val it: Complex = complex(0.346574, -0.785398) { ... } > Complex.Log(Complex(1e300, 1e300));; val it: Complex = complex(691.122, 0.785398) { ... } > let a = Complex(1, 1);; val a: Complex = complex(1, 1) > Complex.Pow(a, Complex(0, 0));; val it: Complex = complex(1, 0) { ... } > Complex.Pow(a, Complex(1, 0));; val it: Complex = complex(1, 1) { ... } > Complex.Pow(a, Complex(2, 0));; val it: Complex = complex(1.22465e-16, 2) { ... } > Complex.Pow(a, Complex(3, 0));; val it: Complex = complex(-2, 2) { ... } > Complex.Pow(a, a);; val it: Complex = complex(0.273957, 0.583701) { ... } > Complex.Pow(Complex(1, 2), Complex(3, 4));; val it: Complex = complex(0.12901, 0.0339241) { ... }
関数 \(\log z \ (z = x + iy)\) は負の実軸 (\(-\infty \lt x \lt 0\)) において、x + 0.0j と x - 0.0j では値が異なります。
> Complex.Log(Complex(-1.0, 0.0));; val it: Complex = complex(0, 3.14159) { ... } > Complex.Log(Complex(-1.0, -0.0));; val it: Complex = complex(0, -3.14159) { ... } > Complex.Log(Complex(-1e300, 0.0));; val it: Complex = complex(690.776, 3.14159) { ... } > Complex.Log(Complex(-1e300, -0.0));; val it: Complex = complex(690.776, -3.14159) { ... }
このように、関数 \(\log z\) は負の実軸上で 2 つの値を持ちます。数学では値を一つ返す関数を「一価関数」、複数の値を返す関数を「多価関数」といいます。ここで、定義域を制限することで多価関数を一価関数にみなすことを考えます。関数 \(\log z\) の場合、負の実軸を定義域から取り除けば、\(\log z\) を一価関数とみなすことができるわけです。
参考 URL 『逆双曲線関数と逆三角関数の branch cut | 雑記帳』によると、この取り除いた領域を branch cut と呼ぶそうです。プログラミングでは branch cut を定義域から取り除くのではなく、その領域では不連続な関数とするそうです。モジュール cmath のドキュメントには、\(\log\) は「分枝切断を一つもち、0 から負の実数軸に沿って \(-\infty\) へと延びており」と記述されています。Python のマニュアルに合わせて、本稿でも branch cut を「分枝切断」と記述することにします。
プログラミング言語の場合、0.0 と -0.0 を区別する処理系であれば、Python のように 2 つの値を区別することができます。0.0 と -0.0 を区別しない処理系では、偏角 \(\theta\) の範囲を \(-\pi \lt \theta \leq \pi\) に制限することで、\(\log z\) の返り値を (\(-\pi\) を取り除いて) 一つにすることができます。
複素数の三角関数の定義は、オイラーの公式から導かれる式の \(\theta\) を複素数 \(z\) に変えたものになります。
\(\sin z, \cos z\) に純虚数 \(ix\) を与えると双曲線関数 (\(\sinh x, \cosh x\)) になります。
これに三角関数の加法定理 [*1] を使うと次の式が導かれます。
これをそのままプログラムすると、次のようになります。
リスト : 三角関数 static member Sin(z: Complex) = let x = z.Real let y = z.Imag Complex(sin(x) * cosh(y), cos(x) * sinh(y)) static member Cos(z: Complex) = let x = z.Real let y = z.Imag Complex(cos(x) * cosh(y), - sin(x) * sinh(y)) static member Tan(z: Complex) = let x = z.Real let y = z.Imag let d = cos(2.0 * x) + cosh(2.0 * y) Complex(sin(2.0 * x) / d, sinh(2.0 * y) / d)
簡単な実行例を示します。
> let a = [ Complex(0, 1); Complex(0, -1); Complex(1, 1); Complex(1, -1)];; val a: Complex list = [complex(0, 1); complex(0, -1); complex(1, 1); complex(1, -1)] > for x in a do printfn "%A" (Complex.Sin(x));; complex(0, 1.1752) complex(0, -1.1752) complex(1.29846, 0.634964) complex(1.29846, -0.634964) val it: unit = () > for x in a do printfn "%g" (Complex.Sin(x).Magnitude());; 1.1752 1.1752 1.4454 1.4454 val it: unit = () > for x in a do printfn "%A" (Complex.Cos(x));; complex(1.54308, -0) complex(1.54308, 0) complex(0.83373, -0.988898) complex(0.83373, 0.988898) val it: unit = () > for x in a do printfn "%g" (Complex.Cos(x).Magnitude());; 1.54308 1.54308 1.29345 1.29345 val it: unit = () > for x in a do printfn "%A" (Complex.Tan(x));; complex(0, 0.761594) complex(0, -0.761594) complex(0.271753, 1.08392) complex(0.271753, -1.08392) val it: unit = () > for x in a do printfn "%g" (Complex.Tan(x).Magnitude());; 0.761594 0.761594 1.11747 1.11747 val it: unit = ()
複素数の双曲線関数の定義は、実数の定義で引数 x を複素数 z に変えたものになります。
sinh z, cosh z に純虚数 ix を与えると三角関数 (sin x, cos x) になります。
これに双曲線関数の加法定理を使うと、次の式が導かれます。
これをそのままプログラムすると、次のようになります。
リスト : 双曲線関数 static member Sinh(z: Complex) = let x = z.Real let y = z.Imag Complex(sinh(x) * cos(y), cosh(x) * sin(y)) static member Cosh(z: Complex) = let x = z.Real let y = z.Imag Complex(cosh(x) * cos(y), sinh(x) * sin(y)) static member Tanh(z: Complex) = let x = z.Real let y = z.Imag let d = cosh(2.0 * x) + cos(2.0 * y) Complex(sinh(2.0 * x) / d, sin(2.0 * y) / d)
簡単な使用例を示します。
> let a = [ Complex(0, 1); Complex(0, -1); Complex(1, 1); Complex(1, -1)];; val a: Complex list = [complex(0, 1); complex(0, -1); complex(1, 1); complex(1, -1)] > for x in a do printfn "%A" (Complex.Sinh(x));; complex(0, 0.841471) complex(0, -0.841471) complex(0.634964, 1.29846) complex(0.634964, -1.29846) val it: unit = () > for x in a do printfn "%A" (Complex.Cosh(x));; complex(0.540302, 0) complex(0.540302, -0) complex(0.83373, 0.988898) complex(0.83373, -0.988898) val it: unit = () > for x in a do printfn "%A" (Complex.Tanh(x));; complex(0, 1.55741) complex(0, -1.55741) complex(1.08392, 0.271753) complex(1.08392, -0.271753) val it: unit = ()
複素数 z の平方根は次の式で求めることができます。
式 (1) を平方根の主値といいます。角度は \(2\pi\) を足すと同じ角度になるので、式 (2) がもう一つの解になります。三角関数の半角の公式を使うと、式 (1) から次の式が導かれます。
これをそのままプログラムすると、次のようになります。
リスト : 平方根 static member Sqrt(z: Complex z) = let x = z.Real; let a = z.Magnitude(); let b = sqrt((a - x) / 2.0); Complex(sqrt((a + x) / 2.0), if System.Double.IsNegative(z.Imag)) then -b else b)
虚部の符号を判定するとき IsNegative を使うと、-0.0 にも対応することができます。
簡単な実行例を示します。
> let one = Complex(1, 0);; val one: Complex = complex(1, 0) > Complex.Sqrt(one);; val it: Complex = complex(1, 0) { ... } > Complex.Sqrt(one * Complex(2, 0));; val it: Complex = complex(1.41421, 0) { ... } > Complex.Sqrt(one * Complex(3, 0));; val it: Complex = complex(1.73205, 0) { ... } > Complex.Sqrt(one * Complex(-1, 0));; val it: Complex = complex(0, 1) { ... } > Complex.Sqrt(one * Complex(-2, 0));; val it: Complex = complex(0, 1.41421) { ... } > Complex.Sqrt(one * Complex(-3, 0));; val it: Complex = complex(0, 1.73205) { ... } > let a = [ Complex(0, 1); Complex(0, -1); Complex(1, 1); Complex(1, -1)];; val a: Complex list = [complex(0, 1); complex(0, -1); complex(1, 1); complex(1, -1)] > for x in a do printfn "%A" (Complex.Sqrt(x));; complex(0.707107, 0.707107) complex(0.707107, -0.707107) complex(1.09868, 0.45509) complex(1.09868, -0.45509) val it: unit = () > for x in a do - let y = Complex.Sqrt(x) - printfn "%A" (y * y);; complex(0, 1) complex(0, -1) complex(1, 1) complex(1, -1) val it: unit = ()
\(\sqrt x\) は \(\log x\) と同じ分枝切断を持っています。x を負の整数とすると sqrt(x) の解は \(i \sqrt x\) になりますが、もうひとつ \(-i \sqrt x\) という解があります。
> Complex.Sqrt(Complex(-2.0, 0.0));; val it: Complex = complex(0, 1.41421) { ... } > Complex.Sqrt(Complex(-2.0, -0.0));; val it: Complex = complex(0, -1.41421) { ... } > Complex.Sqrt(Complex(-1e300, 0.0));; val it: Complex = complex(0, 1e+150) { ... } > Complex.Sqrt(Complex(-1e300, -0.0));; val it: Complex = complex(0, -1e+150) { ... }
三角関数の逆関数を「逆三角関数 (inverse trigonometric function)」といいます。以下に F# に用意されている逆三角関数を示します。
ここでは引数 x を実数とします。asin x は引数 x が与えられたとき sin w = x となる角度 w を求めます。同様に acos x は cos w = x となる角度 w を、atan x は tan x = w となる角度 w を求めます。三角関数には周期性があるので、上式を満たす角度 w は無数に存在します。つまり、逆三角関数の返り値は無数にあることになりますが、通常は一つの値を返すように範囲を制限します。これを「主値」といいます。
逆三角関数の主値を以下に示します。
本ページでは、数式の表示に JavaScript のライブラリ MathJax を使っています。MathJax では、逆三角関数を arcsin, arccos, arctan と表示します。
簡単な実行例を示します。
> for x in [-1.0; -0.5; 0.0; 0.5; 1.0] do printfn "%g" (asin(x));; -1.5708 -0.523599 0 0.523599 1.5708 val it: unit = () > for x in [-1.0; -0.5; 0.0; 0.5; 1.0] do printfn "%g" (acos(x));; 3.14159 2.0944 1.5708 1.0472 0 val it: unit = () > for x in [-1e300; -1.0; 0.0; 1.0; 1e300] do printfn "%g" (atan(x));; -1.5708 -0.785398 0 0.785398 1.5708 val it: unit = ()
複素数の逆三角関数の定義は、複素数の三角関数の定義から導くことができます。\(\arcsin z\) の定義は次のようになります。
\(\arccos z, \arctan z\) は定義だけを示します。
asin, acos, atan は次に示す分枝切断を持っています。
プログラムは次のようになります。
リスト : 逆三角関数 static member Asin(z: Complex) = let z1 = Complex.Sqrt(1.0 + z) let z2 = Complex.Sqrt(1.0 - z) let i = Complex(0, 1) -1.0 * i * Complex.Log(i * z + z1 * z2) static member Acos(z: Complex) = System.Math.PI / 2.0 - Complex.Asin(z) static member Atan(z: Complex) = let i = Complex(0, 1) let y = i * z let x = i * (Complex.Log(1.0 - y) - Complex.Log(1.0 + y)) / 2.0 if z.Real = 0.0 && System.Double.IsNegative(z.Real) && z.Imag > 1.0 then Complex(-x.Real, x.Imag) else x
Atan の場合、公式のままプログラムすると分枝切断 (虚軸上の 1 < ∞ ) で不具合 [*2] が生じます。そのため、z の実部が負のゼロかつ虚部が 1.0 よりも大きいときは、計算結果 x の実部の符号を反転しています。
簡単な実行例を示します。
> let a = [ Complex(0, 1); Complex(0, -1); Complex(1, 1); Complex(1, -1)];; val a: Complex list = [complex(0, 1); complex(0, -1); complex(1, 1); complex(1, -1)] > for x in a do Complex.Asin(x) |> printfn "%A";; complex(0, 0.881374) complex(0, -0.881374) complex(0.666239, 1.06128) complex(0.666239, -1.06128) val it: unit = () > for x in a do Complex.Asin(x) |> Complex.Sin |> printfn "%A";; complex(0, 1) complex(0, -1) complex(1, 1) complex(1, -1) val it: unit = () > Complex.Asin(Complex(2.0, 0.0));; val it: Complex = complex(1.5708, 1.31696) { ... } > Complex.Asin(Complex(2.0, -0.0));; val it: Complex = complex(1.5708, -1.31696) { ... } > Complex.Asin(Complex(-4.0, 0.0));; val it: Complex = complex(-1.5708, 2.06344) { ... } > Complex.Asin(Complex(-4.0, -0.0));; val it: Complex = complex(-1.5708, -2.06344) { ... } > for x in a do Complex.Acos(x) |> printfn "%A";; complex(1.5708, -0.881374) complex(1.5708, 0.881374) complex(0.904557, -1.06128) complex(0.904557, 1.06128) val it: unit = () > for x in a do Complex.Acos(x) |> Complex.Cos |> printfn "%A";; complex(8.65956e-17, 1) complex(8.65956e-17, -1) complex(1, 1) complex(1, -1) val it: unit = () > Complex.Acos(Complex(2.0, 0.0));; val it: Complex = complex(0, -1.31696) { ... } > Complex.Acos(Complex(2.0, -0.0));; val it: Complex = complex(0, 1.31696) { ... } > Complex.Acos(Complex(-4.0, 0.0));; val it: Complex = complex(3.14159, -2.06344) { ... } > Complex.Acos(Complex(-4.0, -0.0));; val it: Complex = complex(3.14159, 2.06344) { ... } > for x in a do Complex.Atan(x) |> printfn "%A";; complex(NaN, Infinity) complex(NaN, -Infinity) complex(1.01722, 0.402359) complex(1.01722, -0.402359) val it: unit = () > for x in a do Complex.Atan(x) |> Complex.Tan |> printfn "%A";; complex(NaN, NaN) complex(NaN, NaN) complex(1, 1) complex(1, -1) val it: unit = () > Complex.Atan(Complex(0.0, 2.0));; val it: Complex = complex(1.5708, 0.549306) { ... } > Complex.Atan(Complex(-0.0, 2.0));; val it: Complex = complex(-1.5708, 0.549306) { ... } > Complex.Atan(Complex(0.0, -4.0));; val it: Complex = complex(1.5708, -0.255413) { ... } > Complex.Atan(Complex(-0.0, -4.0));; val it: Complex = complex(-1.5708, -0.255413) { ... }
(0.0 + 1.0 i) * (0.0 + b i) = 0.0 * 0.0 - 1.0 * b + (0.0 * b + 1.0 * 0.0) = -b + (0.0 + 0.0) i = -b + 0.0i (0.0 + 1.0 i) * (-0.0 + b i) = 0.0 * -0.0 - 1.0 * b + (0.0 * b + 1.0 * -0.0) = -b + (0.0 + -0.0) i = -b + 0.0iIEEE 754 では、-0.0 + -0.0 は -0.0 になるが、0.0 + -0.0 や -0.0 + 0.0 は 0.0 になる。
> 0.0 + 0.0;; val it: float = 0.0 > -0.0 + 0.0;; val it: float = 0.0 > 0.0 + -0.0;; val it: float = 0.0 > -0.0 + -0.0;; val it: float = -0.0
双曲線関数の逆関数を「逆双曲線関数 (inverse hyperbolic function)」といいます。F# には逆双曲線関数が定義されていませんが、.NET のライブラリ System.Math には逆双曲線関数 Asinh, Acosh, Atanh があります。MathJax では逆双曲線関数を \(\sinh^{-1} x, \cosh^{-1} x, \tanh^{-1} x \) で表します。双曲線関数と逆双曲線関数の定義域と値域を示します。
x と y は実数です。x = cosh y の逆関数 y = acosh x を満たす y の値は 2 つありますが、ここでは \(y \geq 0\) を主値として選ぶことにします。
逆双曲線関数の定義は双曲線関数の定義から導くことができます。
関数 \(\tanh^{-1} x\) を使うと、次の関係式から \(\sinh^{-1} x, \ \cosh^{-1} x\) を求めることができます。
簡単な実行例を示します。
> for x in [1.5; 1.0; 0.5; 0.0; -0.5; -1.0; -1.5] do System.Math.Asinh(x) |> printfn "%g";; 1.19476 0.881374 0.481212 0 -0.481212 -0.881374 -1.19476 val it: unit = () > for x in [1.5; 1.0; 0.5; 0.0; -0.5; -1.0; -1.5] do System.Math.Asinh(x) |> sinh |> printfn "%g";; 1.5 1 0.5 0 -0.5 -1 -1.5 val it: unit = () > for x in [1.0; 1.5; 2.0; 2.5; 3.0; 3.5; 4.0] do System.Math.Acosh(x) |> printfn "%g";; 0 0.962424 1.31696 1.5668 1.76275 1.92485 2.06344 val it: unit = () > for x in [1.0; 1.5; 2.0; 2.5; 3.0; 3.5; 4.0] do System.Math.Acosh(x) |> cosh |> printfn "%g";; 1 1.5 2 2.5 3 3.5 4 val it: unit = () > for x in [-0.75; -0.5; -0.25; 0.0; 0.25; 0.5; 0.75] do System.Math.Atanh(x) |> printfn "%g";; -0.972955 -0.549306 -0.255413 0 0.255413 0.549306 0.972955 val it: unit = () > for x in [-0.75; -0.5; -0.25; 0.0; 0.25; 0.5; 0.75] do System.Math.Atanh(x) |> tanh |> printfn "%g";; -0.75 -0.5 -0.25 0 0.25 0.5 0.75 val it: unit = ()
複素数の逆双曲線関数の定義は、実数の定義で引数 x を複素数 z に変えたものになります。
参考 URL 『逆双曲線関数と逆三角関数の branch cut | 雑記帳』によると、\(\cosh^{-1} z\) を式 (1) でプログラムすると「分枝切断線」が複雑になるため、他の式 (たとえば (2) など) でプログラムする処理系が多いようです。ちなみに、ANSI Common Lisp では \(\cosh^{-1} z\) を次の式で定義しています。
asinh, acosh, atanh は次に示す分枝切断を持っています。
プログラムは次のようになります。
リスト : 逆双曲線関数 static member Asinh(z: Complex) = Complex.Log(z + Complex.Sqrt(z * z + 1.0)) static member Acosh(z: Complex) = let z1 = Complex.Sqrt(z + 1.0) let z2 = Complex.Sqrt(z - 1.0) Complex.Log(z + z1 * z2) static member Atanh(z: Complex) = (Complex.Log(1.0 + z) - Complex.Log(1.0 - z)) / 2.0
簡単な実行例を示しましょう。
> let a = [ Complex(0, 1); Complex(0, -1); Complex(1, 1); Complex(1, -1)];; val a: Complex list = [complex(0, 1); complex(0, -1); complex(1, 1); complex(1, -1)] > for x in a do Complex.Asinh(x) |> printfn "%A";; complex(0, 1.5708) complex(0, -1.5708) complex(1.06128, 0.666239) complex(1.06128, -0.666239) val it: unit = () > for x in a do Complex.Asinh(x) |> Complex.Sinh |> printfn "%A";; complex(0, 1) complex(0, -1) complex(1, 1) complex(1, -1) val it: unit = () > Complex.Asinh(Complex(0.0, 2.0));; val it: Complex = complex(1.31696, 1.5708) { ... } > Complex.Asinh(Complex(-0.0, 2.0));; val it: Complex = complex(-1.31696, 1.5708) { ... } > Complex.Asinh(Complex(0.0, -4.0));; val it: Complex = complex(2.06344, -1.5708) { ... } > Complex.Asinh(Complex(-0.0, -4.0));; val it: Complex = complex(-2.06344, -1.5708) { ... } > for x in a do Complex.Acosh(x) |> printfn "%A";; complex(0.881374, 1.5708) complex(0.881374, -1.5708) complex(1.06128, 0.904557) complex(1.06128, -0.904557) val it: unit = () > for x in a do Complex.Acosh(x) |> Complex.Cosh |> printfn "%A";; complex(8.65956e-17, 1) complex(8.65956e-17, -1) complex(1, 1) complex(1, -1) val it: unit = () > Complex.Acosh(Complex(0.0, 0.0));; val it: Complex = complex(0, 1.5708) { ... } > Complex.Acosh(Complex(0.0, -0.0));; val it: Complex = complex(0, -1.5708) { ... } > Complex.Acosh(Complex(-4.0, 0.0));; val it: Complex = complex(2.06344, 3.14159) { ... } > Complex.Acosh(Complex(-4.0, -0.0));; val it: Complex = complex(2.06344, -3.14159) { ... } > for x in a do Complex.Atanh(x) |> printfn "%A";; complex(0, 0.785398) complex(0, -0.785398) complex(0.402359, 1.01722) complex(0.402359, -1.01722) val it: unit = () > for x in a do Complex.Atanh(x) |> Complex.Tanh |> printfn "%A";; complex(0, 1) complex(0, -1) complex(1, 1) complex(1, -1) val it: unit = () > Complex.Atanh(Complex(2.0, 0.0));; val it: Complex = complex(0.549306, 1.5708) { ... } > Complex.Atanh(Complex(2.0, -0.0));; val it: Complex = complex(0.549306, -1.5708) { ... } > Complex.Atanh(Complex(-4.0, 0.0));; val it: Complex = complex(-0.255413, 1.5708) { ... } > Complex.Atanh(Complex(-4.0, -0.0));; val it: Complex = complex(-0.255413, -1.5708) { ... }
// // complex.fsx : 複素数と複素関数 // // Copyright (C) 2022 Makoto Hiroi // [<Struct; CustomEquality; NoComparison>] type Complex = val Real: float val Imag: float // コンストラクタ new(a: float, b: float) = { Real = a; Imag = b } new(a: int, b: int) = { Real = float a; Imag = float b } // メソッド override this.ToString () = sprintf "complex(%g, %g)" this.Real this.Imag // 等値演算子 override this.Equals (other: obj) = let that = other :?> Complex this.Real = that.Real && this.Real = that.Real // とても適当なので実際に使ってはいけない override this.GetHashCode () = (this.Real / this.Imag).GetHashCode() // 複素共役 member this.Conjugate() = Complex(this.Real, -this.Imag) // 偏角 member this.Phase() = atan2 this.Imag this.Real // 絶対値 member this.Magnitude() = if this.Real = 0.0 then abs this.Imag else if this.Imag = 0.0 then abs this.Real else if abs this.Imag > abs this.Real then let temp = this.Real / this.Imag (abs this.Imag) * sqrt(1.0 + temp * temp) else let temp = this.Imag / this.Real (abs this.Real) * sqrt(1.0 + temp * temp) // // 四則演算 // static member ( + ) (x: Complex, y: Complex) = Complex(x.Real + y.Real, x.Imag + y.Imag) static member ( + ) (x: float, y: Complex) = Complex(x + y.Real, y.Imag) static member ( + ) (x: Complex, y: float) = Complex(x.Real + y, x.Imag) static member ( - ) (x: Complex, y: Complex) = Complex(x.Real - y.Real, x.Imag - y.Imag) static member ( - ) (x: float, y: Complex) = Complex(x - y.Real, -y.Imag) static member ( - ) (x: Complex, y: float) = Complex(x.Real - y, x.Imag) static member ( * ) (x: Complex, y: Complex) = let a = x.Real let b = x.Imag let c = y.Real let d = y.Imag Complex(a * c - b * d, b * c + a * d) static member ( * ) (x: float, y: Complex) = Complex(x * y.Real, x * y.Imag) static member ( * ) (x: Complex, y: float) = Complex(x.Real * y, x.Imag * y) static member ( / ) (x: Complex, y: Complex) = let a = x.Real let b = x.Imag let c = y.Real let d = y.Imag if abs c >= abs d then let u = d / c let v = c + d * u Complex((a + b * u) / v, (b - a * u) / v) else let u = c / d let v = c * u + d Complex((a * u + b) / v, (b * u - a) / v) static member ( / ) (x: float, y: Complex) = let c = y.Real let d = y.Imag if abs c >= abs d then let u = d / c let v = c + d * u Complex(x / v, (-x * u) / v) else let u = c / d let v = c * u + d Complex((x * u) / v, (-x) / v) static member ( / ) (x: Complex, y: float) = Complex(x.Real / y, x.Imag / y) // // 複素関数 // // 指数関数 Exp(z) static member Exp(z: Complex) = let a = exp(z.Real) Complex(a * cos(z.Imag), a * sin(z.Imag)) // 対数関数 Log(z) static member Log(z: Complex) = Complex(log(z.Magnitude()), z.Phase()) // べき乗 Pow(x, y) static member Pow(x: Complex, y: Complex) = Complex.Exp(y * Complex.Log(x)) // 三角関数 static member Sin(z: Complex) = let x = z.Real let y = z.Imag Complex(sin(x) * cosh(y), cos(x) * sinh(y)) static member Cos(z: Complex) = let x = z.Real let y = z.Imag Complex(cos(x) * cosh(y), - sin(x) * sinh(y)) static member Tan(z: Complex) = let x = z.Real let y = z.Imag let d = cos(2.0 * x) + cosh(2.0 * y) Complex(sin(2.0 * x) / d, sinh(2.0 * y) / d) // 双曲線関数 static member Sinh(z: Complex) = let x = z.Real let y = z.Imag Complex(sinh(x) * cos(y), cosh(x) * sin(y)) static member Cosh(z: Complex) = let x = z.Real let y = z.Imag Complex(cosh(x) * cos(y), sinh(x) * sin(y)) static member Tanh(z: Complex) = let x = z.Real let y = z.Imag let d = cosh(2.0 * x) + cos(2.0 * y) Complex(sinh(2.0 * x) / d, sin(2.0 * y) / d) // 平方根 static member Sqrt(z: Complex) = let x = z.Real; let a = z.Magnitude(); let b = sqrt((a - x) / 2.0); Complex(sqrt((a + x) / 2.0), if System.Double.IsNegative(z.Imag) then -b else b) // 複素数の逆三角関数 static member Asin(z: Complex) = let z1 = Complex.Sqrt(1.0 + z) let z2 = Complex.Sqrt(1.0 - z) let i = Complex(0, 1) -1.0 * i * Complex.Log(i * z + z1 * z2) static member Acos(z: Complex) = System.Math.PI / 2.0 - Complex.Asin(z) static member Atan(z: Complex) = let i = Complex(0, 1) let y = i * z let x = i * (Complex.Log(1.0 - y) - Complex.Log(1.0 + y)) / 2.0 if z.Real = 0.0 && System.Double.IsNegative(z.Real) && z.Imag > 1.0 then Complex(-x.Real, x.Imag) else x // 複素数の逆双曲線関数 static member Asinh(z: Complex) = Complex.Log(z + Complex.Sqrt(z * z + 1.0)) static member Acosh(z: Complex) = let z1 = Complex.Sqrt(z + 1.0) let z2 = Complex.Sqrt(z - 1.0) Complex.Log(z + z1 * z2) static member Atanh(z: Complex) = (Complex.Log(1.0 + z) - Complex.Log(1.0 - z)) / 2.0