2013年4月6日土曜日

微分解析機

微分解析機の原理は分かっているし, Processingでシミュレータも書いたりしたが, 定量的にはどうか.

佐々木達治郎「計算機械」に解説はあるが, かつて東大理工研や阪大理学部で活躍し, 今は理科大で静態保存されている昭和航空計器製の微分解析機の数値をもとに自分で計算してみた.

下の図は積分機2台で正弦波d2y/dx2=-yを描く場合の機素の接続図である.



I0とI1が積分機で, Axは独立変数の縦軸(バスシャフト), xのA倍で回転する. By,-Dy''はI0の出力軸, I1の被積分関数軸でyのB倍で回転する. Cy'はyの1回微分の軸でI0の被積分関数軸, I1の出力軸だ. ダッシュ'が使えないプログラムではyをy0と, y'をy1 と書くことがある.

図でAx軸のすぐ下の斜め線はその縦軸の回転を, 1:1で積分機I0の変数の横軸(クロスシャフト)へ伝えるはすば歯車. その先の平歯車でm0倍になり円板を回転する. 平歯車の先の軸の回転角と円板の回転角は1:1で固定である.

Cy軸から積分機I0へ向かう横軸, 被積分関数は, n0倍になってから1回転が1ミリのピッチで円板と回転子の相対位置を変える. 図ではI0とI1の円板の位置が上下に揃い, 回転子がずれているが, トルクアンプの関係で, 通常は回転子が上下揃うように構成されている.

回転子は半径rで, j0の平歯車を経て出力軸と結ばれる.

下の積分機I1でも接続は同様だが, Byの軸の斜線の方向が違うのは, 被積分関数軸が縦軸に対して1:-1になっているからだ.

出力卓については次回のブログで説明しよう.

上の積分機I0で考えると, 回転子の円板の中心からの距離lはCy'n0ミリ. その状態でAx軸がdAxだけ回転すると, 円板がdAxm0だけ回転するから, 回転子の位置で円板は円周方向にdAxm0Cy'n0動き, 出力の横軸はdAxm0Cy'n0/rj0(=dBy)だけ回転することになる.

つまり, この値がByに足される.

この微分解析機にはまだ制約がある. rは32ミリ. 円板の半径は約120ミリ. 平歯車のギア比は4, 2, 1, 0.5, 0.25しかない.

そこでA, B, Cやギア比をいろいろ試行錯誤で計算してみると, m0=m1=2, n0=n1=0.25, j0=j1=4, A=256, B=C=448という使えそうな解があることが分かった.

もともとyやy'は-1≤y≤1だから, n0とn1が0.25だとBy,Cy'は, -112ミリから112ミリの間になり, 円板の範囲に収まる.

xを0から2πまで変えるとAxは0から512πまで変わる. yとy'の初期値はそれぞれ0と1だから, ByとCy'の初期値は0と448だ.

この1周の積分を1024分割で実行すると, dAxはπ/2である

これでシミュレートしてみる.
(define m 1) (define n 0.25) (define j 2)
(define a 256) (define b 448) (define c 448)
(define ax 0) (define by0 0) (define cy1 448) 
(define pi (* 4 (atan 1))) (define dax (/ pi 2))
(define r 32)
(define (integrator dx y) (/ (* m dx n y) r j))
(define (step i)
 (if (= (modulo i 16) 0)
  (begin (display (list i ax by0 cy1
    (+ (* (/ by0 b) (/ by0 b)) (* (/ cy1 c) (/ cy1 c)))))
   (newline)))
 (let ((dcy1 (integrator dax (- by0)))
       (dby0 (integrator dax cy1)))
  (set! ax (+ ax dax))
  (set! by0 (+ by0 dby0))
  (set! cy1 (+ cy1 dcy1))))
(do ((i 0 (+ i 1))) ((> i 1024)) (step i))
1024ステップ分出力するのは大変だから1/16つまり65行プリントした. 左からステップ番号, Ax, By, Cy', y2+y'2である.
(   0    0.00000    0.00000  448.00000 1.00000)
(  16   25.13274   43.92436  445.97712 1.00060)
(  32   50.26548   87.45205  439.65678 1.00121)
(  48   75.39822  130.16351  429.09730 1.00181)
(  64  100.53096  171.64681  414.39783 1.00241)
(  80  125.66371  211.50158  395.69748 1.00302)
(  96  150.79645  249.34290  373.17400 1.00362)
( 112  175.92919  284.80504  347.04204 1.00423)
( 128  201.06193  317.54493  317.55119 1.00483)
( 144  226.19467  347.24554  284.98349 1.00544)
( 160  251.32741  373.61893  249.65085 1.00604)
( 176  276.46015  396.40903  211.89195 1.00665)
( 192  301.59289  415.39414  172.06908 1.00725)
( 208  326.72564  430.38907  130.56463 1.00786)
( 224  351.85838  441.24697   87.77739 1.00847)
( 240  376.99112  447.86075   44.11879 1.00908)
( 256  402.12386  450.16415     .00887 1.00968)
( 272  427.25660  448.13236  -44.12771 1.01029)
( 288  452.38934  441.78236  -87.86579 1.01090)
( 304  477.52208  431.17271 -130.78379 1.01151)
( 320  502.65482  416.40305 -172.46777 1.01212)
( 336  527.78757  397.61316 -212.51543 1.01273)
( 352  552.92031  374.98162 -250.54001 1.01334)
( 368  578.05305  348.72414 -286.17396 1.01395)
( 384  603.18579  319.09147 -319.07260 1.01456)
( 400  628.31853  286.36704 -348.91733 1.01517)
( 416  653.45127  250.86423 -375.41882 1.01579)
( 432  678.58401  212.92338 -398.31976 1.01640)
( 448  703.71675  172.90852 -417.39737 1.01701)
( 464  728.84950  131.20387 -432.46556 1.01762)
( 480  753.98224   88.21016 -443.37676 1.01824)
( 496  779.11498   44.34079 -450.02335 1.01885)
( 512  804.24772     .01783 -452.33875 1.01946)
( 528  829.38046  -44.33200 -450.29802 1.02008)
( 544  854.51320  -88.28149 -443.91821 1.02069)
( 560  879.64594 -131.40702 -433.25816 1.02131)
( 576  904.77868 -173.29266 -418.41798 1.02192)
( 592  929.91143 -213.53415 -399.53811 1.02254)
( 608  955.04417 -251.74286 -376.79800 1.02315)
( 624  980.17691 -287.54947 -350.41438 1.02377)
( 640 1005.30965 -320.60761 -320.63922 1.02439)
( 656 1030.44239 -350.59717 -287.75730 1.02501)
( 672 1055.57513 -377.22738 -252.08351 1.02562)
( 688 1080.70787 -400.23970 -213.95984 1.02624)
( 704 1105.84061 -419.41026 -173.75205 1.02686)
( 720 1130.97336 -434.55206 -131.84624 1.02748)
( 736 1156.10610 -445.51682  -88.64505 1.02810)
( 752 1181.23884 -452.19640  -44.56390 1.02872)
( 768 1206.37158 -454.52385    -.02688 1.02934)
( 784 1231.50432 -452.47414   44.53723 1.02996)
( 800 1256.63706 -446.06438   88.69916 1.03058)
( 816 1281.76980 -435.35369  132.03323 1.03120)
( 832 1306.90254 -420.44265  174.12149 1.03182)
( 848 1332.03529 -401.47238  214.55776 1.03244)
( 864 1357.16803 -378.62318  252.95148 1.03306)
( 880 1382.30077 -352.11282  288.93159 1.03369)
( 896 1407.43351 -322.19447  322.15002 1.03431)
( 912 1432.56625 -289.15431  352.28509 1.03493)
( 928 1457.69899 -253.30872  379.04465 1.03556)
( 944 1482.83173 -215.00134  402.16889 1.03618)
( 960 1507.96447 -174.59970  421.43285 1.03680)
( 976 1533.09721 -132.49175  436.64863 1.03743)
( 992 1558.22996  -89.08210  447.66722 1.03805)
(1008 1583.36270  -44.78813  454.37993 1.03868)
(1024 1608.49544    -.03601  456.71951 1.03931)
axの最後の値は (* pi 512) => 1608.495438637974だ.

この出力をもとにサークルテストをしたのが下の図だ.



開始点と終了点である右端がちょっとずれている. 448と456だから2%くらいの誤差である.

0 件のコメント: