werry-chanの日記.料理とエンジニアリング

料理!コーディング!研究!日常!飯!うんち!睡眠!人間の全て!

Pythonでマンデルブロ模様変化の動画作ろうや

久しぶりに何かプログラミングな記事を書こうと思いまして,

そういえば以前,マンデルブロ集合を中途半端に説明して放置してたなぁと

werry-chan.hatenablog.com

↑これですね.ソースコードだけのせて,しかも画像で,これはよろしくない.

今回はちゃんとコード載せますので安心してください.

ほんで今回は,

マンデルブロ集合が動きます.

マンデルブロ集合って複素数座標上の原点スタートが標準なのですが,今回は原点以外から反復関数かけてやると,どのように模様が変化するか見てみます.
www.youtube.com
ゆっくり動くマンデルブロ集合↑

www.youtube.com
早く動くマンデルブロ集合↑


制作過程をtwitterで追っている↑

後半で沢山ウニョウニョしてますね.

M=~~ってのを大きくすると計算するのにかなり時間かかかりますので注意です.

必要なライブラリはnumpyとopencvです.

これがマンデルブロ集合の画像生成プログラムです.↓

import numpy as np
import cv2

def model(X,Y,N,a,b):
     for i in range(N):
         a,b=a**2-b**2+X,2*a*b+Y
         c=a**2+b**2
     c[c<4]=0
     c[c>4]=255
     return(a,b,c)

def mand_plot(a,b,M):
    for N in range(2,20):
     #print("now process is step ",N)
     mdl_tpl=model(X,Y,1,a,b)
     a,b=mdl_tpl[0],mdl_tpl[1]
     if N//5%2==0:
         img[mdl_tpl[2]==255]=N%5*41
     else:
         img[mdl_tpl[2]==255]=255-N%5*41

#1000枚の画像をmand_plptというファイルに保存するのでfilename定義.
filename="./mand_plot/"

#画像の枚数1000枚
img_num=1001
#画像のサイズは1000*1000.後でresizeして良いが,計算中はサイズそのままじゃないとプロットの精度が悪くなるので,保存時にresizeしましょう.
M=1000

#極座標形の長さrの初期値
r=0.0
#極座標形の角度thetaの初期値
theta=0.0
for l in range(1,img_num,1):
    #rとthetaを少しずつ増加させる.thetaは全部で5pi回転する.
    r+=1.0/(img_num-1)
    theta+=5*np.pi/(img_num-1)
    a=r*np.cos(theta)
    b=r*np.sin(theta)
    x,y=[np.linspace(-2,2,M)]*2
    X,Y=np.meshgrid(x,y)
    mdl_tpl=model(X,Y,1,a,b)
    a,b,img=mdl_tpl[0],mdl_tpl[1],mdl_tpl[2]
    mand_plot(a,b,M)
    #画像を保存前にresize
    img=cv2.resize(img,(500,500))
    if l<10:
        cv2.imwrite(filename+"mad"+"00"+str(l)+".png",img)
    elif l<100:
        cv2.imwrite(filename+"mad"+"0"+str(l)+".png",img)
    else:
        cv2.imwrite(filename+"mad"+str(l)+".png",img)
    print("now process is ",l,"\n")

画像を動画に書き出すプログラム.↓

import cv2

fourcc = cv2.VideoWriter_fourcc('m','p','4','v')
#(500,500)は画像のサイズ.
video = cv2.VideoWriter('video_c4_M1000.mov', fourcc, 20.0, (500, 500))

for i in range(1, 1001,1):
    if i<10:
        img = cv2.imread('mad'+"00"+str(i)+'.png')
    elif i<100:
        img = cv2.imread('mad'+"0"+str(i)+'.png')
    else:
        img = cv2.imread('mad'+str(i)+'.png')
    #img = cv2.resize(img, (500,500))
    video.write(img)

video.release()

ますらば東大模試(4)解いてみた

エオ!!

エオ!!

リーロリロリロリロレロ!!

レーロ!!

エーオ!!

ALL RIGHHT!!

遅れ遊ばせながらボエミアンラプソディー見てきました.

フレディ最高ですね.超カッコいい.クイーン大好き.

めっちゃライブ行きたくなりましたわ.


閑話休題


前回までの続きです.

暇な時間の楽しみに少しずつ解いて楽しんでます.

論文書かないとなぁと思いつつ,ぼちぼちだらだらして,高校数学の入試問題解いてるウェリーちゃんです.

前回までの道筋はこちらから辿ってください.

werry-chan.hatenablog.com




問4

辺の長さが で,辺が鏡になっている正方形 ABCD があり,線分 BC 上に
BP : CP = t : (1 −t ),(\frac{1}{2} < t < 1)
となるような点 P を定める。点 A から点 P
に向けて,光線 K を発射し,初めて K により囲まれてできる四角形の面積を S(t)
とする。

(1) K が線分 CD DA 上で最初に反射する点を,それぞれ Q ,R とする。このと
き,\vec{QP}・\vec{QR} を求めよ。

(2) S(t) の最大値と,そのときの t を求めよ。







以下に答えを載せますので,注意です.






































(1)の解答.


図形の設定を紙面左上を点Aとして,時計回りに点B,C,Dが配置する正方形を考える.
また,以下の解答では点Aを原点として\vec{AB}方向をx軸正方向,-\vec{AD}をy軸正方向と考える.

\vec{AP}=(1,-t)

反射の性質を考えると,

\bigtriangleup ABP \sim \bigtriangleup QCP \sim \bigtriangleup QDR  ①

と相似関係が求まる.

相似非は上記①式の三角形の順に

 1:\frac{1-t}{t}:1-\frac{1-t}{t}

と求まる.

以上の関係がつかめれば求める値は容易に計算できる.

 \vec{QP}=(\frac{1-t}{t}, 1-t),  \vec{QR}=(1-\frac{1-t}{t}, 1-2t)

ここで角度関係を求めたいので

\angle BAP=\theta,\angle PQR =\pi -2\thetaとして

 cos(\angle PQR )=-cos2\theta=sin^2\theta -cos^2\theta=\frac{BP}{AP}^2-\frac{AB}{AP}^2=\frac{t^2}{1+t^2}-\frac{1}{1+t^2}

と表せるので

\vec{QP}・\vec{QR}=|QP||QR|cos(\angle PQR )=-\frac{(1+t)(1-t)^2(2t-1)}{t^2}



(2)の解答

点Rで反射した光が線分APと交わる点を,点Sとする.

四角形PQRSは,(1)で求めた相似関係から平行四辺形とわかる.

ゆえに,

 S(t)=\bigtriangleup RQP \times 2 = |\frac{1-t}{t}(1-2t)-(1-t)(2-\frac{1}{t})|=-4t+6-\frac{2}{t}

S'(t)=-4+2\frac{1}{t^2}

これらの式から,t=\frac{1}{\sqrt{2}}で極大値を持ち,

S(\frac{1}{\sqrt{2}})=2(3-2\sqrt{2})



間違ってるとか,別解あったら教えてください.

カヌレ(ホワイトデー用)

Hololensの自分担当のモジュール開発がやっとこさ一段落ついて嬉しいウェリーちゃんです.

実際には,謎のbuildエラー吐かれてて,最後のbuildエラーを処理する必要があるんですがね...

気づけば,もうあと1週間で4月ですよ.春休み消滅ですよ.全く.

しかも4月ってことは,新4年である僕はラボ配属ですよ.タンパク質研の正式配属ですよ.

え?僕のこと情報系と思ってたって?

いやいや,僕は高校までは生物系で二次胚誘導の特殊事例研究してましたよ.

大学での正式な所属は材料開発系学部で量子力学統計力学を主に勉強してますね.

まぁ今の正式所属研究室はタンパク質で非公式所属研究室が情報系ですね.

意味不明な経歴持ちのウェリーちゃんです.

閑話休題


本日は最近ホワイトデー用に作ったお菓子第二弾です.

カヌレ

前回のカトルカール↓に続いて,今回もフランス洋菓子です.
werry-chan.hatenablog.com


f:id:werry-chan:20190323092630p:plain


材料(プリン型6個分)

バター 70g

牛乳 500ml

ラム酒(洋酒で代用可) 35ml

薄力粉 65g

強力粉 35g

砂糖 250g

卵黄(僕は全卵でやります) 4個


①バター50gと牛乳の半量(250ml)を温めて,バターを溶かす.

②手順①にラム酒全量・残りの牛乳を入れて,常温程度まで冷やす.

③薄力粉・強力粉・砂糖をふるっておいたものに手順②を数回に分けて加える.ダマにならないように混ぜる.ゴムベラでなく,泡立て器を使うことを推奨.

④手順③に卵黄を加え,残りのバター20gを焦がしバターにして加える.均一になるよう混ぜる.その後,時間あるなら冷蔵庫で一晩置いておくと良い.

カヌレの型にバター・はちみつを塗布して,手順④を7分目を目処に加える.
(僕はカヌレ型ないのでプリン型でしました.)

⑥220度で45-60分程度焼く.
(僕のオーブン200度までしか上がらないので,オーブントースターで具合見ながら焼きました.上下火1200wで45-60分くらいで焼けました.オーブントースターで焼く場合は,アルミホイルで型の上面に蓋をしておくと消し炭にならない程度に焦げないので良い.)


f:id:werry-chan:20190323092630p:plain


カヌレの完成です.

今回のレシピに書いてある通りに作れば,もっとお手本のようなカヌレが完成するはずです.

僕は手抜きで,手順⑤のバターを勝手にサラダ油にしたり,ハチミツつけ忘れたりしたので微妙な出来でした.あとオーブンの温度が220度まで上がらなかったので,オーブントースターで作ったりしたし...めちゃくちゃですね.

そもそもカヌレって蜜+蝋で作るという洋菓子だという記述が存在するのに,蜜も蝋も適当に扱うなんて言語道断です.

反省します.次はちゃんと作ります.

みんなは上のレシピ守って作って,美味しいカヌレ食べてね.

ますらば東大模試(3)解いてみた

だいぶ半月前の話になるが,大阪に帰省して奈良とか京都とか行って楽しんでました.
つくばに帰ってきたので仕事の続きをしないといけません.あと実験しないとそろそろ次の国際会議が〆切が見えてきて怖い.


暇つぶしに解いてるますらばの問題(3)です.


前回までに解いてきた(1)
ますらば東大模試(1)解いてみた - werry-chanの日記


問いの(2)です.
ますらば東大模試(2)解いてみた - werry-chanの日記




それでは問題文です.


n を正の整数とする.
(1) p を素数とする.このとき,n^p-n は p で割り切れることを示せ.
(2) n^2 + 1 は 4 で割ると 3 余る素因数を持たないことを示せ.



僕,(2)がむずすぎて解けませんでした.
以下に正解を載せますので注意です.































解答です.



まずは(1)です.

帰納法をしたいと思います.

変数がnとpの二つありますね.

これは周知の事かと思われますが,全ての素数pを示す式はありませんので,nを帰納法の変数として動かしてみましょう.


A. n=1で1^p-1=0で割り切れます.


B. kを正の整数として,k^p-kがpで割り切れると仮定する.


C. (k+1)^p-(k+1)について,pで割り切れるか検討する.

(k+1)^pという式の形をみて,一番はじめに思い浮かぶものは,二項定理ですね.二項定理で展開してみましょう.

(k+1)^p=\sum^p_{i=0} {}_pC_i k^i=\sum^{p-1}_{i=1} {}_pC_i k^i +(k^p+1)

すなわち,

(k+1)^p-(k+1)=\sum^{p-1}_{i=1} {}_pC_i k^i +(k^p+1)-(k+1)=\sum^{p-1}_{i=1} {}_pC_i k^i +(k^p-k)-(1+1)


ここでBより,第二項はpで割り切れる事を仮定している.また第三項は0である.

第1項については

{}_pC_i について考えると,{}_pC_i は組み合わせの式なので整数であり,素数pは必ず割り切れない事から{}_pC_i はpの倍数となる.

以上よりCもpで割り切れることがわかった.


ABCより,数学的帰納法より題意は示された.




(2)の解答です.

これ激ムズですね.

正直に言いましょう.解けませんでした.

悔しいですが,公式の模範解答を載せます.

(2) p を 4 で割ると 3 余る素数とする。このとき,n^2 + 1 が p を素因数に持つことは,n^2 + 1 が p で割り切れることと同値である。

(i) n が p で割り切れるとき

n^2 + 1 ≡1 より n + 1 は p の倍数ではない。

(ii) n が p で割り切れないとき明らかに p ≧ 2 である。

x の多項式 x^p − x を x + 1 で割った余りは高々一次式なので ax + b とおけて,その商を Q (x) とすると

x^p − x = (x^2 + 1)Q (x) + ax + b   (2)

x^2 + 1 の各項の係数はすべて であるから,

Q (x) の係数はすべて整数である。

この式に x = i(虚数単位) を代入すると

i^p −i = (i^2 + 1)Q (i) + ai + b

p は 4 で割ると 3 余る素数であるから,各辺を整理すると

−2i = ai + b

両辺の実部と虚部をそれぞれ比較して,a = −2, b = 0 を得る。

したがって,n^p − n = (n^2 + 1)Q (n) − 2n が成り立つ。

(1) より n^p − n は p の倍数である。

また,n は p の倍数でないので 2nも p の倍数でない。

したがって (n^2 + 1)Q (n) は p の倍数でなく,多項式

Q (n) が整数係数であることに注意すると n^2 + 1 も p の倍数ではない。

(i)(ii) より,n^2 + 1 は 4 で割ると 3 余る素数 p で割り切れない。

すなわち,n^2 + 1 はそのような素数を素因数に持たない。

(証明終)



別解や間違ってるなどありましたらコメントくださいな.

むずかったなぁ(2)...解けんかったの悔しいわ.

カスレもどき(カスレはフランス料理だった気がする)

Hololensまじで意味わからんので怒ってます.

WebCam映像取得しようとするとわけわからん挙動して困ってて,てかそもそもC#がよく分からない.

Hololens開発しながら,全く別件の国際会議のジャーナル書いてるとか意味わからんですね.

春休みとは?


f:id:werry-chan:20190322023119p:plain

今回は男料理です.

てか,毎日作るもの思いつかなくて,そういう時にいつも作ってるので詳しい分量はよく分からなくて感覚で作ってます.

まぁ分量で大事なのは玉ねぎ・人参の量が1:1ってことくらいですね.


材料

玉ねぎ (人参と同じくらいの量)

人参 (玉ねぎと同じくらいの量)

ニンニク (明日会う人によって調整)

オリーブオイル 適量(炒めるのに使う程度)

任意の肉(何でも良い) (食べたいだけ)

白ワイン (なければミリンとか日本酒でOK)

パン粉 (オシャレにしたい人か食感楽しみたい人だけ)

出汁 フォンドボー推奨(中華スープの素,鶏ガラでもOK)

①鍋にオリーブオイルにニンニク刻み(チューブニンニクで可)入れて,ニンニクの香りを移すように軽く炒める.
 (いつも5人前以上作るので,最初から最後まで鍋で作業してます.少人数分であれば好みの大きさの鍋使ってください.)

②玉ねぎ,人参を頑張れる限りのミジン切りする.
 (僕はめんどいのでブレンダーでグシャグシャにします.)

③刻んだ玉ねぎ,人参を①へ入れて,炒める.人参と玉ねぎのスープが炒めることで染み出してきたら塩・胡椒・出汁を適量加える.
 (僕は胡椒は断然あらびき派です.あらびき胡椒万歳です.ミルで挽きたてあらびき胡椒は絶品です.)
 この時点で,もう既に美味いので味見してみても良いです.

④肉を食べやすい大きさに切り,塩・胡椒をふっておく.(塩味を馴染ませるために,塩をふってから5-10分程度放置してあると良い.)

⑤肉を鍋に入れて,炒めて,表面に火が通るようにしておく.

⑥白ワインを「ひたひた」まで入れる.
f:id:werry-chan:20190322025030p:plain
これが「ひたひた」です.
(画像はhttps://cookpad.com/cooking_basics/8027より参照)

⑧蓋して煮込む.すぐ食べたければ,グツグツして15分すればOK.
 (鶏肉と豚肉は感染症多いので,芯部温度80℃以上15分以上が原則.)

⑨オシャレに仕上げるには,皿に盛った後に,パン粉ふって,トースターで焦げ目つけて完成.
f:id:werry-chan:20190322023119p:plain


これはかなり美味しいくせに全然手間がなくて良い.

ほぼ玉ねぎと人参切って炒めて煮てるだけ.

アレンジでトマト缶入れても良かった.

これはフランス料理の簡単なスープの作り方って感じで良い.

カトルカール(フランス語で4分の4ケーキ)

最近あんまりアップロードしてませんでしたね.

ダラダラしたり,論文書いたり,Hololens開発したりと忙しかったんですよね.

こないだバレンタインあって,その時彼女にケーキをもらいましたので,そのお返しに作ったケーキが3つあるんですよ.

今回はその一つ目を紹介しましょう.

カトルカール.

f:id:werry-chan:20190321191526p:plain


フランス語で4分の4という意味らしいです.

何となく冲方丁のマルドゥックシリーズを思い出しますね.

マルドゥックシリーズめっちゃ良いのでオススメの小説です.

カトルカールはその名の通り,主な材料を1:1:1:1の比率で混ぜ合わせて作るところから名前がつけられています.

主な材料とは,卵・砂糖・小麦粉・バターです.

これらはケーキの基本ですね.

それではレシピを紹介していきましょう.

材料(18cmホール)

全卵 3個
砂糖(グラニュー糖) 120g
バター 120g
薄力粉 120g
トリモリン(ハチミツで代用) 13g
オレンジジュース 25g
オレンジの皮(あればで良い) 1個分以下
アプリコットジャム 適量


①全卵を全て・砂糖全て・トリモリン全て,をボウルに入れ,軽くときほぐす.大体均一にまざりあえばOK.


②泡立て器で泡立てる.電動ですることを強く奨めます.手でやる奴は筋肉です.ただの筋トレです.電動買って使いましょう.

 大体,「リボンが描ける」状態になるまで撹拌します.

f:id:werry-chan:20190322015651p:plain

このくらいです.(画像は https://www.lotte.co.jp/products/brand/ghana/technique/make3.html より)


③次は別の小さめの鍋を出してください.バター全量を弱火で溶かします.オレンジの皮をすりおろした物をバターへ加えます.


④オレンジの香りがバターに少し移ったかな?と思ったらオレンジジュースをそこへ加えます.

(注意!!バターを熱し過ぎていると突沸爆発事故になるので,オレンジジュースを加える際はバターの温度が100℃以下になっていることを確認.)

オレンジジュースを加えたら,火を止めて,人肌程度の温度に冷やす.


⑤手順②までの生地にふるった小麦粉を全量加えます.混ぜ方は底から返すように,ゴムベラを上手く使って混ぜましょう.


⑥手順④を⑤に加えてバターの油脂性のムラが出ないように(油が浮いてこないように)なるまで混ぜる.


⑦型に流して焼きます.型にバターを塗ったら,完成後外しやすいです.焼く温度は160℃,20分です.
(オーブンの性能によっては温度と時間が変わります.小さいオーブンの人は温度を10℃高め,あるいは時間を長めにするなど工夫してみてください.)
割と膨らみますので,型の深さの7分目までにすることを奨めます.僕は一杯にして焼いて溢れてオーブンめっちゃケーキまみれにしました.掃除しんどかったです.


⑧焼けたら,温めたアプリコットジャムを塗って完成です.


f:id:werry-chan:20190321191526p:plain

オレンジの香りがふわっとして美味しいです.

これお好みでレモンで作っても全然OKですね.

毎回ケーキ作ってる時に砂糖の量入れてみてドン引きしてるんですよね.こんなあるん!!怖って.

ホワイトデー用に作った残り二つのケーキは

マドレーヌ

です.

お楽しみにしててください.

ますらば東大模試(2)解いてみた

筑波山の梅祭りに行ってきて,休憩所で飲んだ梅昆布茶が美味しくて買ってしまったウェリーちゃんです.
3月3日に筑波山梅祭りの梅酒試飲会があるそうですが,その日は大阪に帰ってるかもです.行きたかったので残念です.


さて前回に続き「ますらば東大模試」第2問です.


werry-chan.hatenablog.com
↑問1の僕の解いた解答はこちらです.




それでは問題に入りましょう.

問2


nを正の整数とする.赤色と白色の 2 種類の球が十分な数ある.横一列に並んだn 個の箱それぞれに,赤白どちらかの球を最大1個ずつ入れていく事を考える(球を入れない箱があっても良い).このとき,どの隣り合う 2 つの箱にも同色の球が存在しないような球の入れ方は何通りあるか.


下に解答を載せてますので,注意です.































[解答]

k個の箱が並んでいる時,
k番目が空の場合がa_k通り,赤の場合はb_k通り,白の場合はc_k通り存在すると仮定する.


k+1番目の箱を追加すると,
k番目の空の箱の次に並び得る物は,空or赤or白
k番目の白の箱の次に並び得る物は,空or赤
k番目の赤の箱の次に並び得る物は,空or白

すなわち,

a_{k+1}=a_k+b_k+c_k---①
b_{k+1}=a_k+c_k---②
c_{k+1}=a_k+b_k---③

とk+1番目の箱について考えられる.

ここで②-③式は,

b_{k+1}-c_{k+1}=-(b_k-c_k)

となる.これは左辺と右辺は同型なので

b_{k+1}-c_{k+1}=-(b_k-c_k)=(-1)^2(b_{k-1}-b_{k-1})=...=(-1)^k(b_1-c_1)=0

ここで,b_1=1,c_1=1を利用した.
以上より,b_k=c_kが導けた.(ぶっちゃけ言うと,対称性より明らか.)

b_k=c_kを①, ②式に代入して,

 a_{k+1}=a_k+2b_k---④
 b_{k+1}=a_k+b_k---⑤

④式を

b_k=\frac{a_{k+1}-a_k}{2}

と変形して,⑤式に代入するとa_kに関する式

a_k=\frac{a_{k+2}-a_{k+1}}{2}-\frac{a_{k+1}-a_k}{2}

整頓して

a_{k+2}-2a_{k+1}-a_k=0

となる.

特性方程式を解いてやると,

a_{k+1}=\frac{\beta^k(a_2-\alpha a_1)-\alpha^k(a_2-\beta a_1)}{\beta -\alpha}---⑥

ここで,特性方程式の解を

\alpha =1+\sqrt{2}, \beta=1-\sqrt{2}

とおきました.

最後に,a_1=1, a_2=3を利用して⑥整理すると

a_{k+1}=\frac{(1+\sqrt{2})^{k+1}+(1-\sqrt{2})^{k+1}}{2}

が導かれる.

求めるべき通り数はa_n+b_n+c_nであり,①式を利用すると,

a_n+b_n+c_n=a_{n+1}=\frac{(1+\sqrt{2})^{n+1}+(1-\sqrt{2})^{n+1}}{2}

とわかった.

解答終了.



間違ってるとか,別解あるよーとか,コメントあればよろしくです.