2014年3月2日

黒木先生、scilab でベクトル場をプロット

scilab 使おうと思った。
15
黒木玄 Gen Kuroki @genkuroki

複素有理型函数 f(z) を実数値函数 u, v を使って f(z) = u - iv と表示しておいて、ベクトル場 (u,v) をプロットすれば、二次元理想流体ごっこができる。おそらく探せば画像をすでにたくさん公開している人がいるんじゃないかな?

2014-02-28 23:32:24
黒木玄 Gen Kuroki @genkuroki

複素函数 f(z) を実数値函数u, vでf(z) =u-ivと表示しておいて、ベクトル場(u,v)を複素平面にプロット。ベクトルの長さは色で表示。流れは面積保存で正則な点では渦なし。 添付の画像は f(z)=z http://t.co/ZOu03SBYy2

2014-03-01 10:49:01
拡大
黒木玄 Gen Kuroki @genkuroki

続き。電卓代わりにいつも使っているのがscilabです。添付の図はscilabのchamp1函数によるプロット。 添付の画像は f(z)=z^2 です。 http://t.co/H0QIGRfCbc

2014-03-01 10:50:38
拡大
黒木玄 Gen Kuroki @genkuroki

続き。champ1はベクトルの長さを色で表示、長さの通りにベクトルを表示させたい場合はchampを使う。 添付の画像は f(z)=z^3 です。体積(面積)保存で渦なしの流れ。 http://t.co/NlSBXoH8hx

2014-03-01 10:52:37
拡大
黒木玄 Gen Kuroki @genkuroki

続き。xとyを-1から1まで0.1刻みで動かす。z=x+iyとおいて、f(z)の(実部)と-(虚部)を成分とするベクトルをプロット。ただしベクトルの長さは色で区別。 f(z)=1/z です。湧き出しまたは正電荷に見える。 http://t.co/t7ZQi7DA5X

2014-03-01 10:56:03
拡大
黒木玄 Gen Kuroki @genkuroki

続き。これは f(z)= -1/z です。一つ前の奴の-1倍。これは吸い込み口または負電荷粒子に見える。極でる原点以外では体積保存で渦無しの流れ。原点の極に流体または電気力線が吸い込まれている。 http://t.co/hlfY3ZgW9J

2014-03-01 10:59:13
拡大
黒木玄 Gen Kuroki @genkuroki

続き。ベクトルの長さを色で表示しているところがちょっと分かり難い。極(特異点)以外では体積保存で渦無しの流れ。 f(z)=1/z^2 です。これは正電荷粒子と負電荷粒子がくっついている様子になる。 http://t.co/ADpkvHT8i1

2014-03-01 11:01:24
拡大
黒木玄 Gen Kuroki @genkuroki

続き。f(z)=1/(z-0.2) - 1/(z+0.2) は正電荷粒子がz=0.2に負電荷粒子がz=-0.2にある様子。1/z^2の様子に近い。 http://t.co/uYylpUO9uj

2014-03-01 11:03:50
拡大
黒木玄 Gen Kuroki @genkuroki

続き。回転する渦を見たければ f(z)=i/z の様子をプロットしてみればよい。ベクトルの長さを色で表示しているのでわかりにくいが、原点では無限大の速さで回転している。これも原点以外では体積(面積)保存で渦無しの流れになっている。 http://t.co/BqNA8qP9JT

2014-03-01 11:06:00
拡大
黒木玄 Gen Kuroki @genkuroki

続き。+2の電荷を持つ粒子が一つと-1の電荷を持つ粒子が二つ近くにある様子。 これは f(z) = 2/(z-.4)-1/(z+.1+.6i)-1/(z+.3-.4i) です。 http://t.co/diBHcAmH8G

2014-03-01 11:10:29
拡大
黒木玄 Gen Kuroki @genkuroki

続き。複素函数 f(z) を実数値函数u, vでf(z) =u-ivと表示し、ベクトル場(u,v)を複素平面にプロットした。描かれる流れの図が体積(面積)保存かつ渦無しになることとf(z)が正則函数になることは数学的に同値。正則函数は体積保存で渦無しの流れと同じ。

2014-03-01 11:12:39
黒木玄 Gen Kuroki @genkuroki

続き。以上では簡単な複素函数の流れのみをプロットしてみた。もっと複雑な有理型函数であっても、零点と極の情報さえ得られれば、その近所の流れを全体で適当に繋げれば全体の流れの図が得られる。ただし、体積保存と渦無しという性質に注意する。

2014-03-01 11:14:56
黒木玄 Gen Kuroki @genkuroki

続き。よくある質問【f(z)=u-ivとしていますが、どうしてすなおにf(z) =u+ivとしないのですか?】への回答。虚部を-1倍しているのはコーシー・リーマンの方程式と体積(面積)保存で渦無しという条件が同値になるようにするためです。実際に流れをプロットすればわかる。

2014-03-01 11:18:33
黒木玄 Gen Kuroki @genkuroki

続き。scilabの使い方。(1)scilabを使えるようにする。無料。ググれ。(2)scilabに clear; x=[-1:.1:1]; y=[-1:.1:1]; [xr xc]=size(x); [yr yc]=size(y); と入力。xとyは-1から1まで.1刻み。

2014-03-01 11:23:38
黒木玄 Gen Kuroki @genkuroki

続き for i=1:xc, for j=1:yc, z=complex(x(i),y(j)); if z==0 then v=0; else v=1/z; end;, vx(i,j)=real(v);, vy(i,j)=-imag(v); end;, end; と入力

2014-03-01 11:25:09
黒木玄 Gen Kuroki @genkuroki

続き clf(); champ1(x,y,vx,vy) でベクトルの長さが色で表示されたベクトル場が表示される。色を使わない場合にはchampを使う。たったこれだけです。

2014-03-01 11:26:17
黒木玄 Gen Kuroki @genkuroki

続き。scilabに限らず、多くの人が使っているソフトの使い方はググればわかる。

2014-03-01 11:27:40
黒木玄 Gen Kuroki @genkuroki

以上のscilabで複素函数のベクトル場をプロットする話のまとめ読みは http://t.co/DFMTexNGbb でできます。添付の画像は f(z)=z です。 http://t.co/ZOu03SBYy2

2014-03-01 11:33:12
拡大
黒木玄 Gen Kuroki @genkuroki

scilabでの等高線の描き方 https://t.co/6uFRWtI9oJ たとえば先のツイートの入力の続きで contour2d(x,y,vy,21); とかすればいいのか。最後の21は等高線の本数。

2014-03-01 11:56:21
黒木玄 Gen Kuroki @genkuroki

ツイッターと同じパソコンの裏でやっていたこと。それは複素函数のベクトル場とその不定積分の虚部の等高線を描く遊び。https://t.co/7VA7or8u4f の続き。グラフをどばっと放流します。河合・竹井『特異摂動の代数解析学』p.27を開い待つといいかも。続く

2014-03-02 20:54:06
黒木玄 Gen Kuroki @genkuroki

ツイッターと同じパソコンの裏でやっていたこと。それは複素函数のベクトル場とその不定積分の虚部の等高線を描く遊び。https://t.co/7VA7or8u4f の続き。グラフをどばっと放流します。河合・竹井『特異摂動の代数解析学』p.27を開い待つといいかも。続く

2014-03-02 20:54:06
黒木玄 Gen Kuroki @genkuroki

以下、√(Q(z)) の((実部),-(虚部))をプロットしたベクトル場と∫√(Q(z))dz の虚部の等高線をプロットしたものを放流。平方根の分岐がもろに見えちゃうので注意。 Q(z)=z の場合のベクトル場。 http://t.co/O7mDWSx2VD

2014-03-02 20:58:48
拡大
黒木玄 Gen Kuroki @genkuroki

以下、√(Q(z)) の((実部),-(虚部))をプロットしたベクトル場と∫√(Q(z))dz の虚部の等高線をプロットしたものを放流。平方根の分岐がもろに見えちゃうので注意。 Q(z)=z の場合のベクトル場。 http://t.co/O7mDWSx2VD

2014-03-02 20:58:48
拡大
黒木玄 Gen Kuroki @genkuroki

続き。Q(z)=z の場合の∫√(Q(z))dz の虚部の等高線。 http://t.co/sAEinzxgDO

2014-03-02 21:00:09
拡大
残りを読む(14)

コメント

あびこ🐄 @cowabiko 2014年3月5日
まとめを更新しました。 scalib じゃなくて scilib(サイリブ)だった。恥ずかしい。
0
あびこ🐄 @cowabiko 2014年3月6日
あ、scilab(サイラボ)か。恥ずかしい。
0
黒木玄 Gen Kuroki @genkuroki 2014年3月6日
まず、ごめんなさい。https://twitter.com/genkuroki/status/440098596310896640 2014-03-02 21:17:36の図は間違ってました。 https://twitter.com/genkuroki/status/440276609657212929 に訂正版を投稿しました。
0
黒木玄 Gen Kuroki @genkuroki 2014年3月6日
そして、まとめてくださって、どうもありがとうございます。 https://twitter.com/genkuroki/status/440322068346064896 以降に数学的解説があります。内容は正則函数のベクトル場がその不定積分の等高線に接するということ。
0
あびこ🐄 @cowabiko 2014年3月6日
scilab(スキャリブでもサイリブでもなくサイラボね)、ようはMATLAB(マトラボ)チックなことができるということなんでしょうか(MATLABも触ったことないけど)。熱方程式とか、拡散方程式の勉強に使えるかなーと思って触りはじめたんですが、どうなんでしょう。
0
あびこ🐄 @cowabiko 2014年3月6日
うお、黒木先生ご本人がコメントしてくださってますね。ありがとうございます。
0
黒木玄 Gen Kuroki @genkuroki 2014年3月6日
https://twitter.com/genkuroki/status/441447953299808256 からの連ツイで続きを書きました。 動機は二つ。forループを使うのはできる限り止めたい。0除算関係の処理をもっとサボりたい。
0
黒木玄 Gen Kuroki @genkuroki 2014年3月6日
https://twitter.com/genkuroki/status/441529167205924864 からの連ツイでオープン戸田格子(の非線形常微分方程式)と離散オープン戸田格子をscilabで解く方法(結構簡単)について連ツイしました。
0
黒木玄 Gen Kuroki @genkuroki 2014年3月9日
https://twitter.com/genkuroki/status/442318369643954176 からKdV方程式の数値計算について連ツイしました。scilabを使えば無料でコードはたったの十数行!
0
黒木玄 Gen Kuroki @genkuroki 2014年3月9日
https://twitter.com/genkuroki/status/442321386346786816 から無料で誰でも使えるscilabによる十数行のコードによるKdV方程式の数値計算について追記。
0