2013年10月28日月曜日

数学体験館

この10月にオープンした東京理科大学の数学体験館を見た. パラメトロン計算機や微分解析機や沢山の機械式計算器のある近代科学資料館の地下に ある. 直接地下に行ける入口も作られた. 入場無料. 館長は秋山仁さん. 開館時間は12時(土曜は10時)〜16時. 日曜・月曜・祝日は休み.

実はそれを見るのではなく, 近代科学資料館の所用があっていったら, 体験館はまだ30分は開いていると聞き, いそいで見て回ったので, そのうち改めてゆっくり見学したい.

写真を撮影してはいけないというので, ここに示す図はちらしをスキャンしたものである.





ちらしの絵から面白そうなものを拾うと


これは(1/4)+(1/4)2+(1/4)3+...=1/3の証明である.

正三角形が4つの相似な正三角形に分割され, 上の1つが更に4つに分割され,...た図である. 例えば右下の青の正三角形は全体の1/4であり, 上の同じ大きさの正三角形は全体の1/4であって, その右下の青の正三角形はその1/4だから全体の(1/4)2. このようにして青の正三角形は上の式の左辺になり, 青は赤や黄とともに3色で全体になるから, 青の面積は1/3というのである.

しかし展示は図だけなので, 証明は自分で考えなければならない.



これは最小公倍数・最大公約数算出器. 名前だけみると素晴しい計算機だが, 素因数分解は先に出来ていて, 例になっている90は2の球が1個, 3の球が2個, 5の球が1個で並んでおり, 24の方は2の球3個と3の球1個が並んでいる.

左右に分類器があって, 球をいれるとそれが分類されて出てくるから, それぞれの球の最大個数, 最小個数がわかれば公倍数や公約数が得られる.

これはどうもなぁという装置である.



球の体積の公式を導く話だ. 右にあるスイカをサッカーボールのように筋を入れて, 球の中心との母線で錐体に切り出す. 左にあるようにスイカの錐体が沢山できる. 底面は球の一部だがどんどん細い錐体にすると底面はほとんど平面になり, その錐体の体積は底面積に高さ(つまり球の半径r)を掛けた1/3である. 錐体全部の底面積はスイカの表面積だから4πr2. 従ってスイカの体積は4πr3/3である.

表面積を知っている人なら体積も知っていそうなものではあるが.



三平方スライド. 円のなかに直角三角形と各辺を1辺とする正方形が出来ている. 左の図は小と中の正方形に青と黄の板が入っているところ. 円板を回転すると青と黄の板が分割されスライドして大の正方形にうまく収まることをしめす.

昔ボストンのサイエンスミュージアムに同じ趣向だが水が入っているのを見たことがあった.



二項分布パチンコ. あまりよくは見えぬが右上が坂の上になっていて, パチンコの玉が沢山入っている. 出口を開けると一斉に左下に流れ落ちるが, 途中に通路をランダムにする釘が打ってあり, それによってパチンコ玉は下の区画のどれかに落ちる. その形が二項分布に見えるという実験である.

なにかの本でも見たような気がする.

まぁこういう実験道具が沢山おいてあって, 数学がこういうふうに体験できると分かって面白い.

最初のちらしのフロアーマップにあるように, 立派な数学工房がある. 体験館の実験具の多くはどうもそこで試作されたらしい. お願いすると使わせても貰えるらしい. なにか作ってみるものはないだろうか.

計算機屋からみると, ユークリッドの互除法, エラトステネスのふるい, ハノイの塔のようなアルゴリズム関係の展示もあればいいのにと思うが, すでに展示場は一杯である.

2013年10月17日木曜日

微分解析機

前回の微分解析機のシミュレータの続きである.

実物の微分解析機には出力装置があり, 2つの変数軸をxとyに入れると解の図形が得られる.

MIT Schemeにあるscheme graphicsの最も簡単な使い方は
(define mydevice (make-graphics-device 'x))
で描画の装置を宣言する. 装置の描く場所の正方形は左右がx軸, 上下がy軸であり, 左端はx=-1, 右端はx=1, 下端はy=-1, 上端はy=1である. この座標系で(x0,y0)から(x1,y1) まで直線を引くには
(graphics-draw-line device x0 y0 x1 y1)
とする.

そこで前回のブログのsinz, coszを描くには
(define (draw device step x0 y0 x1 y1)
 (if (> step 0) (begin
  (graphics-draw-line device x0 y0
   (stream-car x1) (stream-car y1))
  (draw device (- step 1) (stream-car x1) (stream-car y1) 
    (stream-cdr x1) (stream-cdr y1)))))
を使う. (draw ..)が描画の関数で, deviceに対してstep回描画する. x0, y0は解の曲線の始点; x1, y1はstreamである.

一方 駆動する方は, 座標を描いた後
(for-each (lambda (y) (graphics-draw-line mydevice -1 y 1 y))
  '(-0.75 0 0.75))
(for-each (lambda (x) (graphics-draw-line mydevice x -1 x 1))
  '(-0.75 0 0.75))

(draw mydevice 6282 0 0.75
 (scale-stream sinz 0.75) (scale-stream cosz 0.75))
のようになっている. stepは積分がzの1/1000で進行し, 2πで一周するから6282ステップ. 次の0と0.75はxとyの初期値(つまり鉛筆を置く位置)で, 1だと画面すれすれになるので0.75倍にしてみた.

次のsinzとcoszが出力軸の変数で, こちらも0.75にギアダウンしている. そういうわけで, この呼び方も出力装置をシミュレートしているといえるであろう.

こうして描いたサークルテストの図が次である.



円の始点と終点の真上でちょっと食い違っているようにも見えるが....

次はezを描く. 描きたい範囲を-1と1の間に変換しなければならない. 下の図を見てほしい. 左と上が実際の座標で右と下がGraphicsの座標である. これによるとx軸は0.4倍すればよい. y軸の方は1が-0.2になっているから0.4倍してから0.6を引く.



従ってプログラムは次のようだ. 始めのfor-eachは座標を描く.

(for-each (lambda (y) (graphics-draw-line mydevice -1 y 1 y))
  '(-0.6 -0.2 0.2 0.6))
(for-each (lambda (x) (graphics-draw-line mydevice x -1 x 1))
  '(-0.8 -0.4 0 0.4 0.8))
 (scale-stream z 0.4) 
  (add-streams (scale-stream ez 0.4)
    (scale-stream ones -0.6)))
回数の1386は丁度上の端に相当するxの値がlog 4=1.386であることによる.



こういう図を描くと, e-zで左半分も描きたくなる.
(define e-z (integrator (delay (scale-stream e-z -1))
  (delay z)
 1.0))
を用意し,
(draw mydevice 2500 0 -0.2
 (scale-stream z -0.4) 
  (add-streams (scale-stream e-z 0.4) 
   (scale-stream ones -0.6)))
を追加すれば, 左も描ける.



2013年10月16日水曜日

微分解析機

8月30日のブログに「これらのシミュレーションにはSchemeのstream処理が適していると思うが, まだ手が付かずにいる.」と書いた.

やっとどうやらシミュレーションが出来るようになったので, 今回はその話題である.

SICPの206ページのintegral関数のように披積分関数を遅延にするのが味噌である. SICPの例題では独立変数がdtだけだが, 微分解析機では「回す」値になにが入るかわからないので, こちらもstreamに対応させなければならない.

独立変数に方は回転角の差が必要になるので, stream-cadrからstream-carを引いたものに, 披積分関数のstream-carを掛けることになる.

まずintegratorは次のようだ.
(define (integrator delayed-integrand delayed-variable 
 initial-value)
 (define (mult-streams a b)
  (cons-stream (* (stream-car a)
            (- (stream-car (stream-cdr b)) (stream-car b)))
    (mult-streams (stream-cdr a) (stream-cdr b))))
 (define int
  (cons-stream initial-value
   (let ((integrand (force delayed-integrand))
         (variable (force delayed-variable)))
   (add-stream (mult-streams integrand variable) int))))
int)
integratorが貰う引数は遅延になっており, intの計算のなかでforceを使って強制する.

次はシミュレーションに必要な補助関数である. これらはSICPにあるものだ.
(define (scale-stream stream factor)
  (stream-map (lambda (x) (* x factor)) stream))
(define (add-streams s1 s2)
 (stream-map + s1 s2))
(define (integers-starting-from n)
 (cons-stream n (integers-starting-from (+ n 1))))
(define integers (integers-starting-from 1))
さてこれらが用意できたら, 8月14日のブログにあった関数からやってみる. ブログの図を見ながらプログラムを眺めてほしい.

最初はsin zとcos z. 独立変数はzである. zは整数のstreamで, 1/1000にして使う. sin zのπ/4, つまり45度での値をstream-refで見ることにする.
(define z (scale-stream integers 0.001))
(define sinz (integrator (delay cosz) (delay z) 0))
(define cosz (integrator (delay (scale-stream sinz -1))
  (delay z) 1.0))
(display (stream-ref sinz 785)) ;=>.7071024791302126
値としてはよさそうである.

次はez.
(define z (scale-stream integers 0.001))
(define ez (integrator (delay ez) (delay z) 1.0))
(stream-ref ez 1000) ;=>2.716923932235898
z2とz3は次のようにする.
(define z (scale-stream integers 0.001))
(define izdz (integrator (delay z) (delay z) 0))
(define z2 (scale-stream izdz 2))
(define iz2dz (integrator (delay z2) (delay z) 0))
(define z3 (scale-stream iz2dz 3))
(stream-ref z3 500) ;=>.12499950000000003
(stream-ref z2 500) ;=>.25049999999999994
0.5の3乗が0.125, 2乗が0.25ならまぁ問題はない.

tan zは多少手強いが,
(define ones (cons-stream 1 ones))
(define z (scale-stream integers 0.001))
(define tanz (integrator (delay 1+tanzsq) (delay z) 0))
(define itanzdtanz (integrator (delay tanz) (delay tanz) 0))
(define tanzsq (scale-stream itanzdtanz 2))
(define 1+tanzsq (add-streams tanzsq ones))
(stream-ref tanz 785) ;=>.9979528633948679
(stream-ref tanz (quotient 3142 6)) ;=>.5761873674775257
テストしたのはtan π/4とtan π/6で, 1/√ 3=.5773502691896258だからこれも合格である.

1/zとlogezはzを1から積分する. z=1での初期値はどちらも0である.
(define z (scale-stream integers 0.001))
(define logez (integrator (delay 1/z) (delay z) 0))
(define 1/z (integrator (delay (scale-stream 1/z -1)) 
   (delay logez) 1))
(stream-ref 1/z 500) ;=>.6664863143078377
(stream-ref logez 1718) ;=>.999948126092027
1/zの500での値は1/1.5なので0.6666.... logezの1718はloge2.718だから1でいいわけである. zの値は1からのはずだが, 回転角の差を使うから, integersで構わない.

微分解析機には結果を描画する出力装置もついているから, 次はscheme-graphicsを駆動するようなプログラムにしたい.