2011年12月27日火曜日

中央値関数の排他的論理和

TAOCPには面白い話が満載である. 今回は多数決関数関連の話だ. この本では, 多数決関数のことを中央値(median)関数というので, そう書くことにしよう.

中央値関数は, 奇数個の引数を大きさの順に並べ, その中央の値を返す(平均値でないことに注意).

    x≤y≤zの時 <xyz>=y.

あるグループで流行っている最中限ゲームみたいなものである.

論理関数としては, x,yを偽か真とする時, 常に反対するか常に賛成する偽か真のサクラの入力を加え, <偽xy>がxとyの論理積, <真xy>が論理和である.

実数値については, <-∞xy>はxとyの最小値, <+∞xy>は最大値になる.

こういう話題がいちおう終わってから,

m>0について,
x0⊕x1⊕...⊕x2m=<¬x0s1s1...s2m>
ここで, sj=<x0xjxj+1...xj+m-1¬xj+m¬xj+m+1...¬xj+2m-1>,
ただし, k≥1について, x2m+k=xkとする;

と書いてあった.

つまり, 何変数であっても, その排他的論理和がサクラ(定数)なしの2段の中央値関数で得られるのである. 美しい.

そもそも排他的論理和は単調関数ではないから, 中央値関数1段だけでは得られない. 3変数の排他的論理和, つまり全加算器もパラメトロンでは2段必要であった.

上の見事な式は, TAOCPによると, 量子化学が専門の佐々木不可止(ふかし)君が考えたという.

佐々木君は東大物理の学生時代に, 高橋研によく出入りしていて, 多数決論理に関心があったのは知っていたが, こういうことも考えたのかと驚いた. 北海道大学の化学科の教授を定年退職し, 苫小牧駒沢大学に移られたが, そこも退職の後, 惜しいことに2007年に他界された.

あの長々とした式は眼がくらくらするので, TAOCPにあったようにm=1とm=2について書き直すと

x0⊕x1⊕x2=<¬x0<x0x1¬x2><x0x2¬x1>>

x0⊕...⊕x4=<¬x0<x0x1x2¬x3¬x4><x0x2x3¬x4¬x1><x0x3x4¬x1¬x2><x0x4x1¬x2¬x3>>

となる.

この式をテストしてみよう. ます多数決関数mを定義する. nは否定をとる関数.

(define (m . xs)
(define (n x) (- 1 x))
(let ((l (length xs)))
(if (> (apply + xs) (/ (- l 1) 2)) 1 0)))
(m 1 0 0) => 0
(m 1 0 1) => 1
(m (n 1)(n 1) 1) => 0

すると5変数(m=2)の排他的論理和はこう書ける.

(define (ex x0 x1 x2 x3 x4)
(let ((s1 (m x0 x1 x2 (n x3) (n x4)))
(s2 (m x0 x2 x3 (n x4) (n x1)))
(s3 (m x0 x3 x4 (n x1) (n x2)))
(s4 (m x0 x4 x1 (n x2) (n x3))))
(m (n x0) s1 s2 s3 s4)))
(map (lambda (x) (apply ex (n2bs x 5))) (a2b 0 32)) =>
(0 1 1 0 1 0 0 1 1 0 0 1 0 1 1 0 1 0 0 1 0 1 1 0 0 1 1 0 1 0 0 1)

でどうやらよさそうだ.

パラメトロンは通常3入力止りだが, 5入力が許されるとしてこの回路はこう描ける.



これをさらに納得すべく, s1,...s4や, その計算の元になっているx1+x2-x3-x4
などを, x0=0である前半の16個について出力してみよう. 自己相対関数なので, 半分やればいいことになっている.

xの否定は1-xと計算し, (l-1)/2の2と比較したから, 閾値と較べる値は, 本来は
    x1+x2+(1-x3)+(1-x4)=2+x1+x2-x3-x4.

従って, 先頭の2を外して,
    x1+x2-x3-x4
の和sumについてs1はsum≤0なら偽(x0=0を除外したから, 等号が要る), sum>0なら真である.

s1と添字がmだけ違うs3は, 足す項と引く項がちょうど反対だから, 和も正負逆で同じ値になる. こういう対を相棒といおう. ところで, 真が奇数個だと, 和が0になることはないので, 相棒同士の真理値は反対になり, x0がキャスティングボートをとって, exは真になる. すなわち排他的論理和である.

ところが, 真が偶数個だと, どこかに相棒同士で和が0になるものがあり, x0は0だから, sは2つとも偽になる. 先頭の¬x0が1を主張しても, 和が0でない他の相棒同士は打ち消し合い, この2個(あるいは4個)の0連合には勝てず, exは0になる.


+ s1 + s2 + s3 + s4 ex
0( 0 0 0 0 0 0 0 0 0)
1(-1 0 -1 0 1 1 1 1 1)
2(-1 0 1 1 1 1 -1 0 1)
3(-2 0 0 0 2 1 0 0 0)
4( 1 1 1 1 -1 0 -1 0 1)
5( 0 0 0 0 0 0 0 0 0)
6( 0 0 2 1 0 0 -2 0 0)
7(-1 0 1 1 1 1 -1 0 1)
8( 1 1 -1 0 -1 0 1 1 1)
9( 0 0 -2 0 0 0 2 1 0)
a( 0 0 0 0 0 0 0 0 0)
b(-1 0 -1 0 1 1 1 1 1)
c( 2 1 0 0 -2 0 0 0 0)
d( 1 1 -1 0 -1 0 1 1 1)
e( 1 1 1 1 -1 0 -1 0 1)
f( 0 0 0 0 0 0 0 0 0)

左端は行数である. 各行のかっこ内の9個の値は, 左からx1+x2-x3-x4, s1, x2+x3-x4-x1, s2, x3+x4-x1-x2, s3, x4+x1-x2-x3, s4, exである. 偶数個の場合, どの対で和が0になるかを調べたのが下だ. 赤の-は, x1とx2の和がx3とx4の和に等しい(つまりs1s3が0になる)場合, 青の-は, x2とx3の和がx4とx1の和に等しい(つまりs2s4が0になる)場合である.


0(0 0-0-0-0-)
1(0 0 0 0 1 )
2(0 0 0 1 0 )
3(0 0 0-1 1-)
4(0 0 1 0 0 )
5(0 0-1-0-1-)
6(0 0-1 1-0 )
7(0 0 1 1 1 )
8(0 1 0 0 0 )
9(0 1-0 0-1 )
a(0 1-0-1-0-)
b(0 1 0 1 1 )
c(0 1 1-0 0-)
d(0 1 1 0 1 )
e(0 1 1 1 0 )
f(0 1-1-1-1-)


和が0になる対は本当にあるか. こういう図を描いてみた.



いま真というか1の数は偶数である. s1とsm+1の範囲で1の数が同じであれば, 一発で和が同じところは見つかったわけでハッピーだ.

一方, s1がオール0でsm+1がオール1なら, 左の密度は0, 右の密度は1. そこで図の下へs2, s3のように加算の範囲をずらしていけば, どこかで2つの範囲の1の数が一致する, 同じ密度になる場所があり, そこで赤や青の-が引ける対が見つかるはずである(平均値の定理を習ったのは高校生の時だっけ).

こういうわけで, 偶数の場合はexが偽になるのである.

「なるほど」であった.

実はサクラを使うと, どんな関数でも2段の中央値関数で作れるが, そのことはいずれ.

2011年12月4日日曜日

平方根の連分数展開

Norman Lehmerの論文に, factor stencilはNの連分数展開でも使えると書いてあった.

その始めの方. まずx1=√N=q1+1/x2とおく. 但し, q1はNの平方根の最大整数. 同様にx2=(√N+P2)/Q2=q2+1/x3のようにする. (ははぁ. TAOCPの連分数展開のUとVがPとQになっているわけだ.)

順に進んで, xn=(√N+Pn)/Qn=qn+1/xn+1になるというまでは明快だが, その後に
「従って
Pn=qn-1Qn-1-Pn-1,
Qn=qn-1(Pn-1-Pn)+Qn-2
がexpediouslyに得られる」
と書いてある.

上の方の式は11月30日のブログの(4)と(5)の間にあるUn=AnVn-Un-1と同じだが, 下の方の式はどう導くか.

今回のブログはその計算法である.

(0)はxnの標準形である. (1)はxn-1で, そのxnに(0)を代入し, 前回のように計算を続けると(1)の最後になる. (1)の最初と分母を等しいとしたのが(2)で, ただちに(3)と書ける. 分子の√Nの後に続く項を等しいとしたのが(4)で, 途中で(3)を利用する. PnとPn-1を交換したのが上の式(5)だ.

一方, (3)の添字を1減らしたのが(6)で, (3)-(6)が(7),(8)である. 1/Qn-1は, (5)を移項した(9)から(10)のように得られ, これを(8)に代入すると(11)を経て, 下の方の式(12)になるのが分かる.

漸化式はP1=0, Q1=1から始まるが, 問題は, Qn-2があるので, Q0が必要なことだ.

実はQ0=Nなのだが, この計算は次の通り.

これで, 連分数のP,Q,qが次々と求まることになる. なんだかまた久しぶりにこのような計算をした.

2011年12月3日土曜日

素因数探し

11月6日の本ブログにある自転車チェーンの篩を作ったのはHenry Lehmerだ. その父親のNorman LehmerもCalifornia大学Berkeley校の数学の教授であった. 父のLehmerが素因数探しの道具として, Factor Stencilを考案したという話がある.

それがどういうものかが分かってきた. 今回はその話をしよう.

aが素数pの平方剰余であるとは, x2=a mod pとなる何かのxがあることなのは周知のとおり. こういうa(≠0)の時, Legendreの記法で(a/p)=+1とする. (htmlなので, 分子分母が横並びだが, 通常は上下に書く.)

xがないときは, aを平方非剰余といい, (a/p)=-1とする.

p=7の時, 剰余は(0,1,2,4), 非剰余は(3,5,6)である. ここで, 0を除外し, 非剰余同士を掛けると, 3*3=2 mod 7, 3*5=1 mod 7, 5*5=4 mod 7のように剰余になり, 剰余と非剰余を掛けると2*3=6 mod 7のように非剰余のまま, 剰余同士を掛けると4*4=2 mod 7で剰余であり, +1と-1で掛算になる.

ところで, ある整数Nを素因数分解する場合, 任意の整数aについて, a2-N=Rは, Nがpで割れるなら, mod pで考えると, Nの項はないのと同じだから, Rはこのpの平方剰余になる.

aを√Nの程度にとると, Rは小さくなり, 平方剰余の乗算に規則でもっと分解すると, さらに小さい剰余か非剰余が得られ, それから素因数の形が決るという.

たとえば, 平方剰余に-1があるなら, すなわちp-1なら, 素因数は4n+1, 2があるなら, 素因数は8n±1の形, というふうに, 素因数はlinear formになっている. それを組み合せて, 素因数の探索の範囲を狭めたい. ところが, 平方剰余が段々大きくなってくると, linear formの形も複雑になり, 式の形を決めにくくなる. Normanの論文によると, Rが113と199の場合, formの形は11088にもなるそうだ.

そこで, 各平方剰余について, 候補になる素数の表をあらかじめ作っておこうというのが, この発想である.

小さい範囲で, テストしてみよう. 下がその剰余(縦方向)と素数(横方向)の表である. 10月26日のブログの図の素数のところを抜き出したものと思えばよい.

3 5 7111317192329313741434753596167717379838997
-1 X X X X X X X X X X X
2 X X X X X X X X X X X
3 X X X X X X X X X X X
5 X X X X X X X X X X
7 X X X X X X X X X
11 X X X X X X X X X X
13 X X X X X X X X
17 X X X X X X X X X
19 X X X X X X X X X X

Normanの素数表には1があるというので, この左に1と2の列があるかも知れないが, 除数としては不要なので, とりあえず3からの表だ. 7の下を見ると, 2と11にXがあり, 2と11=4 mod 7が平方剰余であることを示す. 他の素数についても同じだ. たとえば,

(quadres 97)=>(0 1 2 3 4 6 8 9 11 12 16 18 22 24 25 27 31 32
33 35 36 43 44 47 48 49 50 53 54 61 62 64 65 66 70 72 73 75
79 81 85 86 88 89 91 93 94 95 96)

のうち, 19までの素数, 2,3,11と, 96(= -1 mod 97)があるから-1とにXがある.

また, 横向きに眺めると, -1の行は, 5,13,17,29,37,41,53,61,73,89,97で, すべて4n+1. 2の行の7,23,31,47,71,79は8n-1, 17,41,73,89,97は8n+1である.

さて, これを使い, N=2279(=43x53) の素因数を探してみる. (sqrt N) => 47.738873049120045 だから, aを48から順に増やしながら, a2-Nを計算し, その素因数分解もしておく.

(map (lambda (a) (factorize (- (* a a) n))) (a2b 48 80))
=> ((5 5) (61 2) (17 13) (23 7 2) (17 5 5) (53 5 2) (13 7 7)
(373 2) (857) (97 5 2) (31 7 5) (601 2) (1321) (103 7 2)
(313 5) (13 13 5 2) (79 23) (139 7 2) (67 31) (17 13 5 2)
(67 7 5) (73 17 2) (2621) (1381 2) (83 7 5) (61 5 5 2)
(139 23) (239 7 2) (269 13) (73 5 5 2) (761 5) (283 7 2))

ここに並んでいるのは, 掛けると+1になるもの同士である. (5 5)からは, 5が+1かも, -1かも知れないということが分かるだけだ. (61 2)では, 61が上の表の剰余にないから, 無視. (17 13)は, とりあえず13と17は+1と+1か, -1と-1が分かっただけ. その次の(17 5 5)はラッキーだ. 5と5は打ち消すから, 17が+1と判明した. したがって, (17 13)から, 13も+1と決る.

(13 7 7)からも13が剰余と確信出来る. さらに進むと, (13 13 5 2)や(17 13 5 2)から, 2と5は+1同士か-1同士か, いずれにしろ同じ仲間ということが分かるが, 役にたつかどうかは不明.

しかし, 13と17が得られたので, 13と17の行にXのある素数を探すと, ブラボー! 43と53が見つかった. ついでだが, 43と53の列では, 2と5は空白で, 同じ組なことが示せる. いまは1組だけが見つかったが, 通常は候補がたくさん得られ, 実際に割ってみる必要がある.

Norman Lehmerが1914年に出版したFactor Stencilは, 剰余が±238まで, 素数が48593までの素数表であって, 各剰余ごとに1枚で, Xの代りに素数の位置に孔をあけてある. 今の例では13と17の紙を取り出してきて重ねると, 条件に合う素数のところから光が洩れる仕掛けであった. 当時としては大変な労作であった. いまでは, MacBookで簡単に実験が出来, ありがたい時代だ.