伊勢志摩スカイラインを走ってきた。

日が暮れて見にくいですが。

なかなかきれいな景色でした。しかし鹿さんにいっぱい（4頭）もであいましたよ。

黒木メイサさんのコーラショックCMソングはジミークリフのアイ・キャン・シー・クリアリー・ナウ

最近やってるコーラショックのCM。

http://www.kirin.co.jp/brands/colashock/cm_gallery/index.html

CM曲はJimmy Cliffの I can see clearly nowですね。

オリジナルはこっちのJohnny Nashのやつ。こちらもいいですよ。

二見シーパラダイスへ行って来た（伊勢神宮の帰り）。

あっかんべーをするアザラシで有名ですが、

セイウチのショーも面白かったよ。なにしろでかい。

大阪の女の子がキスされてた。こんな経験めったにできるもんじゃないよね。

夫婦岩を見てきた（伊勢神宮の帰り）。

写真で見るとでかいんだが、結構こじんまりしてました。

かえるが多かったよ。触るとまた帰ってこれるらしい。

伊勢神宮に行って来た（その他）。

超でかい鯉と、

なぜか鶏が放し飼い。

で、おかげ横丁で、赤福の本店へ。お茶と3個の赤福餅セットがあります。

手作りできたての赤福餅はやはり違うなあ。

あとレトロな郵便局とか。

なっちゃん新CMの女の子は桜庭ななみちゃん/鳥取砂丘で撮影だそうです。

砂漠っぽいところを大勢の女子高生が走っているというなっちゃんの新CM。

http://www.suntory.co.jp/softdrink/natchan/cm/index.html

出演は桜庭ななみちゃんですね。撮影場所は鳥取砂丘ですって。

http://www.suntory.co.jp/softdrink/natchan/cm/diary.html

伊勢神宮に行って来た（次は内宮）。

ちょっと歩いてまずは正宮。撮影は階段の下からだけですよ。

で荒祭宮。

風日祈宮。

こちらのお馬さんは白でした。

（続く）

伊勢神宮に行って来た（まずは外宮）。

家族4人で久々に行ってきました。まずは正宮。内部の撮影は禁止ですよ。

近くになぞのパワースポット？

お馬さんもいました。こっちは黒。

で多賀宮。

土宮。

風宮。

 

と拝んだところで次は内宮だ（続く）。

平成22年2月22日22時22分をISC(Inverse Symbolic Calculator)で。

2並びということでこの時間に更新。それだけでは面白くないのでISC (Inverse Symbolic Calculator)で2.222222・・・を推定してみた。

http://oldweb.cecm.sfu.ca/projects/ISC/ISCmain.html

・1/45　 = 0.0222222...

・2/9 = 0.2222222....

がでたね。これでもう決定打ですな。いまいち面白くなかった。

「魍魎の匣」を読んだ。

これは面白かった。前作はちょっと納得できーーん、みたいなとこもあったんですが、これはすごく納得（技術的にできるわけではないけれど。。。）。

鯛を骨だけにして泳がす料理を思い出した。包丁人味平で見たんだが。。。

「虐殺器官」を読んだ。

前から気になっていて、文庫化されていたので読んでみた。日本人が書いた長編SFを読むのはこれがはじめてかもしれない。

大変面白かった。ナノマシンの実用化や脳の機能がある程度解明されたくらいの近未来のお話。雰囲気は「地獄の黙示録」+近未来SFかな。大規模ハンニバルレクターかもしれないが。しかし9.11以降のテロ/内戦をまともに扱ったお話は私初めて読んだ気がする。核戦争とか、情報管理とかものすごくリアルでありそうなお話。

それはそうと解説なし（？私が読み飛ばしただけかも）で出てくるAKはこのカラシニコフ自動小銃ですかね。

http://ja.wikipedia.org/wiki/AK-47

ラストは衝撃って言えば衝撃ですが、アメリカが気に入らない人にはすっとするかも（アメリカ人には書けない話。。。）。

欲を言えば母親のくだりとジョンの元カノへの固執がもう少しかけてたらもっとよかった。あと器官というか、脳のその部分がもっともっともらしい説明が詳しかったらとか。

しかし作者の伊藤さんは若くしてお亡くなりになったということで、今後が楽しみだっただけに残念。海外に翻訳しても受けそうなお話だったのに。。。

いろんなブラウン運動（確率微分方程式）をExcel VBAで計算

前に幾何的ブラウン運動と代数的ブラウン運動の話を書いたけれど、もっと一般にはどうなんだろうと、

dX(t) = μX dt + σ X^n dW

の計算をしてみた。

乱数：メルセンヌ・ツイスタ+Box Muller法でアルゴリズムはミルシュタインスキーム。

いろんなnで計算した結果がこちら。

nがマイナスだとXがでかくなると揺らぎが小さくなるから普通に大きくなりそうなのはわかる。nが大きくなると逆にどっかで押さえ込まれるのかな？もっとばらつきそうなものだが。。。

参考文献はこちら。

https://amzn.to/3H9vMFK

水の上って走れるのか！

東寺の堀にて撮影。ほんとに水の上を走ってたよ。

確率共鳴をExcel VBAで計算してみる。

この「微分方程式による計算科学入門」の確率微分方程式のところのっていた例題。

地球の百万年程度のスパンでの温度変動を表す簡単な数理モデルとして、

dX(t) = [ X (1-X^2) + a cos(ωt)] dt + σ dW(t)

が上げられている。ここでX(t)は地球の気温で、-1が氷期、1が間氷期を表す。

計算はオイラー・丸山スキームで。上からσ=0.1,0.2,0.3,0.4の時の計算結果。

なるほど。こういう風に確率共鳴が起きるのか。

ゼルニケ多項式(Zernike Polynomials)をkeisan.casio.jpにUP!

よく光学の分野で使われるゼルニケ多項式の自作式を作ってkeisan.casio.jpにUPした。

これです。

↓

ゼルニケ多項式(Zernike Polynomials)

式はこんな感じ。

計算結果はこんな感じ。

もうこれで何個だっけ。100個超えてから数えるのがめんどくさくなった。

キリン コクの時間の新CMソングはクワイ河マーチですね。

最近流れている松たか子さんのCM。映画「戦場にかける橋」のテーマでもあるクワイ河マーチですね。しかもミッチ・ミラーって確かオリジナルを演奏してる。

http://www.kirin.co.jp/about/toku/ad/kokunojikan/index.html

これって、よく替え歌で

「さる、ごりら、チンパンジー」

と歌われるやつですが、私が感心したのはドランクドラゴンの塚地氏が

「さる、ごりら、一般人」

という替え歌にしてたこと。すごいアイデア。

ソフトバンクCM（旅行編）の高校生は福井で選ばれたエキストラですって。

夫婦で旅行中に3人の女子高生に会うというCMですが、

http://mb.softbank.jp/mb/campaign/3G/

（CMソングは金子由香利の「時は過ぎてゆく」）

ロケ地は福井の一乗谷朝倉氏遺跡だそうです。でかわいい高校生は福井で選ばれたエキストラさんだって。

http://www.fukuishimbun.co.jp/modules/news2/article.php?storyid=9792

純朴でかわいらしいですね。ちょっと上戸さんにも似てる人が選ばれたんでしょうね。

「3点を通る円」「4点を通る球」をkeisan.casio.jpにUP!

このサイトに4点を通る球の話で検索されてきた方がいるのですが、そういやあんまりそういう計算みないな、と思って自作式を作成しました。

リンクはこれです。

4点を通る球

球　(x-a)^2 + (y-b)^2 + (z-c)^2 = r^2　で、球上の4点を入力すると(a,b,c,r)を計算します。

ついでに3点を通る円も。

3点を通る円

円　(x-a)^2 + (y-b)^2 = r^2　で、円上の3点を入力すると(a,b,r)を計算します。

プログラムは出すほどのものでもないですが、要望があれば（5次方程式の解は要望があったので、ソース公開した）。

東寺に行って来た。

近鉄の東寺駅で降りてすぐです。大人一人800円なり。五重塔の中の金の像もいいんだけれど、金堂と講堂の中の巨大金の像、特に金堂の薬師如来像には感激。全然どんなものかしらなくて入ったらすぐ巨大なのがありました。写真撮影、スケッチなどは禁止ですよ。

正規乱数をkeisan.casio.jpにUP!

最近確率の計算を続けているのでおまけとしてこんなのつくってみた。

↓

正規乱数

線形合同法+Box-Muller法を使っているだけ。

プログラムはこれ↓

av = 平均;
sig = 標準偏差;
n = 乱数の個数;
seed=乱数の種;
for (j = 1 ; j <= n; j=j+1) {
    seed=mod(1566083941*seed+1,2^32);
    r=seed / 2^32;
    t = sqrt(-2.*ln(1.-r));
    seed=mod(1566083941*seed+1,2^32);
    r=seed / 2^32;
    u = 2.*pi*r;
    println(j, t*cos(u)*sig + av, t*sin(u)*sig + av);
}

まあ乱数のアルゴリズムがこれなんで、実用性はないかもしれないが。

山崎駅でオーバーラン/いつまで復旧にかかるのか。

昨日、JR京都線に乗っていたら、山崎駅をとまるときに列車がオーバーラン。ちょっと後戻りするのに５分かかった。やっちまったなー、また日勤教育かな、とか思ってたら、いつまでも動き出さない。放送では運転手が連絡をしてるので待ってくださいとかあって、これも10分くらいかかった。何を連絡してるのかとか思ってたらようやくのろのろ運転。近くの踏み切りがずっとしまっていたらしく待っている人が大勢いた。と思ったら！急停止。本部から急停止信号がでたらしく、どうも踏み切りの確認ができてない、とか言う話。もう通り過ぎたよ!

ということで、安全性の確認はぜひ徹底的にしてもらいたい（私はJR福知山線で田舎に帰るので、あの痛ましい事故の場所をずっと通っていたので人事じゃない）ですが、もう少し効率てきにできないですかねえ。

確率ダフィング方程式をExcelで計算

下の参考文献ではなぜかドゥフィンク振動子になってますが（本当の発音はわからん）。

dX(t) = Y dt

dY(t) = ( X (a-X^2) - Y) dt + σX dW

3種類の初期条件で計算。σ=a=1。

 

平衡点がふらふら動いてる。

参考文献はこちら。

LismoのCMの女子高生は川口春奈ちゃん。

ポカリスエットのCMのときはほんとに子供という感じでしたが、成長されましたねー。

最近のLismo Fesのはこちら。

http://www.kddi.com/lismo/fes/?did=ad_lismo-to-fes_00001#/tvcm

だいぶ前の学生服のはこちら。

http://www.kddi.com/lismo/#/cp/homeroom/h_cm

ビッグ・ナラはそれほどビッグでない。

あちこちに看板があってさぞでかいのかな、と思ったらそうでもないなあ。

でもきれいなお店ですよ。

奈良でスガキヤさんを見つけた。

大安寺近くのイオンにて。関西にもあるのか、知らなかった。ラーメンフォークで食べるんでしたよね。

http://www.sugakico.co.jp/

ってwebサイトみたら京都にもあるじゃないの。今度行って見よ。

ちなみにここで食べてないのはその前にちかくの「がっつり亭」で黄身丸弁当を食べたからです。一見、敷き詰められたご飯の上に目玉焼きが乗っているだけなんだが、そのご飯の下におかずがあるというもの。ただ、今回入っていたのは豚のしょうが焼き、から揚げが半分とあとは焼きそば！炭水化物だらけで正直ご飯とのバランスが取れなかった。。。

確率微分方程式をExcel VBAで（幾何的ブラウン運動と代数的ブラウン運動）。

前から疑問に思っていたこと。物理で出てくるのはたいてい、ランジュバン方程式みたいに

dX(t) = μX dt + σ dW

のようにノイズがXの関数でないもの。これは加法的ノイズとか、代数的ブラウン運動というみたい。ところが金融工学とかブラックショールズで出てくるのは

dX(t) = μX dt + σX dW

のようにノイズがXの関数になっている。これは乗法的ノイズとか、幾何的ブラウン運動とかいうみたい。

で、すごく金融工学が人工的に見えて私があまりこの手の話を好きでない理由のような気がする。

Xにさえ依存しなければ伊藤微分もストラトノビッチ微分も一緒で、めんどくさいことが一切なくなるというのに。。。でも数学的にはXに依存しないと面白そうなことはあまりなさそうかな。

とりあえず計算結果はこちら。

乱数：メルセンヌ・ツイスタ+Box Muller法でアルゴリズムはミルシュタインスキーム。

うーん、しかしｄX/X　（変化率）　へのノイズの影響一定という物理的な状況もないわけではないだろう。なにかあるだろうか。

参考文献はこちら。

https://amzn.to/3H9vMFK

dBの平均と平均のdB

先日も書いたけれど、再度詳しく書いてみる。dBで表される量、たとえばインサーションロスとか考えてみよう。

PdB = 10*log P

とか。で、これを測定値として平均を取るとき、なにも考えずに

(PdB1 + PdB2 + .... + PdBn) / n

を計算する人がいる。でもこれは間違いで、実際は真値のPに直して

10*log ((P1 + P2 + .... + Pn)/n )

とするのが正解。それは

(PdB1 + PdB2 + .... + PdBn) / n = 10 * log{ (P1*P2*...*Pn)^(1/n) }

になってしまって相加相乗の定理により、

(P1 + P2 + .... + Pn)/n ≧ (P1*P2*...*Pn)^(1/n)

で同じ値になるのは測定値が全部同じ場合のみだから。

しかしなぜこういう間違いが多いかというと、これが大きくずれるのは、ばらつきが大きい場合で、ばらつき少なかったら大して変わらない。実際に計算してみた結果がこちら。真値の平均mと標準偏差σの比で平均値がどう変わるかを計算してみたもの↓。

10%くらいだったらどっちの平均もかわらない。てことで注意しない人が多いのでしょう。しかし、ちょっとでもロスを少なくしたい人が間違った計算してたら是非直しましょう。

もちろん、数十%とかでかいばらつきがある測定は論外。

でもこれちゃんと計算で示せるかと思ったらできなかったのでシミュレーションに逃げてしまった。計算力低下もはなはだしい。。。

「眩暈」を読んだ。

北半球と南半球で洗面所で水を流したときの回り方の向きが違うというのは迷信なわけで。。。（コリオリ力はそんなミクロなサイズでは効かない。台風レベルじゃないと。。。）　でも赤道をまたぐ観光地で、洗面器の水の回り方が違うのを見せるパフォーマンスをする人がいるらしいですが、それは見てみたい（もちろんトリック）。というのはネタバレなのだろうか。黒猫館の殺人とどっちが先だったんだろう、というのもネタバレなのだろうか。でもそこらへんはもう本の前半で解決してるからいいよね。

確率ロジスティック方程式をExcelで（アルゴリズム比較）。

先日はオイラー・丸山スキームをつかってみて計算したロジスティック方程式の確率版。

dX(t) = X (α-βX）dt + σX dW

乱数：メルセンヌ・ツイスタ+Box Muller法。

今度はミルシュタインスキームだ。

σ=1.0

で計算。

あんまりかわらんなあ。パラメータが差が出ないようなものなのかな。

参考文献はこちら。

RFワールドNo.9を買った。

なんだかんだと毎号買ってるなあ。今月の特集は電波伝搬の基礎。

まあそれはいいんだけれど、TWTってまだつかわれてるのか。知らなかった。あとケーブルテレビの技術って今回初めて知った。

それよりなにより、コラムであったんだが、dB単位の平均を単純にとるな（和をとって数で割るな）という話があった。これはよく間違えている人がうちの会社でもいてよく注意することがある。

たとえば、デシベル単位でないP1, P2, ... , Pnの平均は

Pav = (P1 + P2 + ... Pn)/n

でいいんだけど、これをdBでやっちゃうと

(10logP1 + 10logP2 + ... +10logPn) /n= 10log(P1*P2*...*Pn)^(1/n)

ってこと。で、10logPavとこれを比べると、

相加相乗の定理:

(P1 + P2 + .... +Pn)/n   >=   (P1*P2*...*Pn)^(1/n)

だからdBで平均とると必ず真値に直したものより小さくなる。でここで終わればいいんだが、まだ先があって、それは

「P1,...Pnのばらつきが小さかったら問題になるほど差がでない！」

（しかも全部同じ値なら＝がなりたつ）

ということ。なのでいつまでも平均のとり方をこれでやっている人が減らない。。。

PLANT3に行って来た。

福知山のかき末に行くついでによってきた。

http://www.plant-co.jp/

なかなか大きな１フロアのスーパーで、アメリカにいたころのKマートとかウォルマートを思い出したよん。

リアル「漁師力学」をはじめて見た。

研究者になりたくて博士課程までいったけれど全然就職なくて、行方不明になったりフリーターになったり。。。という話はたまに聞く。

で、物理で博士課程いって、漁師さんになった人がさっきテレビ見てたらNTVの「満天！青空レストラン」にでてた。ほんとにみたのははじめてだ！漁師力学というおやじギャグがリアルにできる！

ブログもみつけたよ（リンクは略）。

何があったんだろうね。。。

確率ロジスティック方程式をExcel VBAで計算。

先日、確率微分方程式の数値計算について調べたのだが、とりあえずオイラー・丸山スキームをつかってみることにする。

方程式はロジスティック方程式の確率版、

dX(t) = X (α-βX）dt + σX dW

だ。乱数はここは当然メルセンヌ・ツイスタ+Box Muller法で（このブログの過去ログ参照）。

σ=0, 0.1, 0.5, 1.0, 1.5, 2.0

で計算したものがこちら。

σが大きくなるともう絶滅してるなあ。これはこれで面白い。

参考文献はこちら。

PS3/Wii/DSi LL/Wii Fit plusを一度に買った。物欲爆発。

42型のプラズマテレビを買ったら、つなぐものがほしくなって一度に買ってしまった。やりすぎか。。。

しかし金を持ってる人が使わないと景気は回復しない！これでいいのだ。

確率微分方程式をExcelで数値的に計算する（準備編）。

年初の予告通り、確率微分方程式

dX = f(t,X)dt + g(t,X)dW

を数値的に計算しようと思う。しかし、ほとんど文献がないなあ。日本語で読めるものはこの「微分方程式による計算科学入門」くらいかな。

https://amzn.to/3H9vMFK

当初考えていたのは、dWのない

dX=f(t,X)dt

はこのブログの看板ともいえる、Dormand-Prince(Runge-Kutta8次）で計算できるから、それにちょこっとノイズを加えりゃいいんじゃないの？とか思っていたが、上の本読んで気づいたのは、「全然ちがーーーーう！！！」

考えれば当たり前で、そもそも dW っていわば√(dt)　みたいなもんだから、そもそも次数の勘定が全然違う。

で、まず普通のオイラー法に相当するオイラー・丸山スキームは

X_n+1 = X_n + f(t_n, X_n) Δt + g(t_n, X_n) ΔWn

です。ここでΔW_n=ξ*√Δt　で、ξは平均0、分散1のガウス乱数。

√がでてることでもわかるように、これでも収束次数は1/2！！！

もう少し精度がいいものとしてミルシュテイン･スキームも紹介されている。

X_n+1 = X_n + f(t_n, X_n) Δt + g(t_n, X_n) ΔWn + [∂g/∂x g](t_n,X_n)*(ΔWn^2-Δt)/2

だが、これでやっと収束次数1！！！！

こりゃ大変だ。しかしfに関する高次のものは長時間積分するわけじゃないからべつにいらないのか。

ではこれを使っていろいろ計算してみよう（続く）。

春日神社で初詣

だいぶ前に撮った写真のお蔵だし。

結構有名な能舞台があるよ。

「数学を生み出す魔法のるつぼ」を読んだ/ISCをはじめて知った。

表題の本を読んだのだが、この本全般にわたってなんか日本語訳が変じゃないか？

まあそれはいいとして、この本に出てくるISC (Inverse Symbol Calculator)が面白い。

リンクはこちら。

http://oldweb.cecm.sfu.ca/projects/ISC/ISCmain.html

ISC+ (ISC plus)もある。ただし本に出てくるリンクは古くてこっちになったみたい。

http://glooscap.cs.dal.ca:8087/

何ができるかというと、たとえば

19.999099979

を入力すると、

Best guess: exp(Pi)-Pi

と返してくれる。数値データから、推測してくれるルックアップテーブルですね。

これを使って自分の今年の誕生日（2.0100614）を推測してもらおう。そのままぴったりというのはなかったが、

Best guess: Im((20+13*I)^(19/11))

これが近いとでた。201.0061347...ってことらしい。

今年のバレンタインデー2.0100214で推測すると

Best guess: 3+14*x-4*x^2+3*x^3

の解だって。

本当か調べると、-0.2010021411101だ。これの一部ってことか。すごいのかなんなのかわからなくなってきた。

適当な数をだすからいまいちなのかな。では勝手な関数を作って逆算してもらおう。

たとえば、

e^π+2*π-eなんてどうかな。→　(ln(5)+5)/(-Pi+2/3) と出た。微妙だ。

どうも本に載ってるみたいにばっちりとした答えがない。なんかやり方間違ってるのかな。

capsule/Joe Pass/Casiopea/SantanaのCDを買ってきた。

Capsuleは初めて聴いた。こういうのはやっぱり好き。Perfumeが歌ってもよさそうな曲もあって、逆にPerfumeでも手を抜いてないってことですね。

フォー・ジャンゴは昔探してて全然見つからなかったのがようやく買えた。タイトル曲は私が昔買ったジャズギターの教則本で楽譜だけあってオリジナルを聴いたことがなかったという。。。

カシオペアは黄金時代の名曲がずらり。ドライブのお供に。

サンタナの邦題「哀愁のヨーロッパ」は私がギター初めて一曲通して弾けた思い出の曲ですよ。