黒木玄 Gen Kuroki
- いいね数 389,756/311,170
- フォロー 995 フォロワー 14,556 ツイート 293,980
- 現在地 (^-^)/
- Web https://genkuroki.github.io/documents/
- 自己紹介 私については https://twilog.org/genkuroki と https://genkuroki.github.io と https://github.com/genkuroki と https://github.com/genkuroki/public を見て下さい。
2018年02月03日(土)
twitter.com/wkbme/status/9...
#Julia言語 ではmatplotlib.animationをほぼそのまま使えます。私のGIF動画の大部分はJuliaからmatplot.animationを使って作ったものです。
タグ: Julia言語
posted at 00:10:02
#Julia言語 でPlotsパッケージ+GRバックエンドを使えばmatplotlib.animationよりもずっと高速にGIF動画を作れます。
nbviewer.jupyter.org/gist/genkuroki...
の1.7.5節を見て下さい。
matplotlib.animationで70秒以上かかった動画作成が10秒で終わっています!
タグ: Julia言語
posted at 00:16:18
twitter.com/ceptree/status...
Plotsが速いのではなく、GRバックエンドが速い。
Python側でもGRバックエンドを使えれば速いはず。
gr-framework.org/python.html#
タグ:
posted at 00:30:09
「ニセ科学」と呼ばれているものはどれもあからさまに非科学なので、科学と非科学の境界問題とかそういう微妙な話を持ち出す必要はありません。そういうのは科学哲学の仕事。ニセ科学はもっとずっと簡単な話。むしろ、「なぜそんなあからさまな非科学が広まるのか」という問題だと思ってもらえば
タグ:
posted at 00:44:58
matplotlibでもgr使えそうだけど、通常のmatplobliのバックエンドを変えるやり方とは違ってめんどそう。
gr-framework.org/tutorials/matp...
gr-framework.org/tutorials/mpl_...
タグ:
posted at 00:47:04
でもアニメーションやインタラクティブプロットにいいらしい
Starting with version 0.6.0 of the GR framework Matplotlib’s capabilities can be greatly expanded, especially for dynamic or interactive plots.
タグ:
posted at 00:48:09
ちなみにpyenv使ってるからなのか、MacOSXだからなのか、このやり方でできなかった。
gr-framework.org/tutorials/matp...
タグ:
posted at 00:52:40
@genkuroki おそらくPyPlot_plot_Assetの中でsubplotを使って毎回Axesオブジェクトを新規作成してるのがボトルネックです。matplotlib 2.2で修正されるかもしれませんがAxes作成には数100ms無駄に時間がかかります。figと同時にaxを作ってフレーム毎の内容の更新もclfを使わなければだいぶマシになるかと。
タグ:
posted at 02:04:40
™ (blueskyに同アカウント名で避 @tmaehara
数値的存在証明です。ある方程式が解をもつかの判定を数値計算を使いつつ「数学的に厳密に」行なえます。 twitter.com/jaialkdanel/st...
タグ:
posted at 05:41:38
™ (blueskyに同アカウント名で避 @tmaehara
精度保証付き数値計算の最も有名な利用例は「ローレンツ方程式はストレンジアトラクタを解にもつか」というSmalの14th問題のTuckerによる解決 mathworld.wolfram.com/news/2002-02-1... だと思います。「21世紀の数学の難問」が「数値計算を使う計算機支援証明」で解決されるのは良いストーリーですね。
タグ:
posted at 05:52:05
™ (blueskyに同アカウント名で避 @tmaehara
精度保証のほうが計算速度が早くなるタイプの問題もあります。計算幾何の問題で適当にやると幾何的に破綻してしまう問題について、ふつうは厳密計算法と呼ばれる任意桁演算(結構重い)をするところを精度保証(任意桁よりは軽い)で済ますテクニックがあります www.eng.nus.edu.sg/civil/REC2010/... twitter.com/jaialkdanel/st...
タグ:
posted at 06:05:47
Dr. Chris Rackauckas @ChrisRackauckas
I accepted an applied mathematics instructorship @MIT. I will continue to work on methods for efficiently solving differential equations and estimating dynamical models from biological/pharmaceutical data in #julialang. Thank you #julialang community for your support.
タグ: julialang
posted at 06:22:05
非公開
タグ:
posted at xx:xx:xx
@genkuroki 気になって色々試してみたんですが意外とsubplotは大したことなくてhistに200msかかってます。確かに128本の長方形をAxesに追加するための細々した設定を含めて毎回やってれば遅くなりそうです。
タグ:
posted at 07:31:12
@genkuroki matplotlibでヒストグラムのアニメーション作成時間を短くするには、histメソッドを使わずにPatchesと呼ばれる長方形を直接描画するかなり面倒な方法を使うようです。 matplotlib.org/gallery/animat...
タグ:
posted at 07:36:12
@genkuroki アニメーションは全然使わないのであまり詳しくはないのですが、おそらくFuncAnimationに渡す関数の中では各コマで更新したい部分をピンポイントでいじるのがコツなんだと思います。sinカーブを伸ばす例だとset_dataでLine2Dの点を一つずつ増やしていくアレです。
タグ:
posted at 07:45:15
@genkuroki Artistをピンポイントで更新していくにはplt.***で済ませる方式ではきついのでオブジェクト指向方式にしたほうがいいです。ax.plotやhistなどが作ったArtistオブジェクトはAxesが持つ”箱”に入っています。cf. qiita.com/skotaro/items/...
タグ:
posted at 07:51:38
@genkuroki 先にglobalで宣言しておいたfig, axesを使って直前のコマの長方形と線を消してhist&plotメソッドで毎回書き直す方針のままOO方式で書き換えてみましたが70秒から60秒になっただけでした。
タグ:
posted at 08:05:25
数学的にはODEはやり尽くされていて、これからはPDEになるという話を聞いたことがあるが、実際今どうなっているのか定かではない。
実装力雑魚勢としては、触れたくないものの一つではある。優秀で正しいと認められているライブラリが存在していればあるいは・・・
タグ:
posted at 08:17:11
@SING_A_WELL なるほど、PDEだと計算量的にも、数値計算の正当性を示すというよりは、数値的に解の存在を示すことが目的ということですね。一方でODEの場合は、前者の使い方を現実的に出来る程度には高速化しているということですかね。
タグ:
posted at 08:21:09
@SING_A_WELL PDEに関しては今後のブレイクスルー次第ということですかね。昨日精度保証の研究者に伺ったところによると、ODEは現在は多様な手法のライブラリが混在しているような状況だと仰ってました。そこからデファクトスタンダードが生まれてくれると使いやすいんですがね。
タグ:
posted at 08:33:34
#Julia言語 小ネタ
2次元ランダムウォークのプロットもこんな感じで簡単
using Plots
w = cumsum(randn(10000,2))
plot(w[:,1], w[:,2], legend=false, size=(400,400), lw=0.5) pic.twitter.com/cF8I6nCkk1
タグ: Julia言語
posted at 09:23:16
#Julia言語 小ネタ
3次元ランダムウォークのプロットもこんな感じで簡単
using Plots; pyplot()
w = cumsum(randn(10000,3))
plot(w[:,1], w[:,2], w[:,3], legend=false, size=(500,500), lw=0.5) pic.twitter.com/nImrJVwfue
タグ: Julia言語
posted at 09:25:07
バーゼル問題は2003年日本女子大の自己推薦入試で出たのも有名
examist.jp/legendexam/200...
バーゼル問題でウォリス使う証明は定番だが、sin^n xの定積分は、数研、東書の数3教科書の発展などの項目に紹介されている(ウォリスまではやってない)。高校の教科書の「発展」項目まで理解していると大学で楽
タグ:
posted at 09:38:12
高校数学の教科書は、数研(大島利雄ほかツイッタラー)・東京書籍(俣野博・河野俊丈ほか)のシェアが高く、あと実教(岡本和夫ほか)を含めた3社で寡占だったはず。難度の違う数種類が出てるが、各社最高難度の高校の教科書の内容を理解しておれば、大学1年の微積分の期末試験はもっとできたはずだ
タグ:
posted at 09:45:04
はてなブログに投稿しました #はてなブログ
Julia入門 - Strings - St_Hakky’s blog
st-hakky.hatenablog.com/entry/2018/02/...
タグ: はてなブログ
posted at 10:00:23
非公開
タグ:
posted at xx:xx:xx
#Julia言語 ニュース
twitter.com/search?f=tweet...
ValidatedNumerics.jlの作者のJuliaおじさんが大人気。
github.com/JuliaIntervals...
タグ: Julia言語
posted at 10:12:09
Wallis積は「偶数全部の積の二乗÷奇数全部の積の二乗=半円周の長さ」
ja.wikipedia.org/wiki/%E3%82%A6...
と覚えやすい。スターリングの公式、Γ(1/2)、ζ'(0)もWallis積がキーになって計算できる基本的な極限。原典は
John Wallis, Arithmetica Infinitorum 1656
archive.org/stream/Arithme...
タグ:
posted at 10:14:21
Arithmetica Infinitorumの頃は「π」もまだ使われてなかったし、微積分もNewton以前なので読みにくい。Wallisの原典の解読は
GA.Dickinson, Wallis's Product for π/2 参照
www.jstor.org/stable/3607585
Wallisの本自体の英訳が出版されている
www.springer.com/gp/book/978038...
タグ:
posted at 10:18:28
twitter.com/Paul_Painleve/...
に補足。以下のリンク先で実教出版の高校数学IIIの教科書でどのように
∫_0^{π/2} sin^n x dx
を扱っているかの写真があります。
twitter.com/genkuroki/stat...
教科書の著者達の意図は明らか。意図をくんだ教え方がされると素晴らしいと思う。
タグ:
posted at 10:21:45
ここまで教科書に載せることに成功したなら、ウォリスの公式の証明に使えるとか、ウォリスの公式はランダムウォークの逆正弦法則の証明に使えるとか、∫_0^{π/2} sin^n x dx の公式はベータ函数の特殊値になっているとか、……についてもついでに紹介してしまえばいいのに!
タグ:
posted at 10:23:46
∫_0^{π/2} sin^n x dx の計算のような素朴な事柄はガチ本物の数学の世界の入り口の一つであり、高校の数学の教科書を見るとそういう計算が大量に載っている。しかし、その計算のどこがどう面白いかに関する説明は皆無。
数学で最も重要なのは「何がどう面白いか」の部分だと思う。
タグ:
posted at 10:25:47
非公開
タグ:
posted at xx:xx:xx
非公開
タグ:
posted at xx:xx:xx
Wallis積の証明は∫_0^{π/2} sin^n x dx(これもWallis積分と呼ばれる)を求めて、nが偶数と奇数の場合から挟み込めばすぐ出る。前掲Dickinsonを参照。非自明な極限の例として微分の定義に必要な lim[x→0]sin x/x=1, lim[n→∞](1+1/n)^n=eの次に習うWallis積は汎用性が高く、微積分の華の一つ(終
タグ:
posted at 10:33:25
GR、さすが中性子散乱の研究所で作られてるだけあって二次元スピン配列の三次元画像を作るdrawspins関数がある pic.twitter.com/fYgVwZxKQX
タグ:
posted at 10:51:38
初期値で10^-12の摂動を与えているが、このルンゲ・クッタ法が6桁の精度なので、最初からルンゲ・クッタ法の誤差の方が支配的ではないのかという記事。
初期値鋭敏性と計算方式そのものの誤差拡大の切り分けをきっちりしないと、何を計算しているか分からなくなるという話。
verifiedby.me/adiary/0113
タグ:
posted at 11:53:25
こんな記事がホットエントリに上がってくるのがすごい/いま組んでるプログラムも全部倍精度で定義すればええやろ!って乱暴にもほどがある実装してるので,この辺きちんとやらないと / “二重振り子の精度保証付き数値計算 - kashiの…” htn.to/P1cgzF
タグ:
posted at 12:17:03
adhara_mathphys @adhara_mathphys
私は以前、象の方(tkmtsoさん)とおそらく同じ方法(6次のシンプレクティック陽解法、とても硬いバネで棒を表す、200桁くらいの多倍長精度計算)で三重振り子をときました。
twitter.com/adhara_mathphy...
タグ:
posted at 12:19:14
adhara_mathphys @adhara_mathphys
精度保証ではないので誤差を見積もったことにはなりませんが、初期状態が近いものについての差(確かこれは先っぽの位置についてだと思います)の拡大はこんな具合でした。
タグ:
posted at 12:21:07
adhara_mathphys @adhara_mathphys
真の解との誤差の拡大の仕方の目安にはなると考えています。
真の解とのズレはもっと大きいと思いますが拡大の仕方は近いのではないでしょうか。(1step後の誤差が大きいので上に平行移動した感じになると思われます)
タグ:
posted at 12:25:27
adhara_mathphys @adhara_mathphys
シンプレクティック法のおかげでエネルギーの変化が抑えられている、というのは本当です。
twitter.com/adhara_mathphy...
タグ:
posted at 12:27:50
adhara_mathphys @adhara_mathphys
そもそも剛体でやっていないのでその時点で真の振り子とは違いますが、それは触れないでくださいすみませんですね。
剛体でシンプレクティック法を入れるのは難儀だと思います。
タグ:
posted at 12:32:22
adhara_mathphys @adhara_mathphys
拘束条件つきのハミルトン系ということです。
剛体に対するラグランジアンからオイラー・ラグランジュ方程式は書けますが、これをハミルトン方程式に変えるのは厄介です。
ぜひ教えていただきたいところです。
タグ:
posted at 12:40:48
adhara_mathphys @adhara_mathphys
いや、ハミルトン方程式にはできた気がしますね。
その後、シンプレクティック陽解法で必要な、dp/dt=qの関数、dq/dt=pの関数のような変数分離ができていなかったことが障害だった気がします。
タグ:
posted at 12:47:24
adhara_mathphys @adhara_mathphys
最も簡単なシンプレクティック陰解法はシンプレクティック・オイラー法でしょう。 pic.twitter.com/ebnSwHsdHc
タグ:
posted at 12:58:31
twitter.com/Dsuke_KATO/sta...
#Julia言語 数学の教科書などでは行列を書くときにはコンマを使わずにスペースで区切って書く習慣になっています。Juliaはそれと同じです。
たぶん、Juliaは数学を理解している人にとって違和感がないスタイルを採用しているのだと思います。
タグ: Julia言語
posted at 13:00:09
adhara_mathphys @adhara_mathphys
Julia DifferentialEquations.jl のシンプレクティック法は全て陰解法ですね。
docs.juliadiffeq.org/latest/solvers... pic.twitter.com/esj6CETGC3
タグ:
posted at 13:04:23
@genkuroki (型を気にせずに)[1,2,3]のコマンドが通ったので、[1,2,3;4,5,6]は二次元配列になるだろうという予想が外れたので少し驚きました。
タグ:
posted at 13:06:12
非公開
タグ:
posted at xx:xx:xx
非公開
タグ:
posted at xx:xx:xx
@adhara_mathphys わかりきってるのに間違えるやつですね。笑
それはそうとこの実装みると面白いと思いますよ。シンプレクティック法ではないですが、陽的ルンゲ=クッタの実装みたらお世辞にも綺麗とは言えないコードでした。
タグ:
posted at 13:17:06
あとから何かを追加すればどんなプログラミング言語でもできるのは当たり前のことではあるのですが、ユーザー側がそういう面倒なことをしなくてもできるようになっていることが大事。
タグ:
posted at 13:31:18
@adhara_mathphys 自分で陽的ルンゲ=クッタ法を実装した時は外部ファイルから読み込むようにしたんですけど、気になってライブラリを見にいったら直書きだったんで、驚きました。直接書くメリットもあるので、こちらの方がいいのかも知れませんが、ちょっとデバックする気がおきないですよね。
twitter.com/ceptree/status...
タグ:
posted at 13:50:43
#Julia言語 縦ベクトル
1
2
と縦ベクトル
3
4
を横に並べて2×2行列を作りたいときには
A = [[1,2] [3,4]]
のような書き方ができます。空白は横に並べること。
a = [1,2]
b = [3,4]
A = [a b]
としても同じ。
vcatとhcatのドキュメントはあんまり詳しくない。
docs.julialang.org/en/v0.6.1/stdl...
タグ: Julia言語
posted at 14:09:20
非公開
タグ:
posted at xx:xx:xx
非公開
タグ:
posted at xx:xx:xx
#Julia言語 最新の DifferentialEquations.jl を使いたい人は色々注意が必要なようです。できるだけ最初にPkg.add("DifferentialEquations")して、ダウングレードに気をつける。
以上の記録を誰でも
nbviewer.jupyter.org/gist/genkuroki...
で閲覧できるようにしておきました。
タグ: Julia言語
posted at 17:11:11
はてなブログに投稿しました #はてなブログ, #Julia言語, #モデリング
Julia言語を使って常微分方程式を解く(捕食者-非捕食者モデル) - システムとモデリング
otepipi.hatenablog.com/entry/2018/02/...
posted at 17:14:49
Julia言語を使って常微分方程式を解く(捕食者-非捕食者モデル) - システムとモデリング otepipi.hatenablog.com/entry/2018/02/...
タグ:
posted at 17:46:24
#Julia言語 DIfferentialEquations.jl 4.0.0 を試したくて、こんなことをしてみたわけですが、よくわからない。
nbviewer.jupyter.org/gist/genkuroki...
を4.0.0で実行してみました。ハミルトニアンを H(p,q) から H(p,q,params=Float64[])に変更するだけで動いた。しかし、~続く
タグ: Julia言語
posted at 18:41:05
#Julia言語 続き~、色々面倒なので、
juliadiffeq.org/2018/01/24/Par...
に書いてある通り、
Pkg.pin("DifferentialEquations",v"3.1.0")
を実行して、一つ前のバージョンを使うようにしたいと思います。
4.0.0での実行結果は
nbviewer.jupyter.org/gist/genkuroki...
で公開してあります。
タグ: Julia言語
posted at 18:44:20
#Julia言語
DifferentialEquations.jl v3.1.0での実行結果
nbviewer.jupyter.org/gist/genkuroki...
v4.0.0での実行結果
nbviewer.jupyter.org/gist/genkuroki...
後者はかなり変です。仮に私のソース変更が足りないことが原因であっても、色々面倒なので一つ前のバージョンに戻るのが正解な感じ。
タグ: Julia言語
posted at 18:47:23
@genkuroki なるほどというか開発途上でめっちゃ都合上等主義で期待バンバンなんで当然だとは思うのですが、情報としては「Diff...の4.0」はまだ実用にはヤバいことが判明したので「目的があって使いたい人」はアップデートやPkg.addはしないのが吉という情報を端的に流して欲しいと。そのあとで様々な経験を語る
タグ:
posted at 19:00:12
#Julia言語 twitter.com/genkuroki/stat...
で逆正弦法則を5行のコードで数値的に確認できることを説明しました。
そこで扱った逆正弦法則は「公平なギャンブルにおいてトータルで得した状態になっている時間の割合が逆正弦分布に従う」というタイプの逆正弦法則です。続く
タグ: Julia言語
posted at 19:14:39
#Julia言語 トータルで最も勝っていた時刻の分布を5行で確認
using Plots
N, Niters = 10^4, 10^4
X = [findmax(cumsum(randn(N)))[2] for i in 1:Niters]/N
histogram(X, normed=true, bins=100, alpha=0.5, legend=false)
plot!(x -> 1/(π*sqrt(x*(1-x))), 0.003, 0.997, lw=3, color="red") pic.twitter.com/eZpDQDdN3A
タグ: Julia言語
posted at 19:24:09
#Julia言語 最後に勝ち負けがひっくり変える時刻の分布
添付画像を見ればわかるように6行で数値的に確認可能。 pic.twitter.com/1fYL5kHmuO
タグ: Julia言語
posted at 19:27:23
Juliaさん少しずつ分かってきた。module, export, using, importあたりとpush!(LOAD_PATH, “.”)で自作モジュールは管理できるのね。代数も直感的にコーディングできそうなんでちょっと楽しみになってきた。
タグ:
posted at 21:56:42
ちょっと前の #julialang v0.7 から Pkg.update() が depwarn になった。
リリースノートにはまだ書いていないっぽい。
v1.0 で初めに出会う depwarn はきっとこいつだな。 pic.twitter.com/m3MIa10sU3
タグ: julialang
posted at 22:21:28
#Julia言語
nbviewer.jupyter.org/gist/genkuroki...
植物・草食動物・肉食動物
1=植物は周囲に草食動物が存在すると草食動物に置換.
2=草食動物は周囲に肉食動物が2匹以上いると肉食動物に置換.
2=草食動物は周囲に植物が1つもないと植物に置換.
3=肉食動物は周囲に草食動物が5匹以上いないと植物に置換. pic.twitter.com/ebtwZh45Q9
タグ: Julia言語
posted at 22:34:39
#Julia言語 ツイッターに投稿した動画の画質が悪い。
nbviewer.jupyter.org/gist/genkuroki...
にはもっときれいなやつがあります。
Plotsパッケージ+GRバックエンド。動画作成に10秒かかりません!
タグ: Julia言語
posted at 22:36:35
@readonly_true 遅くなりましたがとりあえずカルマン渦発生するcsv吐くjuliaコードGitHubに挙げたのでご自由にお使いください。ただjuliaにしてはめっちゃ遅いのでかなり処理的に効率悪いことをやっている感があります。
github.com/koukyo1994/kar...
タグ:
posted at 22:51:53
@readonly_true juliaコード自体はjulia_scriptsに置いてあります。
必要環境はreadmeに書いておきます。
可視化もjuliaでやりたかったのですが、上手くいかないので一応matplotlibで可視化するコードもnbfiles/pyplot.ipynb内に入れてあります
タグ:
posted at 22:55:09