2016年5月1日日曜日

Prefix adder

60年くらい前, 東大物理学科高橋研究室でパラメトロン計算機PC-1を作っていたころ, 36ビットの並列加算器で, 高速に繰上げ(carry)したいという問題があった.

後藤さんはパラメトロンの多数決素子を巧妙に使い, 桁数に対して対数時間で繰上げを検出する回路を考案して, それを取り入れた.

先日出版されたDigital Design and Computer Architecture Arm Editionを眺めていたら, 同じような考えによる, Prefix Adderという回路の話があった. 今回は それが話題である.





この上の図が16ビットの二進数Ai, Bi (i=0,1,...,15)と下から来る繰上げCinを足すときのCarry-Lookaheadの回路だ.

上の図の上段の白い四角, 中頃の黒い四角, 下段の角のとれた白い四角の説明が下の図である.

ある桁 i が, 下の桁 i-1 からの繰上げの有無に関係なく, 上の桁 i+1 へ繰上げを生じるとき, 繰上げをgenerateするということで, Gi=1, そうでないとき, Gi=0とする. AiとBiが共に1なら繰上げがあるから, Gi=AiBi.

ある桁 i が, 下の桁 i-1 からの繰上げを上の桁 i+1 へ伝えるとき, 繰上げをpropagateするということで, Pi=1, そうでないとき, Pi=0とする. AiとBiのいずれかが1なら繰上げがあるから, Pi=Ai+Bi.

そうすると, i からi+1へ出る繰上げCi

Ci=Gi+PiCi-1

になる. このGやPをPrefixという.

これは単一の桁についての話であったが, iからjの桁のブロックについても同じことがいえる. ただGi, Piの代わりにGi:j, Pi:jのように書く.

そうだと思って上の図を見てみると, 左下のあたりに黒四角で左下隅に14:-1と書いたものがある. つまり15の桁への下全体のブロックからの繰上げはG14:-1 (と, P14:-1 と, Cin) で決るということを示す.

GiとPiがAiとBiとで書けたように, Gi:jとPi:jがそのブロックを途中で分けたGi:kとGk+1:j, Pi:kとPk+1:jで同じように書くと, Gi:jはi,k間で繰上げが出るか, k+1,j間で繰上げが出て, i,k間で繰上げが伝わるかだから, Gi:j=Gi:k+Pi,kGk-1,j.

Pi:jはi,k間でもk-1,j間でも繰上げが伝わるかだから, Pi:j=Pi:kPk-1:j.

ここまではDigital Designの本に書いてあることだ. 上の図のように半分, 半分と区間を分るにはi,jからkをどう決めればよいかは, この本には書いてない. そこがこの本とThe Art of Computer Programmingの違うところである. TAOCPならk(i,j)が詳しく記述してある筈だ.

しかしそう難しいことではなくて, λx(x>0)をxを二進法で表したとき, 最も左にある1のビット番号とすると, (λ1=0, λ2=1,λ3=1, λ4=2,...)

k=2λ(i-j)+j

と簡単に得れれる. Schemeでシミュレートしてみたが, これで合っているようであった.

2016年3月1日火曜日

菱形六十面体

菱形六十面体についてはこれまで何度もこのブログに書いた. 夫々の菱形が黄金比であることは, MathWorldにそう書いてあったのを鵜呑みにしたきらいがある.

今回はそれを確かめてみたい.

下の図は正十二面体である. 中心を通る赤, 緑, 青の線はx, y, z軸を示す. それぞれの軸は正十二面体の30本ある稜線のうち6本の中央を通過する. 次の図との関係で偶数だけがついているが, その頂点を持つ正五角形で座標を考える.


頂点間の距離, 稜線0, 2の長さを2とすると, 正五角形の対角線8, 4の長さは2φ. 従って4の頂点の座標は(φ, φ, φ). それを基準に他の座標も分かって, 1は(φ2, -1, 0);2は(φ2, 1, 0); 6は(1, 0, φ2); 8は(φ, -φ, φ)である. (1のx座標は中心から0までの距離と, 中心から4までの距離が等しいことから計算できる.) この正五角形に5個の菱形を作ったのが下の図の菱形六十面体だ.



これをx軸の方向から撮影した写真がこれである.



写真では視差があるから, 注目している正五角形のついて正確な図にしたのが次の図だ.



ちょうど1のところをx軸がこちらへ向かっている. その両側へ1離れたところに02が来る. また1から真上に行ったところに中心10があるわけだが, その高さは現時点では未定.

さて菱形は真上から見れば菱形であるが, 少し横から見れば菱形に見える方向も皆無ではないものの, 一般には菱形に見えない. しかし平行四辺形であるのは確かである.

従って0,1,10,9の菱形, 1,2,3,10の菱形を考えると, 0,91,3の辺は1,10の辺と平行である. また9,100,1と平行なども分かる.

すると8,75,40,1などに平行である. 一方 10,5, 6,7の菱形は対角線の方向から見ているから菱形のままであり, 10,55,6の長さは等しい. 5, 6の高さの差(φ2-φ=1)は既知だから, 10の高さも確定し, 計算するとφ-1=1/φである. 5,4の幅が1だから, 5のy座標も1/φと決まる.

ところで残りは奥行きである. 2, 3 の辺は, 横から見ると, 正面図の6, 5の辺に相当していることが分かる. 従って3の奥行きが分かり, 3, 4 がz軸に平行なことも分かる. 従って10, 5, さらに5, 6もz軸に平行だったことも分かった.

つまり正面から見た菱形10,5, 6,7は真の形であった.

という次第で菱形の縦横比は1:φであったわけである.

2016年2月24日水曜日

MathematicaのCal

いろいろなプログラム言語でCalのプログラムを書いてみた. 今回はMathematicaである.

Mathematicaというと, 一発処理のような印象である. 例えば Prime[25] と入力すると, 25番目の素数 97が出力される. しかしなんとかプログラムを書いてみようと思い, 例によってCalを書いてみた. 結局いろいろなノウハウを覚えつつ試行錯誤の結果完成したのが次のプログラムである. その見どころを説明しよう.



0行目. 仮引数はyr_のように後に下線をつける. 本体で使うときはその下線はいらない. システムの関数手続きは大文字で始めるが, ユーザの関数は小文字で始める. 本体の定義は:=の後に書く.

1行目. Module[{ の次にローカル変数を書く. ローカルな関数名も書く.

2,3行目. printd[d_]は, その月の最後の日付をzとして, 1≤d≤zの時は出力し, それ以外は3文字の空白にする. さらにdが1桁なら前に2文字の空白, 2桁なら1文字の空白をつける.

4,5行目. prline[d_]は, dからd+6までをprintdで文字列にしたものを, StringJoin(<>を使う)し, Printで出力する.

6,7行目. 年と月を出力. また曜日名を出力.

8,9行目. うるう年なら, ローカル変数lを1に, そうでなければ0にする.

10行目. 各月の朔日の曜日を決めるため, 次の年の初めまでの日数を計算.

11行目. 1,2月の日数と, 月初めから年末までの日数(mod 7)を設定.

12,13行目. 3月以降の月に上のような値を計算する.

14行目. 朔日の曜日wを確定する.

15行目. dを第1週の日曜に相当する日(1-w)から月末日まで, 7ずつ増やしてprlineする.

こうして出来たのが次だ.



なれない言語でプログラムを書くのはしんどい.

2016年2月15日月曜日

月齢カレンダー

このブログの2008年12月に月齢計算のことを書いた. それは計算機を使って月齢を計算する自己流アルゴリズムであったが, その昔, 私は高校生の頃, 先輩に教えられて毎年月齢定数というのを使って月齢を概算していた.

それではその年の月齢定数をCとすると, その年のmd日の近似的な月齢を (C + m + d) mod 30 で得るものであった.

それで時々年の初めころ, 今年の月齢定数はいくつかなとほどほどに計算して記憶していたこともあった.

ところでインターネットであちこち見ていたら堀さんの式というのに出会った. http://www.asj.or.jp/geppou/archive_open/1968/pdf/19680704.pdf

詳しくは元記事を見ていただくことにして, 大体は月齢定数と同様の趣旨であるが, mの代わりに f(m)を足すのである. この値が1月から順に0,2,0,2,2,4,5,6,7,8,9,10なのである. 6月から後はf(m) = m - 2なので, その辺は昔の月齢定数と同様だが, 年初がすこし違っていて, この式の方がより正確になるらしい. 元記事のタイトル「鬼鬼西」は1月から6月までの値の覚え方である.

Cの方は, (((y - 11) mod 19) × 11) mod 30 だそうで, 本年2016年なら 2005 mod 19 = 10なので, 110 mod 30 = 20が今年の月齢定数になる. 元記事の掲載された1968年はこの定数がちょうど0だったのでタイムリーな発表であった.

f(2)=2だから (20+2+8) mod 30 = 0 となり, 2月8日が新月つまり春節であった. 3月の満月をxとすると, 20 + 0 + x = 15 (mod 30) だからx = 25となって, 今年の復活祭は3月27日である. また3月9日には部分日蝕があるが, その月齢は 20 + 0 + 9 = 29 だ.

なおこの定数は 2005 mod 19 = 10 から後8年はこの値が1ずつ増え, 従って定数も毎年11ずつ増えるから初めから計算しなおすことはなく, 来年2017年は1, 2018年は12, 2019年は23,...と分かる. こうして..., 7, 18と来たら, その次は29ではなく0とする(なる?)のだ. 次回0なのは2025年である. 1968年, 1987年, 2006年, 2025年と19年ごとに定数が0の年が回って来る.

2016年2月14日日曜日

菱形六十面体

2010年5月17日のこのブログの続きである. そこにはWolfram MathWorldのRhombic Hexecontahedronの項に 20 golden rhombohedra can be combined to form a solid rhombic hexecontahedron. と書いてあるのが気になっていて, そのrhombohedronの頂点の立体角を計算し, たしかに20個あると全立体角, 4πステラジアンになると書いた.

golden rhombohedra は菱形六面体とでも訳そうか, 6つの面が黄金菱形になっているものである.

Wolfram MathWorldのGolden Rhombohedronの項は, これには鋭角のものと鈍角のものがあるというが, ここで使うのは鋭角の方である. この項には1辺の長さaの鋭角黄金六面体の 体積は1/5√(10-2√5) a3とあり, また菱形六十面体の項にはその体積が4√(2(5+√5))a3とあるから, 体積的にも確かにそうである.

ざっと考えてみると, 菱形六十面体の母体である正十二面体には頂点が20あるから, 多分その頂点の数と菱形六面体が対応しているであろうと推察できるが, あまり自信はない.

昨今, 3Dプリンタが使えるようになったので, ちょっと作ってみようという気になった.



上の写真の上は20個の菱形六面体が出来たところ, 下はその5個を使い, 菱形六十面体の内, 正十二面体の一面に相当する面が完成した所である.

この写真では分らないが, 3Dプリンタで作ったものは, 存外 外面がざらざらで, 接着するのが大変であった. 次回はもっと時間がかかっても, もう少し精密な面にしてみたい.

さて, 次の5個で同様に正十二面体の反対側の面に相当する面も作る. この両方の5個組を北極と南極とすれば, 残りの10個で赤道に相当する部分を構成することになる.

それには10個の内2個ずつを接着して2個組を5組作り, それを一周するように繋げるのである. これは仲々面倒で, どの面とどの面が合さるのか, いろいろこねくり回さなければならない. まぁそうこうする内に固定出来た. 下の写真は接着中の写真で, 小型のクランプやしゃこ万力で押えているところである. 左に転がっている2個は北極と南極の5個組である.



出来上がった赤道部分のこのようである.



次の図はPostScriptで描いた, 菱形六十面体と上でいう北極, 赤道, 南極の図である. 手前右上から光が当っているように影を付けてみたが, 効果はいまいちであった.

全体
北極
赤道
南極

相対的な位置が分るように, 外接する正十二面体と, 中心を通るx,y,z軸も描いてある.

これで見ると20個すべての菱形六面体が中心に接しているのが分る. つまり赤道の組は中央が厚さ0である.

3Dプリンタで作ったモデルにはこのような精度は残念ながらない. でもこんな実験が出来るのも3Dプリンタのお蔭である.

2016年2月3日水曜日

菱形六十面体

菱形六十面体についてこのブログで何度も取り上げた. その最初, 2010年5月16日のに「角D'B'E'は黄金菱形なので, 黄金比φ=1.61に対して2cot-1(φ)=63.435°である. その一端を持ち上げると, 角DBEは少しずつ広がり, 72°になるわけだ. その時の傾き角を計算すると, 驚いたことに31.7175°であった. 」と書いている.


長い間 気になっていたが, 先頃 証明した. 図を見てほしい.

図の下は破線のAD'B'E'が縦1, 横φの黄金比の菱形である. 実線のADBEは角Aが72度で縦が1の菱形である. この横の対角線の長さを以下eとする. 黄金比の方の角Aの半分は31.7175°である.

図の上の直角三角形は底辺が下の実線の菱形の横の対角線, 斜線が黄金比の対角線で, これを描いてみると角Aがやはり31.7175°になるのである. この角を以下θとする.

証明したいのは, 上の31.7175°, つまりcos θ=e/φが下の31.7175°, つまりtan 2/φに等しいということである. 一連の計算は下の通り.


φの値はよく知られている. tan 36°はcos 36°がφ/2 から計算できる.

(0)上の図のθが欲しいので, atanの式を書く. 以下atanの引数を計算する.

(1) φとeの値を代入する. (2)根号の中, 第1項を展開し, 第2項の分母を有理化する. (3)根号の中を通分し, 1/eの値を書く.

(4)根号の中の分子を計算する.

(5)根号の中, 分母分子を2倍すると, 分母分子とも開平出来た.

一方 1/φも計算すると(6)(7)同じ値になった. めでたしめでたし.

2016年1月25日月曜日

曜日の計算

久し振りに曜日の計算の話だ. 2011年5月25日の私のブログ(曜日の計算)に

「21世紀中の 2000+y年m月d日の曜日は, (y+floor(y/4)+h(m)+d+4) mod 7 を計算するとそれが曜日になる.」

ただし「うるう年の1月と2月はh(m)の値から1を引く必要がある.」

とあり,
m 1月 2月 3月 4月 5月 6月 7月 8月 9月10月11月12月
h   2   5   5   8   3   6   1   4   7   2   5   0
の表がある.

今日の話題は上の式の「+4」と直ぐ上のh(m)である. 21世紀は+4だと書いたが, では他の世紀ではどうなるか. 20世紀に当て嵌めて計算してみると, その時代には+5であった.

我々の生きている20, 21世紀中は, この+5と+4を覚えておけばいいが, もう少し一般的に出来ないかと考えてみた. グレゴリオ暦は4世紀で曜日が元に戻るから, 年を100で割った整数部をさらに4で割った剰余に対してこの値を決めればよいわけだ.

1900年からは(剰余が3で)+5
2000年からは(剰余が0で)+4
もっとやってみると
2100年からは(剰余が1で)+2
2200年からは(剰余が2で)+0
だから, この剰余の0,1,2,3に対応して +4,+2, 0,-2と覚えておけばよいことが分った.

剰余  0  1  2  3
b    +4 +2  0 -2
上の式は (y+floor(y/4)+h(m)+d+b(c)) mod 7 になる.

h(m)の表は, 1月の2, 2月の5は覚えるとして, 後は奇数の月も偶数の月も5回ずつあり, 偶数月は12-m, 奇数月は3,5,7月は8-m, 9,11月は16-mと覚えるのがよさそうだとこの頃は思っている.

使ってみよう. グレゴリオ暦に改暦された1582年10月15日の曜日は (82+20+2+15-2) mod 7=117 mod 7=5 だから金曜. 日本でグレゴリオ暦に改暦された1873年1月1日の曜日は (73+18+2+1+0) mod 7=94 mod 7=3 だから水曜であった.