LoginSignup
4
3

More than 1 year has passed since last update.

@Helmet(へるめっと)

python使って複素積分を可視化・イメージ化して理解してみるチャレンジ

Last updated at Posted at 2022-07-19

概要

複素積分をpython使ってベクトルで表してみました。
イメージの一つとして誰かの参考になれば幸いです。

matplotlibの練習にもなります。クリック検出とか、リアルタイム更新とかします。

イメージ

1. 複素関数をベクトル場で表現してみる

複素積分を可視化する前に、積分する対象である複素関数を可視化できないとキツイので、やってみます。
今回は以下の手順でベクトル場にしてみました。(f(z)=1/zで考えてます。)

  1. z平面、w平面を用意　（横軸=実数軸、縦軸=虚数軸です。）
  2. z平面のベクトルの先端にw平面のベクトルをくっつける
  3. 上記を繰り返してベクトルをたくさん作る。
    image.png
    こんな感じでpythonでもベクトル場を作っていきます。
    ただ、ベクトルの大きさを表すことが出来ないので、色で表すことにしています。
    image.png

2. 複素積分をベクトル表現にしてみる

経路をC:|z|=1の単一閉曲線に沿う複素関数f(z)=1/z の複素積分

\int_C\frac{1}{z}dz = 2\pi i

をやってみると下のような感じになります。
画像1.png
それぞれの矢印は以下の項目を表してます。
fukuso_2.png
左は　　　ベクトル場と経路Ｃを表すグラフで、
真ん中は　その時のf(z),dzベクトルとその積のグラフ
右は　　　これまで計算してきた積のベクトルとその和のグラフ

結果は、
sum f(z)*dz = -1.22 + 6.13 i
なので、まあ2πiと言えなくもないかなと。
もうちょっと円の精度良くすれば値は近づきます。

コーシーの積分定理はこんな感じ。
f(z)が単一閉曲線C及びその内部で正則なとき、

\int_Cf(z)dz = 0

でやってみると、
プレゼンテーション1.gif

って感じで、割とゼロに近づいてく様子が分かるんじゃないかと思います。
値はsum f(z)*dz = 0.08 + 0.02 iでした。

注意：

製作者は複素積分についてほぼ理解していません。
学術的に正しいかは分かりません。
てかグリーンの定理から既にわかりません。教えてえrい人。

環境

  • windows10
  • python 3.7.5

全体スクリプト

描画スクリプトの全容です。サイト1)2)を参考にしました。
コピペして実行して、グラフ上でマウスのホイールを押したら経路Cを描けると思います。

complex_plot.py
import numpy as np
import matplotlib.pyplot as plt

LX, LY=2,2 #描画範囲
gridwidth=0.01 #矢印を描く密度

# 画面準備(1*3個)
fig, ax = plt.subplots(1, 3)
fig.set_size_inches(20, 7) #画像サイズ

# タイトルを付ける
ax[0].set_title('C & f(z)')
ax[1].set_title('f(z) & dz & f(z)*dz')
ax[2].set_title('f(z)*dz & ∫f(z)dz')
# 格子を加える
ax[0].grid()
ax[1].grid()
ax[2].grid()
# 範囲指定
ax[0].set_xlim([-LX,LX])
ax[0].set_ylim([-LY,LY])
ax[1].set_xlim([-1,1])
ax[1].set_ylim([-1,1])
ax[2].set_xlim([-1,1])
ax[2].set_ylim([-1,1])


# 被積分複素関数の定義
def f(Z):
    # return Z
    # return np.sin(Z)
    return 1/Z
    # return np.sin(Z)/(Z**2+1)
#零点(特異点)は手打ちする
Xzero, Yzero = 0, 0

"ベクトル場作成"
X, Y= np.meshgrid(np.arange(-LX, LX, gridwidth), np.arange(-LY, LY,gridwidth))
Z = (X + Y * 1j)
W = f(Z)
U = W.real
V = W.imag

#正規化したかったら以下もつける
R = np.sqrt(U**2 + V**2)
U = U / R
V = V / R

ax[0].plot(Xzero,Yzero,'o',color='black')
strm = ax[0].streamplot(X, Y, U, V, color=R, linewidth=2, cmap='summer')
plt.colorbar(strm.lines)

"クリックイベント"
C_start = []
C_end = []
dz = []
save_fz_dz = []
sum_fz_dz = 0+0j
def draw_C_vec(event):
    global C_start
    global C_end
    global dz
    global save_fz_dz
    global sum_fz_dz
    
    if event.button == 2: # ホイールクリックなら
        # 始点取得
        if C_start == []: 
            C_start = [event.xdata, event.ydata]
            ax[0].plot(C_start[0], C_start[1],
                       "o", color="blue")
            # 更新
            fig.canvas.draw()
            print("get_start")

        # 終点取得＋処理
        else:
            C_end = [event.xdata, event.ydata]
            # dzベクトル取得
            dz = (C_end[0]-C_start[0]) + (C_end[1]-C_start[1]) * 1j

            # 要素の計算 (z_i, f(z), f(z)*dz)
            z = C_start[0]+C_start[1]*1j
            f_z = f(z)
            fz_dz = f_z * dz

            " 以下グラフ描画 "
            # 経路C(dz_i)をプロット
            ax[0].arrow(C_start[0],C_start[1],
                        dz.real, dz.imag,
                        width=0.01,color="blue")
            # クリック点プロット
            ax[0].plot(C_end[0], C_end[1],
                        "o", color="blue")
            # f(z)ベクトル
            ax[0].arrow(C_start[0],C_start[1],
                        f_z.real,  f_z.imag,
                        width=0.01,color="cyan")

            # 描画し直し
            ax[1].clear()
            ax[2].clear()
            # タイトルを付ける
            ax[1].set_title('f(z) & dz & f(z)*dz')
            ax[2].set_title('f(z)*dz & ∫f(z)dz')
            # 格子を加える
            ax[1].grid()
            ax[2].grid()
            # 範囲指定
            ax[1].set_xlim([-1,1])
            ax[1].set_ylim([-1,1])
            ax[2].set_xlim([-1,1])
            ax[2].set_ylim([-1,1])

            # f(z)ベクトル
            ax[1].arrow(0,0,
                        f_z.real, f_z.imag,
                        width=0.01,color="cyan")
            # dzベクトル
            ax[1].arrow(0,0,
                        dz.real, dz.imag,
                        width=0.01,color="blue")
            # f(z)*dzベクトル
            ax[1].arrow(0,0,
                        fz_dz.real, fz_dz.imag,
                        width=0.01,color = "red")
            
            # 積分値グラフ
            sum_fz_dz += fz_dz
            save_fz_dz.append(fz_dz)

            # これまでのf(z)*dzを描画
            for one_fz_dz in save_fz_dz:
                ax[2].arrow(0,0,
                            one_fz_dz.real, one_fz_dz.imag,
                            width=0.01,color = "red")
            # 積分値を描画
            ax[2].arrow(0,0,
                        sum_fz_dz.real, sum_fz_dz.imag,
                        width=0.01,color = "gold")
            # 3桁まで表示
            print("sum f(z)*dz = {z.real:.2f} + {z.imag:.2f} i".format(z=sum_fz_dz))
            
            # グラフ更新
            fig.canvas.draw()
            # 値の更新
            C_start = C_end.copy()
            
# クリックに反応するよう設定
cid = fig.canvas.mpl_connect('button_press_event', draw_C_vec)

# 描写開始
plt.draw()
plt.show()

説明

データ・グラフの用意

対応スクリプト（クリックすると開きます）
グラフ準備
import numpy as np
import matplotlib.pyplot as plt

LX, LY=2,2 #描画範囲
gridwidth=0.01 #矢印を描く密度

# 画面準備(1*3個)
fig, ax = plt.subplots(1, 3)
fig.set_size_inches(20, 7) #画像サイズ

# タイトルを付ける
ax[0].set_title('C & f(z)')
ax[1].set_title('f(z) & dz & f(z)*dz')
ax[2].set_title('f(z)*dz & ∫f(z)dz')
# 格子を加える
ax[0].grid()
ax[1].grid()
ax[2].grid()
# 範囲指定
ax[0].set_xlim([-LX,LX])
ax[0].set_ylim([-LY,LY])
ax[1].set_xlim([-1,1])
ax[1].set_ylim([-1,1])
ax[2].set_xlim([-1,1])
ax[2].set_ylim([-1,1])

色々グラフの準備をします。

被積分関数の用意

被積分関数
# 被積分複素関数の定義
def f(Z):
    # return Z
    # return np.sin(Z)
    return 1/Z
    # return np.sin(Z)/(Z**2+1)
#零点(特異点)は手打ちする
Xzero, Yzero = 0, 0

関数に入れるだけにしました。
ただ、ゼロ点だけは手打ちのままです。

ベクトル場の作成

対応スクリプト（クリックすると開きます）
ベクトル場作成
"ベクトル場作成"
X, Y= np.meshgrid(np.arange(-LX, LX, gridwidth), np.arange(-LY, LY,gridwidth))
Z = (X + Y * 1j)
W = f(Z)
U = W.real
V = W.imag

#正規化したかったら以下もつける
R = np.sqrt(U**2 + V**2)
U = U / R
V = V / R

ax[0].plot(Xzero,Yzero,'o',color='black')
strm = ax[0].streamplot(X, Y, U, V, color=R, linewidth=2, cmap='summer')
plt.colorbar(strm.lines)

ここでは、まず格子上の各座標をX,Yで作って、Zに複素数に変換して格納してます。
そのZを先ほどの関数で変換してWにします。
そのあとUとかVとか、streamplot()関数の引数を用意していく感じです。

クリック反応

対応スクリプト（クリックすると開きます）
クリック反応
"クリックイベント"
C_start = []
C_end = []
dz = []
save_fz_dz = []
sum_fz_dz = 0+0j
def draw_C_vec(event):
    global C_start
    global C_end
    global dz
    global save_fz_dz
    global sum_fz_dz
    
    if event.button == 2: # ホイールクリックなら
        # 始点取得
        if C_start == []: 
            C_start = [event.xdata, event.ydata]
            ax[0].plot(C_start[0], C_start[1],
                       "o", color="blue")
            # 更新
            fig.canvas.draw()
            print("get_start")

        # 終点取得＋処理
        else:
            C_end = [event.xdata, event.ydata]
            # dzベクトル取得
            dz = (C_end[0]-C_start[0]) + (C_end[1]-C_start[1]) * 1j

            # 要素の計算 (z_i, f(z), f(z)*dz)
            z = C_start[0]+C_start[1]*1j
            f_z = f(z)
            fz_dz = f_z * dz

            " 以下グラフ描画 "
            # 経路C(dz_i)をプロット
            ax[0].arrow(C_start[0],C_start[1],
                        dz.real, dz.imag,
                        width=0.01,color="blue")
            # クリック点プロット
            ax[0].plot(C_end[0], C_end[1],
                        "o", color="blue")
            # f(z)ベクトル
            ax[0].arrow(C_start[0],C_start[1],
                        f_z.real,  f_z.imag,
                        width=0.01,color="cyan")

            # 描画し直し
            ax[1].clear()
            ax[2].clear()
            # タイトルを付ける
            ax[1].set_title('f(z) & dz & f(z)*dz')
            ax[2].set_title('f(z)*dz & ∫f(z)dz')
            # 格子を加える
            ax[1].grid()
            ax[2].grid()
            # 範囲指定
            ax[1].set_xlim([-1,1])
            ax[1].set_ylim([-1,1])
            ax[2].set_xlim([-1,1])
            ax[2].set_ylim([-1,1])

            # f(z)ベクトル
            ax[1].arrow(0,0,
                        f_z.real, f_z.imag,
                        width=0.01,color="cyan")
            # dzベクトル
            ax[1].arrow(0,0,
                        dz.real, dz.imag,
                        width=0.01,color="blue")
            # f(z)*dzベクトル
            ax[1].arrow(0,0,
                        fz_dz.real, fz_dz.imag,
                        width=0.01,color = "red")
            
            # 積分値グラフ
            sum_fz_dz += fz_dz
            save_fz_dz.append(fz_dz)

            # これまでのf(z)*dzを描画
            for one_fz_dz in save_fz_dz:
                ax[2].arrow(0,0,
                            one_fz_dz.real, one_fz_dz.imag,
                            width=0.01,color = "red")
            # 積分値を描画
            ax[2].arrow(0,0,
                        sum_fz_dz.real, sum_fz_dz.imag,
                        width=0.01,color = "gold")
            # 3桁まで表示
            print("sum f(z)*dz = {z.real:.2f} + {z.imag:.2f} i".format(z=sum_fz_dz))
            
            # グラフ更新
            fig.canvas.draw()
            # 値の更新
            C_start = C_end.copy()
            
# クリックに反応するよう設定
cid = fig.canvas.mpl_connect('button_press_event', draw_C_vec)

最後の行を用意すれば、クリック時に関数が呼べるので、その関数に機能を入れていきます。
ベクトルの値を保存したり、足し算したりと色々してるのでごちゃごちゃになってます。お許しください！博士！

ホイールをクリックすることで、event.xdata,event.ydataでカーソル位置から始点or終点を取得して、そこからベクトルを書いていきます。

ベクトル書く部分は、

  1. 左グラフに　経路、点、f(z)　を描く
  2. 真ん中グラフに　f(z)、dz、f(z)*dzのベクトル　を描く
  3. 右グラフに　積分結果と、これまでのf(z)*dzベクトル　を描く
    という順番になってます。

クリックについては、左右のクリックは拡大縮小移動などに使うので、今回はホイールをクリックするようにしました。
左右クリックにしたい場合は、if event.button == 2のとこを1とか3にすると使えます。

図から分かる複素積分の用語・性質の考察

例

f(z)=sin(Z)
詳細（クリックすると開きます） pythonでは return np.sin(Z) のように定義

範囲=[-2,2]
image.png

image.png
範囲=[-6,6]
image.png
(0,0)の点から放射状に出てそうです。
波っぽいのが出てきました。
範囲を広げると、周期関数であることが分かりやすいですね。

f(z)=\frac{1}{Z^n}
詳細（クリックすると開きます） pythonでは return 1/Z**2 のように定義

範囲=[-2,2]
左からn=1,2,3
image.png
これを見ると、-n乗のZという関数は、

  • 特異点から、2*n本の出ていく部分と入ってくる部分がある。
  • 単一閉曲線で積分しようとすると、水色のf(z)ベクトルがn周する
    などが言えそうです。
    文献1)によると、位数は回転数を表しているのだそう。
    となると、
    位数m = 分母の(z-a)をm乗してる = 積分時のf(z)がm周回転
    ってコト？！
f(z)=\frac{e^{iZ}}{Z^2+1}
詳細（クリックすると開きます） pythonでは return 1/Z**2 のように定義

僕の持ってる教科書に出てきたやつです。
範囲=[-2,2]
image.png
ここから教科書ではcosx/(x^2+1)の実積分を導いてるんですが、うーん。わからぬ。

正則かどうか

定義としては、「コーシー・リーマンの関係式が成り立つか？」です。

\frac{\partial u}{\partial x}=\frac{\partial v}{\partial y}\\
\frac{\partial v}{\partial x}=-\frac{\partial u}{\partial y}

グラフから直感的なことを言うなら、
f(z)=∞となる場所があるかないか（矢印の出発点・終着点があるかないか）

  • 出発点 = f(z) → +∞
  • 終着点 = f(z) → -∞
    という事になりそうです。

数学的・物理的背景

この動画を参考にしています。（https://www.youtube.com/watch?v=EyBDtUtyshk&t=2164s）
グラフを見てお分かりかもしれませんが、電磁気学の電場や磁場、流体力学の流れ場などの流束に関係してきます。
この動画では、コーシー・リーマンの関係式とは、複素共役ベクトル(u-iv)を考えた時、
発散が無い

\nabla \cdot w^* = 0\\
\nabla \cdot\left(\begin{matrix} 
u \\ 
-v
\end{matrix} 
\right)
= \frac{\partial u}{\partial x} - \frac{\partial v}{\partial y} = 0\\
∴\frac{\partial u}{\partial x}=\frac{\partial v}{\partial y}

回転が無い

\nabla \times w^* = 0\\
\nabla \times \left(\begin{matrix} 
u \\ 
-v
\end{matrix} 
\right)
= \frac{\partial v}{\partial x} - ( -\frac{\partial u}{\partial y}) = 0\\
∴\frac{\partial v}{\partial x}=-\frac{\partial u}{\partial y}

という事をあらわあしています。
また、動画では、複素共役を扱っていますが、複素共役を用いると、

  • 上記の関係式が成り立つ=正則となる
  • 流束として考えることが出来るようになる
    というメリットがあるため、複素共役を取っているんだと思います。

正則ならば積分は経路によらない

  • 正則ならば、f(z)ベクトル場の矢印の方向が大体一緒なので、経路がちょっと変わっても、結果は一緒となる？
  • f(z)ベクトルが変わっても、それと対応して経路C(dz)の向きも変わるから一緒になる。
  • 非正則な点を横切ると値が∞に吹っ飛ぶので一緒も何もなくなる。

留数・留数定理

グラフ上では何を意味してるのかわかりませんでした。誰か教えてください。

孤立特異点の種類（除去可能な特異点・位数nの極・真性特異点）

除去可能な特異点はグラフ上で特異点っぽくならない
位数については下のやつ。
真性特異点はグラフ上でも特異点

位数

水色のf(z)ベクトルを見ていくと、水色ベクトルの回転数と位数が一致してそう
参考はこちら1)

動画では、ベクトルが出てくる点をsource, 入ってくる点をsinkと呼び、それを合成する数が位数であるとも考えられる…？

まとめ

複素積分の視覚化にチャレンジしました。
ここまでやっても複素関数論分からないんですがどうしたらいいでしょうか…
ご閲覧ありがとうございました。ご意見・ご指摘等、よろしくお願いいたします。

参考サイト

1)複素関数を見てみよう
https://monaqua.netlify.app/Speak_20210613_sundaymath.pdf
2)Matplotlibの使い方｜Pythonによる可視化入門 #3
https://lib-arts.hatenablog.com/entry/Python_visualize3
3)Complex integration, Cauchy and residue theorems | Essence of Complex Analysis #6
https://www.youtube.com/watch?v=EyBDtUtyshk&t=2164s
めちゃわかりやすいです。サイト製作者様にこの場をお借りして感謝申し上げます。

4
3
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
4
3