2011年9月25日日曜日

再帰曲線

このブログにdragon曲線のことを書いたのは, 2年以上も前であった.

dragon曲線は, 紙を同じ向きに半分半分と折り, 折り目を90度になるように広げたものであった.

最近, この折る向きを毎回逆にしたらどうなるかと思った. 紙を折ってやってみるのは, やはり3,4回が限度であり, 手元の計算機の威力を借りたくなった.

下の図を見て欲しい.



左端のAは折る前である. 下の黒丸は出発点を示す. その上の方に右向きの矢印が示すように, 最初は右に折る. するとBのようになる. 紙の長さが倍にになっているが, 今は折り方が問題なので, 分かりやすく描いてある.

出発点から途中右折したので, その角には+が付けてある.

Bを今度は左へ折る. 最初の折り目の+は左の下へ移動する. Cの出発点から辿ると, 最初は左折, 次がBで出来た右折, さらに右折する.

Cを次は右へ折るとDになる. 上の方の+や-は, 今回出来たのもで, 下の方のは左がCの下にあった折り目, 右がCの上にあったものだ.

出発点からの右折左折を抜き出すと

B +
C - + +
D + - - + + + -
この伝でつづけると次は
E - + + - - - + + - + + + - - +
となる.

一段上の列の要素を間に+と-を交互に入れたものになっているが, dragon曲線の場合と同じで, このパターンに着目する.

EをX3, DをX3とすると, X3の左半分はX2の+-を反対にしたものだ. 右半分は中央だけが違って右半分を殆んど同じである.

そこで, 'で+-の反転を表わすとすると,
X3=X2'+Y2'.

結局
X0=+, Y0=-,
Xn=Xn-1'+Yn-1',
Yn=Xn-1'-Yn-1'
となることが分かる.

Schemeでプログラムしてみる. この引数のnは上の漸化式の2n-1になっている.

(define (r s) (if (eq? s '+) '- '+))

(define (x n i) (cond ((= i n) '+)
((< i n) (r (x (/ n 2) i)))
((> i n) (r (y (/ n 2) (- i n))))))
(define (y n i) (cond ((= i n) '-)
((< i n) (r (x (/ n 2) i)))
((> i n) (r (y (/ n 2) (- i n))))))

(map (lambda (i) (x 8 i)) (a2b 1 16))
=> (- + + - - - + + - + + + - - +)


これで右折左折の情報が得られたから, いよいよ交互折りの曲線を描くことする.

そして出来たのが次の図である. X1からX8までが, 色を変えながら描いてある. 左中ほど上の, 直角に下の曲がった青の線がB(X1)である.



左から右向きに出発するのが, 出発点である. その下の緑の線がC(X2)に対応する.

折角描いてみたが, 竜にも鳳凰にもならず, 単に三角形が出来ただけであった. 思うにdragon曲線を考えた人も, これもやっては見たが, つまらない結果だったので, このことは書いて置かなかったのかもしれない.

つまらない絵しか描けないことが分かっただけでも, 一応の知見が得られたというべきか.

2011年9月2日金曜日

整数立方根

平方根に較べ, 立方根はおおいに疎外されている感があるが, 立方根を計算したいこともある.

その場合, 実用的には
(expt 9270 (/ 1 3)) => 21.00680051861007
なことで済ませてしまう.

さて, Warrenの「ハッカーの楽しみ」を眺めていたら, 整数立方根という話題があった(訳書の226ページ). こういうプログラムが書いてある(32ビット用だ).



同書の巻頭の「推薦の辞」に私が書いたように, 本書にはそのアルゴリズムでよい理由が殆んど書いてない.

このアルゴリズムは下の説明のように出来ているのだ.

例によってSchemeに書き直すと




最初のsの設定が30でなく18なのは, そんなに大きい数でテストすることもないからだ. また最後にyだけでなく, xも返すのは, 余りも欲しいからである.

途中の経過を見るdisplayが3行ある. 9270の立方根を計算してみる.



すなわち, 立方根は21で, 剰余が9. 整数の範囲で開立をやめ, 剰余が得られるから, 整数立方根というわけである. 同じものを十進法で書くと,



つまりbの始めの方は, xより小さくなるまで218, 215, 212,... と減っていき, 4096がxより初めて小さくなった時点でyを1とする. yは最後は二進で10101だが, その最初の1(=24)が決ったのである.

xの方はそれを引いて, 9270-4096=5174になった.

次はyのその下のビットを1にするか, 0にするかを決めることである. y=16の時, その下のビットは8だから(16+8)3=163+3*162*8+3*16*82+83=13824

しかし, 最初の4096はすでに引いてあるので, 比較するbは, 上の4段目の13824-4096=9728である. こうするよりは, y=2, (y+1)3-y3=3y(y+1)+1を作り, 左シフトする方が楽というのが, プログラムの趣旨であろう.

8のビットは立たなかった. 次は4のビットで, (16+4)3-163=3904である. これはxより小さいからyは1増えて101になった.

このようにしてyのビットを上から次々を立てていき, 1ビットごとにきめていく. 最後に残ったxが開立の剰余である.

このプログラムの場合は, 方針は最初から見え見えであったが.

2011年8月28日日曜日

多面体描画道楽

Archimedes多面体13種のうち, rhomb(斜方とか菱形)という形容詞のつくものが2つある. 斜方立方八面体と斜方二十十二面体である. 何ゆえに斜方といわれるのか.





これらの図に示すように, この2種類の多面体には, それぞれ3種類の面がある.

斜方立方八面体では, 青は元々の立方体の面, 緑は正八面体の面である. 斜方二十十二面体では, 青は元々の正十二面体の面, 緑は正二十面体の面である.

では, それぞれの赤の面はなにか.

赤の面を延長して, 緑の面の上方での交点を決め, 赤の正方形に外接する菱形を考えることが出来る. 赤は元々そういう面の多面体の面である.

この多面体がどちらも菱形になるので, 斜方といわれる所以である. つまり, 斜方(赤)立方(青)八(緑)面体, 斜方(赤)二十(緑)十二(青)面体という命名であったわけだ.

ではその斜方多面体を描いてみよう.


上は斜方立方八面体である. x, y, zの対称軸も描いてある. z軸(青)のまわりに少し回転し, y軸(緑)のまわりにも僅かに回転したようになっている.


そのz軸, y軸まわりの回転をやめ, x軸方向から眺めたのが, この上の図だ. そこで, 最初の赤い正方形を含む面を考えてみると, この図の赤線で示す菱形が得られる. その座標が分かれば, 次の図が描けるわけだ.



これは, 対角線の比が1:√2の菱形が12個で構成されているので, 斜方十二面体(rhombic dodecahedron)という.


斜方二十十二面体の方は, 面が多いせいか, 多少手ごわい.



これが元の斜方二十十二面体で, 上と同様にz軸, y軸について回転してある.



回転角を0とし, x方向から眺めると, このように見えるはずだ. 赤線はそこにかぶせた菱形である. この方は, 菱形30枚で構成され, 斜方三十面体(triacontahedron)という.



菱形の対角線の比は1:φ(黄金比)である.

同一の菱形だけで出来る多面体はこの2つだけらしい. またKeplerはこの2つの多面体のあることを知っていたといわれる.

2011年8月27日土曜日

3シリンダ機関車

ずいぶん前のこのブログに3シリンダ機関車のことを書いた. そこで引用しているURLが, 最近なくなっているのに気づいたので, もう一度3シリンダ機関車の弁装置について書くことにする.

今でも蒸気機関車は人気があり, 各地で復活されているが, その構造はあまり理解されていないに違いない.

簡単にいってしまえば, 機関車の左右の前方につけた2個のシリンダ(気筒)の中にピストンを入れ, ピストンの両側から交互に蒸気を補給してピストンを前後に動かす. ピストンの往復運動を主連棒で繋ぎ, 動輪を回転して機関車を進める. ピストンがどちらかの端に寄っていると, 蒸気の補給が出来ず, 死点といって力が出せないが, 反対側のピストンが最も力の出る中央にいるように, 左右のピストンの位相を90度づらせて配置してある. このことは, 鉄道博物館で実物の機関車を見ると確認出来る.

ピストンの前後に交互に蒸気の補給をしているのが, シリンダと並んだ蒸気室で, その中に滑り弁があり, これが蒸気室の中を前後に滑りながら, 蒸気をピストンの一方に送り込む.

Walschaertの走り装置の図はここにあるので参照されたい.

滑り弁は, 動輪に動力を伝える主連棒と90度づれたエキセン棒から, 逆転装置を経由した弁棒で押し引きされて運動する. この辺の装置は大変に込み入っているから, 説明は省略!

ところで, 表題の3シリンダ機関車は, 2つのシリンダでは力不足の時, 第3のシリンダを両シリンダの中央に搭載し, 合わせて3本の主連棒で動輪を駆動する方式である. 動輪軸の中央はクランク状になっていて, そこに中央の主連棒が接続されている.

問題は, 左右の動輪のエキセン棒のようなものが, 中央のクランクからも出ているかということだが, そういう複雑な走り装置にはなっておらず, 中央のシリンダの弁の制御は両側の弁の制御から構成していたというのが前回のブログであった.

下の図を見て欲しい. Cylinder AとCylinder Bと書いたのが, 両端のシリンダの蒸気室で, 左(機関車の後方)から弁棒が入り, 途中に滑り弁があり, 弁棒の延長が右方に突き出している. 機関車の先頭にAOD, BDCのような梃子があり, 支点AはCylinder Aの弁棒の先に, 支点BはCylinder Bの弁棒の先に固定されている. 梃子だから, 図で左右に揺れると支点の上下移動もあるはずだが, 今は各支点は左右にだけ動くとしている.

支点Oは固定点である. またOD=OA/2, つまりAが動くとDはAの半分の振幅で逆位相に動く. そのDを支点として, DB=DCの梃子があり, CがCylinder Cの弁棒になっている. そうすると, CはAとBが120度の位相差で動くと, さらに120度の位相差で動き, 3番目の滑り弁が望み通りに動くのである.

このアニメーションがここにある.

アニメの下でくるくる回る3芒星は, 120度の位相で回転することを示す.

2011年8月26日金曜日

ベクトルの和

「とっておきの数学パズル」に, 多面体の各面の, 直交外向きで, 面積に比例する大きさのベクトルの和が0になることを証明せよというのがあった. Q5.13

多分そうなりそうだが, 証明となると....

「いかにして問題を解くか」には, 「一様な四面体の重心を求めよ」という問題が解けなければ, 平面の類似の問題「一様な三角形の重心を求めよ」を考えてみるという教訓がある.

その伝でいくなら, 多角形の各辺の, 直交外向きで, 辺の長さに比例する大きさのベクトルの和が0になることを証明せよとなる.

たとえば, 下の図の左は, 最も簡単な多角形の三角形の場合で, 辺の長さをa,b,cとする. (5:3:4になっている.) 従って, 外向きのベクトルは, 黒い矢印のように出来る. ベクトルの和を計算するのだから, aの先端までbの付け根を移動し(青矢印), その先端までcの付け根を移動する(青矢印)と, 3つのベクトルで構成される三角形は, もとの三角形と相似になり, cの先端はaの付け根に一致し, 和が0なことが分かる.



しかし, こういう方法はどうだろうか. (図の右)

元の三角形の左右方向をx軸. 上下方向をy軸とする. 各ベクトルのx, y成分を考える. bはx成分しか, cはy成分しか持たない. aは, 図にax, ayと書いて示したx成分とy成分がある. ベクトルaの大きさが, 斜面aの長さに比例するなら, axはaの斜めの面の, x方向からみた断面(青線で示す)の長さに比例するわけで, この断面の長さはbの辺の長さに等しく, 従ってaxの大きさはbの大きさを同じで方向が反対である.従って, ベクトルのx成分の和は0.

同様にして, y方向の成分の和も0. 従って, ベクトル全体の和も0になる.

これを3次元の問題でも使ってみよう. ある面の, その面積に比例した直交外向きベクトルaのx成分は, ベクトルとx軸のなす角をθとすると, a cosθ, この面のx方向から見た断面積もa cosθのはずである.



したがって, この多面体の直交外向きのベクトルのx+方向の和は, x+方向から見た断面積になり, x-方向の和は同じ断面積で反対向きになり, x方向全体は打ち消す.

y,z方向も同じわけだから, 結局和は0というので証明出来たのではないか.

ところで, ある時気づいた, 多面体を水の中に置くという考えも棄て難い. ある面に立てた, 面積に比例する直交外向きのベクトルは, その面に働く水圧に相当する. 水中の多面体が流れないのは, その水圧の総和が0だからだと思うと, この証明は終りらしいが.

2011年8月12日金曜日

3での整除

TAOCPに, レジスタ上の64ビットの数が3で整除出来るかどうかの判定プログラムをMMIXで書けという問題7.1.3-215がある.

十進の場合, 各桁の和が3で整除出来るなら, 元の数も3で整除出来るのは知られている.

理由は1も10も100も..., つまり10nが3を法として1であるからだ.

二進の場合も同様に考えられる. 二進の1, 10, 100,...は, 交互に3を法として1と-1だから, ビットに右から0,1,2,...と番号をつければ, 偶数番目の1の数と奇数番目の1の数の差が3の倍数の時に3で整除出来る.

これは十進法で11で整除出来る判定と同じだ. 3003や1023や5060は11で整除出来る.

だから, (nu x)をxの横方向の和をとる関数, (band x y)をxとyとのビットごとのandをとる関数とし, 32ビットでお許し頂くと

(= (modulo (- (nu (band x #x55555555)) (nu (band x #xaaaaaaaa))) 3) 0)

ということになる. 最後にまた3での剰余をとるのは違反っぽいが, 十進の場合の各桁の和の3の整除と同じで, 大きい数をいきなり割るのよりは簡単であり, とりあえずはこういう考えでよかろう.

この辺でTAOCPの解答をみると, 一見妙なプログラムである. MMIXのアセンブリコードに慣れていない向きには申し訳ないが, コメントをつけてある.



最後の3で割った剰余をみるのに, 001001...を左シフトし, 符号ビットが1になるのを利用している. これはうまい.

ではその前は何をしたか. レジスタ0は偶数ビットと奇数ビットの横方向の和(これをr0とする). レジスタ1は偶数ビットの横方向の和(r1とする). ところで, 奇数ビットの横方向の和(r2とする)も欲しい. r2=r0-r1のはず. レジスタ0はその後32-r0になり, 2ADDUでレジスタ1の2倍に足される.

つまりレジスタ1はr1+r1+(32-r0)=32+r1-r2になり, #x249...をこれだけ左シフトすることになる. bの定数は32ビット目, 8文字目の次が#x9なので, r1=r2のとき符号ビットが1, 負数になり, 3で整除出来ることが分かる.

巧妙なプログラムであった. 斯くの如きプログラムの解読は楽しい.

2011年8月11日木曜日

多面体描画道楽

立方体は正方形6枚で出来ている. 正方形は1つの内角が90度. 従って正方形全体で内角の和は360度. それがF=6枚あるから, 2160度, つまり12π(パイ)になる.



一方頂点の数Vは8. それに対して内角の総和は(2V-4)π=12πという関係をEulerは知っていたらしい.

別の例. 正十二面体は正五角形がF=12枚. 頂点はHamilton世界巡礼でお馴染のV=20. 正五角形の内角は108度. 従って内角の総和は, 108*5*12=6480. 6480/180=36(π). 2V-4も36だ.

さらに2種類の面のある, 切頭二十面体(サッカーボール)は, 正五角形が12, 正六角形が20, 頂点は60. 内角の総和は, 108*5*12+120*6*20=20880, 20880/180=116. 2*60-4=116となる.

こう考えてみた.

Eulerの式に, 面の数F, 辺の数E, 頂点の数Vとすると, F+V=E+2という有名なのがある.

いま, ある多面体がn1角形, n2角形,...,nF角形のFの面で出来ていたとする.

n角形の内角の和は(n-2)πであったから, 計算したい内角の総和Sは
S=(n1-2)π+(n2-2)π+...+(nF-2)π=(n1+n2+...+nF)π-2Fπ.

一方, この多面体の辺の数を計算すると, 各多角形の辺は2回ずつ寄与しているから, E=(n1+n2+...+nF)/2.

従って, n1+n2+...+nF=2E. これをSの式に入れると,

S=2Eπ-2Fπ. Eulerの式から, E-F=V-2.

ゆえにS=2(V-2)π.

なるほど.

また, Descartesの定理というのもある.



左下は, 正四面体の頂点は, 正三角形が3個集まっていて, その内角の和が, 360度より足りない度合い(deficiency)が180度であることを示す. 正四面体には頂点が4個あるから,deficiencyの総和は720度. 4πである.

これを他の正多面体でやって見ると, やはり, どれも720度になるのである.

サッカーボールではどうか.

この図のように, 120度, 120度, 108度に挟まれた角は12度. それが60個あるから, やはり720度である.




Eulerの方は, サッカーボールの図でいうと, 348度の方を総和している. 60倍すれば先ほどの20880度になる.

V個のradianで表わしたdeficiencyをd1, d2,...,dVとする.

Eulerの方の内角の和は, (2π-d1)+(2π-d2)+...(2π-dV)
=2πV-(d1+d2+...+dV) と足す.

これが2(V-2)πだったから, d1+d2+...+dV=2πV-2(V-2)π=4π.

これもなるほどであった.

正多面体では, 面の数の順に正四面体, 立方体(正六面体), 正八面体, 正十二面体, 正二十面体と並べる. 正四面体の頂点は尖っているから, あとに行くほど丸くなるかと思うのは誤解であって, 正八面体も正二十面体も結構尖っている. それは正八面体は立方体より, 正二十面体は正十二面体よりも, deficiencyが大きく, 頂点の数が少ないためであろう.

サッカーボールがよく転がるのは, deficiencyが12度で, V=60もあるからである.