第13回：複素数

■ 複素数を作る

以下では，x, y は Int64 型または Float64 型の数とする．

関数 complex(x) は，実数部が $x$ の複素数を作る．
関数 complex(x,y) は，実数部（実部）が $x$ ,虚数部（虚部）が $y$ の複素数を作る．
複素数の型は，実数部と虚数部の型に合わせて Complex{Int64} または Complex{Float64} となる．$x, y$ で Int64 と Float64 型が混在した場合には．後者になる．

julia> complex(1)1 + 0im
julia> complex(1.0)1.0 + 0.0im
julia> complex(1.0, -1.0)1.0 - 1.0im

定数 im は虚数単位である．これは，complex(0,1) または complex(0.0,1.0) と同じ意味である．複素数 $x+iy$ を作るのに，以下のように書いてもよい（が，乗算と加算の演算を含むので，関数 complex を使うほうが好ましい )．

1 - 1im
1.0 - 1.0im

▶ 複素数と整数・浮動小数点数との四則演算

複素数と整数または浮動小数点数との四則演算は，演算子 +, -, *, / を用いる．

julia> z = complex(1, -1)1 - 1im
julia> z + 23 - 1im
julia> z - 2-1 - 1im
julia> z * 22 - 2im
julia> z / 20.5 - 0.5im

零 0 による除算は，実数部と虚数部の各々で行われ，Inf ないし NaN の値となる．

julia> complex(1, 0) / 0Inf + NaN*im

▶ 複素数同士の四則演算

複素数同士の四則演算にも，演算子 +, -, *, / を用いる．

julia> z = complex(1,-1)1 - 1im
julia> w = complex(2,-2)2 - 2im
julia> z + w3 - 3im
julia> z - w-1 + 1im
julia> z * w0 - 4im
julia> z / w0.5 - 0.0im

複素数を値とする変数に対して，更新演算子 +=, -=, *=, /= も使える．

julia> z *= w0 - 4im

▶ 複素数のベクトル

整数または浮動小数点数のベクトルを作るのと同様な方法で，複素数のベクトルを作ることができる．

julia> [complex(0, 0), complex(1, 0), complex(1, 1)]3-element Vector{Complex{Int64}}:
 0 + 0im
 1 + 0im
 1 + 1im
julia> [complex(i, 2i) for i = -2:2]5-element Vector{Complex{Int64}}:
 -2 - 4im
 -1 - 2im
  0 + 0im
  1 + 2im
  2 + 4im

複素数の零 complex(0.0,0.0) を $5$ 個含むベクトルを作るには，以下のように書けばよい．

julia> zeros(Complex{Float64}, 5)5-element Vector{ComplexF64}:
 0.0 + 0.0im
 0.0 + 0.0im
 0.0 + 0.0im
 0.0 + 0.0im
 0.0 + 0.0im

▶ 複素数の実数部・虚数部・共役複素数

関数 real(z) は複素数 z の実数部（実部，real part）を返す．
関数 imag(z) は複素数 z の虚数部（虚部，imaginary part）を返す．

julia> z = complex(1, -1)1 - 1im
julia> real(z)1
julia> imag(z)-1

複素数ベクトルの各要素の実数部ないし虚数部を計算するには，dot記法を用いる．

julia> zs = [complex(0, 0), complex(1, 0), complex(1, 1), complex(0, 1), complex(0, 0)]5-element Vector{Complex{Int64}}:
 0 + 0im
 1 + 0im
 1 + 1im
 0 + 1im
 0 + 0im
julia> real.(zs)5-element Vector{Int64}:
 0
 1
 1
 0
 0
julia> imag.(zs)5-element Vector{Int64}:
 0
 0
 1
 1
 0

関数 conj(z) は複素数 z の「共役（きょうえき）複素数（conjugate complex number）」を返す．

julia> z = complex(1, -1)1 - 1im
julia> conj(z)1 + 1im

▶ 負の数に対する平方根

Base.sqrt - Function

平方根 sqrt(x)は，負の実数 x に対して例外を出すが，引数を複素数の引数を与えれば計算できる．

julia> # 例外を出す
       sqrt(-1)ERROR: DomainError with -1.0:
sqrt will only return a complex result if called with a complex argument. Try sqrt(Complex(x)).
julia> # 計算できる
       sqrt(complex(-1))0.0 + 1.0im

正負の実数に対して，平方根の実数部と虚数部を描け．

using PyPlot
plt.axes().set_aspect("equal")

xs = -1:0.05:1
zs = complex.(xs)
sqzs = sqrt.(zs)
plt.plot(xs, real.(sqzs), label = "real part")
plt.plot(xs, imag.(sqzs), label = "imaginary part")
plt.xlabel("x")
plt.ylabel("sqrt of x")
plt.legend()
plt.xlim(-1.5, 1.5)
plt.ylim(-1.5, 1.5)
plt.axhline(0, lw = 0.5, color = "k")
plt.axvline(0, lw = 0.5, color = "k")

◀ 2次方程式の解：一般の場合

2次方程式 $x^{2} -bx + c = 0$ の解を求める方法について，▼ 2次方程式で紹介したが，実数解のみに留めていた．以下では，複素数解を含めて求めてみる．一方の解 $x_{1}$ を，2次方程式の解の公式で求め，他方の解を，係数と解の関係により求める．

\[\begin{aligned} x_{1} &= \dfrac{b+\sqrt{b^2-4c}}{2}, \\ x_{2} &= \dfrac{c}{x_1} \end{aligned}\]

係数 $b=1$ を一定とし，係数 $c$ を $-1$ から $1$ の範囲で動かす．

実数部を実線で，虚数部を点線で表した．

using PyPlot
plt.axes().set_aspect("equal")

b = 1
cs = range(-2, 2, length = 200)
ds = complex(b .* b .- 4cs)
x1s = (b .+ sqrt.(ds)) / 2;
x2v = cs ./ x1s;
plt.plot(cs, real.(x1s), "b-", label = "x1, real part")
plt.plot(cs, imag.(x1s), "b:", label = "x1, imag part")
plt.plot(cs, real.(x2v), "g-", label = "x2, real part")
plt.plot(cs, imag.(x2v), "g:", label = "x2, imag part")
plt.legend()
plt.xlabel("c")
plt.xlim(-2.5, 2.5)
plt.ylim(-2.5, 2.5)
plt.axhline(0, lw = 0.5, color = "k")
plt.axvline(0, lw = 0.5, color = "k")

定数 $c$ の値により，解の形が異なる様子が観察できる．すなわち，

範囲 $c \lt \dfrac{1}{4}$ では実数解（虚数部は零 $0$ ）
範囲 $c = \dfrac{1}{4}$ では重解 $\dfrac{1}{2}$
範囲 $c \gt \dfrac{1}{4}$ では複素解（共役複素数）

▶ ガウス平面

複素数 $x + i y$ は，平面上の点 $(x,y)$ と一対一に対応する．複素数を平面に表したものをガウス平面という．

using PyPlot
plt.axes().set_aspect("equal")
z = complex(1, 1)
w = complex(-2, 1)

plt.plot(real(z), imag(z), "b.", label = "z")
plt.plot(real(w), imag(w), "r.", label = "w")
plt.legend()
plt.xlabel("real part")
plt.ylabel("imag part")
plt.xlim(-3, 3)
plt.ylim(-3, 3)
plt.axhline(0, lw = 0.5, color = "k")
plt.axvline(0, lw = 0.5, color = "k")

複素数ベクトルで表された図形を，ガウス平面に描こう．

using PyPlot
plt.axes().set_aspect("equal")

zs = [complex(0, 0), complex(1, 0), complex(1, 1), complex(0, 1), complex(0, 0)]

plt.plot(real.(zs), imag.(zs), ".-")
plt.xlabel("real part")
plt.ylabel("imag part")
plt.xlim(-2, 2)
plt.ylim(-2, 2)
plt.axhline(0, lw = 0.5, color = "k")
plt.axvline(0, lw = 0.5, color = "k")

▶ 複素数の絶対値と偏角

関数 abs(z) は，複素数 z の絶対値 $\left\vert{z}\right\vert$ を返す．
関数 abs2(z) は，複素数 z の絶対値の2乗 $\left\vert{z}\right\vert|^{2}$ を返す．
関数 angle(z) は，複素数 z の偏角 $\angle{z}$ をラジアンで返す．

結果（値域）は $-\pi$ から $\pi$ である．

これは，ガウス平面上で，複素数を極座標で表示したものに対応している．次の ▶ オイラーの公式も参照．

julia> for z in [
          0,   complex(1, 1), 1im, complex(-1, 1),
          -1, complex(-1, -1), -im, complex(1, -1), ]
          @show z, abs(z), angle(z)
       end(z, abs(z), angle(z)) = (0 + 0im, 0.0, 0.0)
(z, abs(z), angle(z)) = (1 + 1im, 1.4142135623730951, 0.7853981633974483)
(z, abs(z), angle(z)) = (0 + 1im, 1.0, 1.5707963267948966)
(z, abs(z), angle(z)) = (-1 + 1im, 1.4142135623730951, 2.356194490192345)
(z, abs(z), angle(z)) = (-1 + 0im, 1.0, 3.141592653589793)
(z, abs(z), angle(z)) = (-1 - 1im, 1.4142135623730951, -2.356194490192345)
(z, abs(z), angle(z)) = (0 - 1im, 1.0, -1.5707963267948966)
(z, abs(z), angle(z)) = (1 - 1im, 1.4142135623730951, -0.7853981633974483)

▶ オイラーの公式

実数 $\theta$ に対して，指数関数 $\exp{i\theta}$ は，以下のように書き表される．これをオイラーの公式という．

\[\exp{i\theta} = \cos\theta+ i \sin\theta\]

指数関数は，複素数を引数とするように拡張されている． $\exp{i\theta}$ をガウス平面上に描く．これは，単位円（半径 $1$ ）の円を描く．

plt.axes().set_aspect("equal")

zs = [exp(im * t) for t = 0:pi/18:2pi]

plt.plot(real.(zs), imag.(zs), ".")
plt.xlabel("real part")
plt.ylabel("imag part")
plt.axhline(0, lw = 0.5, color = "k")
plt.axvline(0, lw = 0.5, color = "k")

関数 cis(x) は，数 x に対して $\exp{i x}$ を計算する．同じ結果が得られることを確認しよう．

using PyPlot
plt.axes().set_aspect("equal")

zs = [cis(t) for t = 0:pi/18:2pi]
plt.plot(real.(zs), imag.(zs), "o")

zs = [exp(im * t) for t = 0:pi/18:2pi]
plt.plot(real.(zs), imag.(zs), ".")

plt.xlabel("real part")
plt.ylabel("imag part")

plt.axhline(0, lw = 0.5, color = "k")
plt.axvline(0, lw = 0.5, color = "k")

◀ ガウス平面内で回転する

複素数に $\exp(i \theta)$ を乗ずることは，ガウス平面上で，原点に対して，反時計方向に角度 $\theta$ だけ回転することに相当する．

using PyPlot
plt.axes().set_aspect("equal")

r15 = cis(pi * 15 / 180)
zs = complex(1, 1)

for i = 1:10
   global zs
   plt.plot([0, real.(zs)], [0, imag.(zs)], ".-")
   zs *= r15
end

plt.xlabel("real part")
plt.ylabel("imag part")
plt.axhline(0, lw = 0.5, color = "k")
plt.axvline(0, lw = 0.5, color = "k")

plt.xlim(-2, 2)
plt.ylim(-2, 2)

◀ アルキメデスの渦を描く（複素数版）

実数 $t$ に対して，複素数 $z=t \exp(i t)$ の軌跡を，ガウス平面上に描いてみよう．これは，アルキメデスの渦である．参考 → ▼ アルキメデスの渦を描く

using PyPlot
plt.axes().set_aspect("equal")

ts = range(0, 4pi, length = 200)
zs = ts .* exp.(im * ts)
plt.plot(real.(zs), imag.(zs))

●▼ 花曲線（複素数版）

実数 $t$ に対して，複素数 $z= \cos(nt) \exp(i t)$ の軌跡をガウス平面上に描いてみよ．これは，花曲線となる．参考 → ▼ 花曲線を描く

◀ 複素数同士の乗算の意味

複素数同士の乗算は，極座標で表示すると，その意味が明らかになる．

オイラーの公式を使うと，複素数 $z, w$ は，その絶対値と偏角を用いて，以下のように書いて，

\[\begin{aligned} z &= \left\vert{z}\right\vert \exp(i \angle{z}), \\ w &= \left\vert{w}\right\vert \exp(i \angle{w}) \end{aligned}\]

複素数 $z$ と $w$ の積を求めると，次のようになる．

\[zw = \left\vert{z}\right\vert \left\vert{w}\right\vert \exp\left( i \left(\angle{z}+\angle{w}\right)\right)\]

つまり，積 $zw$ の絶対値は，2つの複素数の絶対値の積である．また，積 $zw$ の偏角は，2つの複素数の偏角の和である．

まとめると，複素数同士の積は，ガウス平面上で拡大縮小と回転を同時に行う演算である．

using PyPlot
plt.axes().set_aspect("equal")

z = 1 * cis(pi / 3)
w = 2 * cis(pi / 4)
zw = z * w

plt.plot([0, real(z)], [0, imag(z)], "r.-", label = "z")
plt.plot([0, real(w)], [0, imag(w)], "b.-", label = "w")
plt.plot([0, real(zw)], [0, imag(zw)], "g.-", label = "z*w")

zs = [cis(t) for t = 0:pi/18:2pi]
plt.plot(real.(zs), imag.(zs), "r", lw = 0.5)
plt.plot(real.(2zs), imag.(2zs), "g", lw = 0.5)

plt.legend()
plt.xlabel("real part")
plt.ylabel("imag part")

plt.axhline(0, lw = 0.5, color = "k")
plt.axvline(0, lw = 0.5, color = "k")
plt.xlim(-3, 3)
plt.ylim(-3, 3)

◀ 複素数の平方根とは

上の特別な場合として，$z$ の二乗を検討する． $w=z$ として，以下を得る．

\[z^2 = {\left\vert{z}\right\vert}^2 \exp\left( i 2 \angle{z} \right)\]

これから，$z$ の平方根 $\sqrt{z}$ は，以下のように求められる．

\[z = \sqrt{\left\vert{z}\right\vert} \exp\left( i \dfrac{\angle{z}}{2} \right)\]

すなわち，複素数の平方根は，

絶対値が，元の複素数の絶対値の平方根
偏角が，元の複素数の偏角の半分

となる．

以下に，数値例を示す．

\[\begin{aligned} z & = -1 + i \sqrt{3} = 2 \left( -\dfrac{1}{2} + i \dfrac{\sqrt{3}}{2} \right) = 2 \exp\left(i \dfrac{\pi}{3}\right), \\ \sqrt{z} & = \sqrt{2} \exp\left(i \dfrac{\pi}{6} \right) = \sqrt{2} \left( \dfrac{\sqrt{3}}{2} + i \dfrac{1}{2} \right) \end{aligned}\]

z = complex(-1, sqrt(3));
@show z, abs(z), angle(z);
w = sqrt(z);
@show w, abs(w), angle(w);

(0.7071067811865476 + 1.224744871391589im, 1.414213562373095, 1.0471975511965976)

using PyPlot
plt.axes().set_aspect("equal")

z = complex(-1, sqrt(3))
w = sqrt(z)

plt.plot([0, real(z)], [0, imag(z)], "r.-", label = "z")
plt.plot([0, real(w)], [0, imag(w)], "b.-", label = "sqrt(z)")

zs = [cis(t) for t = 0:pi/18:2pi]
plt.plot(real.(zs) * 2, imag.(zs) * 2, "r", lw = 0.5)
plt.plot(real.(zs) * sqrt(2), imag.(zs) * sqrt(2), "b", lw = 0.5)

plt.legend()
plt.xlabel("real part")
plt.ylabel("imag part")

plt.axhline(0, lw = 0.5, color = "k")
plt.axvline(0, lw = 0.5, color = "k")
plt.xlim(-3, 3)
plt.ylim(-3, 3)

◀ 図形を回転する

複素数に $\exp(i \theta)$ を乗ずることは，ガウス平面上で，原点に対して，反時計方向に角度 $\theta$ だけ回転することに相当する．平面図形を複素数ベクトルとして表して，図形を回転しよう．

using PyPlot
plt.axes().set_aspect("equal")

r15 = cis(pi * 15 / 180)
zs = [complex(0, 0), complex(1, 0), complex(1, 2), complex(0, 0)]

for i = 1:10
   global zs
   plt.plot(real.(zs), imag.(zs), ".-")
   zs *= r15
end

plt.xlabel("real part")
plt.ylabel("imag part")

plt.xlim(-3, 3)
plt.ylim(-3, 3)
plt.axhline(0, lw = 0.5, color = "k")
plt.axvline(0, lw = 0.5, color = "k")

◀ 伝達関数

制御工学や電気電子回路では「線形システム」に着目する．これは，入力と出力が，時刻 $t$ に対する定数係数の微積分方程式で表されるようなものである．線形システムでは，

角振動数 $\omega=2\pi{f}$ の正弦波 $x(t) = x_0 \cos \left(\omega{t}+\phi \angle{X}\right)$ を入力に与えると，
同じ角振動数の正弦波 $y(t) = y_0 \cos \left(\omega{t}+\angle{Y} \right)$ が出力として得られること

が知られている．

そこで，

複素数 $X = x_0 \exp(i\angle{X})$ を入力の正弦波 $x(t)$ と同一視し（＝同じものと考え），
複素数 $Y = y_0 \exp(i\angle{Y})$ を出力の正弦波 $y(t)$ と同一視すると，

入力と出力との比は，複素数 $H$ になる．

\[Y = H X\]

ここで，$X, Y, H$ は，一般に，角振動数 $\omega$ または周波数 $f$ の関数である．

\[Y({\omega}) = H({\omega}) X({\omega})\]

複素数 $H(\omega)$ を，「伝達関数（transfer function）」という．

伝達関数の絶対値は，入力と出力の振幅の比（振幅比）を与える．
伝達関数の偏角は，入力の位相と出力の位相の「ズレ」（位相差）を与える．

\[\begin{aligned} \left\vert{H({\omega})}\right\vert &= \dfrac{\left\vert{Y({\omega})}\right\vert}{\left\vert{X({\omega})}\right\vert} = \dfrac{y_0}{x_0}, \\ \angle{H(\omega)} &= \angle{Y(\omega)} - \angle{X(\omega)} \end{aligned}\]

「周波数応答（frequency response）」とは，角振動数 $\omega$ または周波数 $f$ の関数として伝達関数の振幅 $\left\vert{H({\omega})}\right\vert$ と位相 $\angle{H(\omega)}$ を描いた図である．線形システムの振る舞いを観察するのに便利である．

◀ 共振回路

自己インダクタンス（コイル） $L$, 電気容量（コンデンサ） $C$, 電気抵抗 $R$ を直列に接続した $LCR$ 直列回路の電源に正弦波電圧 $e(t)$ を加える．抵抗 $R$ の両端子間の電圧 $v(t)$ は，同じ角振動数を持つ正弦波 $v(t)$ となる．正弦波 $e(t)$ と $v(t)$ を複素数 $E(\omega)$ と $V(\omega)$ で表したとき，両者の比は，以下の複素数 $H(\omega)$ で表される．

\[H(\omega) = \dfrac{V(\omega)}{E(\omega)} = \dfrac{1}{1 + i Q \left(\dfrac{\omega}{\omega_{0}} - \dfrac{\omega_{0}}{\omega}\right)}\]

ここで，$\omega_{0} = 2\pi{f_0}$ と $Q$ は，回路素子 $L, C, R$ の値から決まる正の定数である．$f_{0}$ は周波数の次元，$\omega_0$ は角振動数の次元を持ち，$Q$ は無次元である．

まず，$f_0= 1\;\mathrm{kHz}, Q = 1$ として，$H(\omega)$ の振幅 $\left\vert{H(\omega)}\right\vert$ を描く．横軸周波数は対数で表示する．

f0 = 1e3
w0 = 2pi * f0

fs = exp10.(range(1, 5, length = 200))

ws = 2pi * fs

q = 1
h1 = 1 ./ (1 .+ im * q * (ws / w0 - w0 ./ ws))

using PyPlot
plt.plot(fs, abs.(h1))
plt.xscale("log")
plt.ylabel("Amplitude")
plt.xlabel("f / Hz")
plt.axvline(f0, lw = 0.5, color = "k")

振幅は，単峰性の極大値 $1$ をとる．

極大となる周波数 $f_{0}$ を「共振周波数（resonance frequency）」，これに対応する角振動数 $\omega_{0}=2\pi{f_0}$ を「共振角振動数（resonance angular frequency）」という．

今度は，振幅の二乗 $\left\vert{H(\omega)}\right\vert^2$ と位相 $\angle{H(\omega)}$ を同時に描こう．

f0 = 1e3
w0 = 2pi * f0

fs = exp10.(range(1, 5, length = 200))

ws = 2pi * fs

q = 1
h1 = 1 ./ (1 .+ im * q * (ws / w0 - w0 ./ ws))

using PyPlot
fig, axs = plt.subplots(2, 1)
ax1 = axs[1]
ax1.plot(fs, abs2.(h1))
ax1.set_ylabel("Amplitude")
ax1.set_xscale("log")
ax1.axvline(f0, lw = 0.5, color = "k")
ax1.axhline(1 / sqrt(2), lw = 0.5, color = "k")

ax2 = axs[2]
ax2.plot(fs, angle.(h1) * 180 / pi)

ax2.set_xscale("log")
ax2.set_xlabel("f / Hz")
ax2.set_ylabel("Phase")
ax2.set_ylim(-100, 100)
ytics = [-90, -45, 0, 45, 90]
ax2.set_yticks(ytics)
for y in ytics
   ax2.axhline(y, lw = 0.5, color = "k")
end
ax2.axvline(f0, lw = 0.5, color = "k")

位相は $90^{\circ}$ から始まり $-90^{\circ}$ に単調減少する．振幅が極大となる周波数 $f_0$ で位相は $0$ となる．

さらに，$\left\vert{H(\omega)}\right\vert^2 = \dfrac{1}{2}$ となる周波数（２つある）で，位相は $\pm 45^{\circ}$ をとる．

次に，$f_0, \omega_0$ を変えずに，$Q$ の値を増減して，振幅の二乗 $\left\vert{H(\omega)}\right\vert^2$ を描く．

f0 = 1e3
w0 = 2pi * f0

fs = exp10.(range(1, 5, length = 200))

ws = 2pi * fs

using PyPlot
for q in [0.5, 1, 2]
   h1 = 1 ./ (1 .+ im * q * (ws / w0 - w0 ./ ws))
   plt.plot(fs, abs2.(h1), label = "Q=" * string(q))
end
plt.legend()
plt.xscale("log")
plt.xlabel("f / Hz")

plt.axhline(1 / sqrt(2), lw = 0.5, color = "k")
plt.axvline(f0, lw = 0.5, color = "k")

┌ Warning: Assignment to `h1` in soft scope is ambiguous because a global variable by the same name exists: `h1` will be treated as a new local. Disambiguate by using `local h1` to suppress this warning or `global h1` to assign to the existing global variable.
└ @ ch13.md:691

定数 $Q$ が小さくなると，峰が鋭くなる（幅が狭くなる）様子が観察される．

更に，位相のグラフを加えよう．

f0 = 1e3
w0 = 2pi * f0

fs = exp10.(range(1, 5, length = 200))

ws = 2pi * fs

using PyPlot
fig, axs = plt.subplots(2, 1)
ax1 = axs[1]
ax2 = axs[2]

for q in [0.5, 1, 2]
   h1 = 1 ./ (1 .+ im * q * (ws / w0 - w0 ./ ws))
   ax1.plot(fs, abs.(h1), label = "Q=" * string(q))
   ax2.plot(fs, angle.(h1) * 180 / pi)
end
ax1.legend()
ax1.set_xscale("log")
ax1.set_ylabel("Amplitude")
ax1.axvline(f0, lw = 0.5, color = "k")
ax1.axhline(1 / sqrt(2), lw = 0.5, color = "k")

ax2.set_xscale("log")
ax2.set_ylabel("Phase")
ax2.set_xlabel("f / Hz")
ax2.set_ylim(-100, 100)
ytics = [-90, -45, 0, 45, 90]
ax2.set_yticks(ytics)
for y in ytics
   ax2.axhline(y, lw = 0.5, color = "k")
end
ax2.axvline(f0, lw = 0.5, color = "k")

┌ Warning: Assignment to `h1` in soft scope is ambiguous because a global variable by the same name exists: `h1` will be treated as a new local. Disambiguate by using `local h1` to suppress this warning or `global h1` to assign to the existing global variable.
└ @ ch13.md:725

●◀ 練習：共振回路

複素数 $H(\omega)$ の分母の虚数部分

\[Q \left( \dfrac{\omega}{\omega_0} - \dfrac{\omega_0}{\omega} \right)\]

を，角振動数 $\omega$ の関数として描け．その特徴を記せ．

振幅の二乗が最大値の半分 $\left\vert{H(\omega)}\right\vert^2 = \dfrac{1}{2}$ になる角振動数は2つある．その周波数 $\omega_1, \omega_2$ を数値的に求めてみよ．参考 → ▼ 「はさみうち」法による，方程式の求解

２つの角振動数の差 $\Delta\omega = \omega_2 - \omega_1$ を，角振動数の半値全幅（FWHM; full width of half maximum）という．

定数 $Q$ は，半値全幅 $\Delta{\omega}$ と共振角振動数 $\omega_{0}$ の比に，ほぼ等しい．

\[Q \simeq \dfrac{\omega_0}{\Delta\omega} = \dfrac{f_0}{\Delta{f}}\]

これを，数値的に確かめてみよ．

ヒント：半値全幅を与える角振動数 $\omega_{1}, \omega_2$ では，以下が成り立つ．複号 $\pm$ が，$\omega_{1}, \omega_{2}$ のどちらかに対応するかを考えよ．

\[H(\omega_{1,2} ) = \dfrac{1}{1 \pm i}\]

◀ 低域通過フィルタ

自己インダクタンス L と電気抵抗 R を直列に接続したLR直列回路の電源に正弦波電圧 $e(t)$ を加える．抵抗 $R$ の両端子間の電圧 $v(t)$ は，同じ角振動数を持つ正弦波 $v(t)$ となる．正弦波 $e(t)$ と $v(t)$ を複素数 $E(\omega)$ と $V(\omega)$ で表したとき，両者の比は，以下の複素数 $H(\omega)$ で表される．

\[H(\omega) = \dfrac{V(\omega)}{E(\omega)} = \dfrac{1}{1+i\dfrac{\omega}{\omega_0}}\]

ここで，$\omega_0 = 2\pi{f_0}$ は，回路素子 $L, R$ の値から決まる正の定数である．$f_0$ は周波数の次元，$\omega_{0}$ は角振動数の次元を持つ．

まず，$f_{0} = 1\;\mathrm{kHz}$ として，$H(\omega)$ の振幅 $\left\vert{H(\omega)}\right\vert$ を，描く．横軸の周波数は対数で表示する．

f0 = 1e3
w0 = 2pi * f0

fs = exp10.(range(1, 5, length = 100))

ws = 2pi * fs
h1 = 1 ./ (1 .+ im * ws / w0)

using PyPlot
plt.plot(fs, abs.(h1))
plt.xscale("log")
plt.xlabel("f / Hz")
plt.ylabel("Amplitude")

plt.axhline(1, lw = 0.5, color = "k")
plt.axhline(1 / sqrt(2), lw = 0.5, color = "k")
plt.axvline(f0, lw = 0.5, color = "k")

振幅は，周波数が低いとき $1$ にほぼ等しく，周波数の増加に伴い，単調減少する．周波数が低い正弦波をそのまま通し，周波数の高い正弦波を減衰させるので，「低域通過フィルタ（Low Pass Filter; LPF）」と呼ばれる．

周波数 $f_0$ における振幅は $\dfrac{1}{\sqrt{2}}$ である．周波数 $f_0$ は，振幅が減衰する周波数の目安であり，「遮断（しゃだん）周波数（cut-off frequecy）」と呼ばれる．

振幅 $A$ の常用対数（底が $10$ の対数）をとり，$20$ 倍したものを「利得（gain）」という．利得の「単位」はデシベル（dB）である．

\[G = 20 \log_{10} A\]

上のグラフを，縦軸を利得に変換して描く．

f0 = 1e3
w0 = 2pi * f0

fs = exp10.(range(1, 5, length = 100))

ws = 2pi * fs
h1 = 1 ./ (1 .+ im * ws / w0)

using PyPlot
plt.plot(fs, 20 * log10.(abs.(h1)))
plt.xscale("log")
plt.xlabel("f / Hz")
plt.ylabel("gain / dB")

plt.axhline(0, lw = 0.5, color = "k")
plt.axhline(-3, lw = 0.5, color = "k")
plt.axvline(f0, lw = 0.5, color = "k")

振幅 $1$ は利得 $0\;\mathrm{dB}$ である．減衰動作の周波数領域では，周波数が $10$ 倍になると利得は $-20\;\mathrm{dB}$ 減少する．この傾きを $-20\;\mathrm{dB}/\mathrm{decade}$ と称する（decade は $10$ 倍の意味）．

更に，位相を加えて描く．

f0 = 1e3
w0 = 2pi * f0

fs = exp10.(range(1, 5, length = 100))

ws = 2pi * fs
h1 = 1 ./ (1 .+ im * ws / w0)

using PyPlot
fig, axs = plt.subplots(2, 1)
#
ax1 = axs[1]
ax1.plot(fs, 20 * log10.(abs.(h1)))
ax1.set_xscale("log")
ax1.set_xlabel("f / Hz")
ax1.set_ylabel("gain / dB")

ax1.axvline(f0, lw = 0.5, color = "k")
ax1.axhline(0, lw = 0.5, color = "k")
ax1.axhline(-3, lw = 0.5, color = "k")
#
ax2 = axs[2]
ax2.plot(fs, rad2deg.(angle.(h1)))
ax2.set_ylim(-100, 10)
ax2.set_yticks([-90, -75, -60, -45, -30, -15, 0]) # y軸の刻みを指定する
ax2.set_xscale("log")
ax2.set_xlabel("f / Hz")
ax2.set_ylabel("phase / degree")
ax2.axvline(f0, lw = 0.5, color = "k")
ax2.axhline(-45, lw = 0.5, color = "k")
ax2.axhline(0, lw = 0.5, color = "k")
ax2.axhline(-90, lw = 0.5, color = "k")

位相は周波数の増加に伴い単調減少する．低い周波数では $0^{\circ}$ に，高い周波数では $-90^{\circ}$ にそれぞれ漸近する．

遮断周波数 $f_0$ では伝達関数は $H(\omega_0) = \dfrac{1}{1+i}$ である．したがって，遮断周波数 $f_0$ での利得は $-3\;\mathrm{dB}$ ,位相は $-45^{\circ}$ である（註: $\log_{10} 2 \simeq 0.3$ ）．

以下の伝達関数に対して，周波数応答を描け．

\[H(\omega) = \dfrac{-A}{1+i\dfrac{\omega}{\omega_0}}\]

まず，$A = 1$ を保ったまま $f_0 = 100, 1000, 10000\;\mathrm{Hz}$ と増やしてみよ．

次に，$A = 1, 10, 100$ と増やしてみよ．

それぞれ，どのように変化するか，言葉で記述してみよ．

■ 複素数に拡張された関数

平方根や指数関数以外でも，実数を引数とする関数の多くが，複素数を引数とするように拡張されている．

▶ 対数関数（複素数）

複素数 $z = r\exp(i\theta)$ と極座標表示したとき，その自然対数は，以下のように計算できる．すなわち，実数部は絶対値の自然対数，虚数部は偏角である．

\[\log{z} = \log r\exp(i\theta) = \log{r} + i \theta\]

複素数 $z = 1 + i y$ の自然対数 $\log{z}$ の実数部と虚数部を描く．

using PyPlot
plt.axes().set_aspect("equal")

ys = range(-2pi, 2pi, length = 101)
zs = complex.(1, ys)
cs = log.(zs)
plt.plot(ys, real.(cs), label = "real log (1+i*y)")
plt.plot(ys, imag.(cs), label = "imag log (1+i*y)")
plt.legend()

plt.xlabel("y")
plt.legend(loc = 4)
plt.ylim(-3, 3)
plt.xlim(-3, 3)

plt.axhline(0, lw = 0.5, color = "k")
plt.axvline(0, lw = 0.5, color = "k")
plt.axhline(pi / 2, lw = 0.5, color = "c")
plt.axhline(-pi / 2, lw = 0.5, color = "c")
#
plt.axhline(log(sqrt(2)), lw = 0.5, color = "m")
plt.axhline(pi / 4, 0.6, 1, lw = 0.5, color = "m")
plt.axhline(-pi / 4, 0, 0.4, lw = 0.5, color = "m")
plt.axvline(1, lw = 0.5, color = "m")
plt.axvline(-1, lw = 0.5, color = "m")

変数 $y$ の増加に伴い，$\log{z}$ の虚数部は $-\dfrac{\pi}{2}$ から $\dfrac{\pi}{2}$ へ単調に増加する（シアン色の補助線）．実数部は下に凸で，$y=0$ で極小値 $\log{1} = 0$ をとる．

特に，$y = \pm{1}$ において，虚数部は $\pm\dfrac{\pi}{4}$ ,実数部は $\log\sqrt{2}$ をとる（マゼンダ色の補助線）．

上のグラフは，ガウス平面上の $z = 1 + i y$ の軌跡から理解できるであろう．

using PyPlot
plt.axes().set_aspect("equal")

plt.xlim(-3, 3)
plt.ylim(-3, 3)
plt.xlabel("real part")
plt.ylabel("imag part")
plt.axhline(0, lw = 0.5, c = "k")
plt.axvline(0, lw = 0.5, c = "k")
plt.axvline(1, lw = 0.5, c = "k")
for y in [-2, -1, 1, 2]
   plt.plot([0, 1], [0, y], "b.-")
   text(1.1, y, "1+i " * string(y))
end

◀ 練習：対数関数（複素数）

複素数 $z = x + i$ の自然対数 $\log{z}$ の実数部と虚数部を描け．その結果を，ガウス平面上の $z$ の軌跡を描いて，考察せよ．

▶ 三角関数と双曲線関数

双曲線関数 $\cosh{t}, \sinh{t}$ や三角関数 $\cos{t}, \sin{t}$ は，指数関数 $\exp(z)$ を用いて，定義されることもある．

\[\begin{aligned} \cosh x &= \dfrac{\exp(x)+\exp(-x)}{2}, \\ \sinh x &= \dfrac{\exp(x)-\exp(-x)}{2}, \\ \cos x &= \dfrac{\exp(ix)+\exp(-ix)}{2}, \\ \sin x &= \dfrac{\exp(ix)-\exp(-ix)}{2i} \end{aligned}\]

したがって，三角関数に純虚数を与えると，双曲線関数となる．

\[\begin{aligned} \cos ix & = \cosh x, \\ \sin ix & = i \sinh x \end{aligned}\]

上の等式がなりたつことを，グラフで観察しよう．

using PyPlot
plt.axes().set_aspect("equal")

xs = range(-1, 1, length = 21)
cz = cos.(im * xs)
plt.plot(xs, real.(cz), "r-", label = "real cos(ix)")
plt.plot(xs, cosh.(xs), "ro", label = "cosh(x)")
#
sz = sin.(im * xs)
plt.plot(xs, imag.(sz), "b-", label = "imag sin(ix)")
plt.plot(xs, sinh.(xs), "bo", label = "sinh(x)")
#
plt.xlabel("x")
plt.xlim(-1.8, 1.8)
plt.ylim(-1.8, 1.8)
plt.legend()
plt.axhline(0, lw = 0.5, color = "k")
plt.axvline(0, lw = 0.5, color = "k")

★ 今回のまとめ

複素数
複素数のベクトル・行列
複素数の加減乗除
ガウス平面
オイラーの公式
極座標表示
複素数に拡張された関数
応用：伝達関数