単純な素因数分解アルゴリズムを実装してみる
「OpenSSLとPythonでRSA暗号の原理を知る」で書いたように、RSA暗号は公開鍵である整数modulusを二つの素数に素因数分解することができれば破ることができる。 ここでは、いくつかの単純な素因数分解アルゴリズムをPythonで実装し、一般的なPCで解くことができるbit数などを調べてみる。 なお、正確には素因数分解ではなく、因数をひとつ見つけるアルゴリズムについて考えるものとする。
環境
Ubuntu 14.04.3 LTS 64bit版、Intel Core i5-4200U (1.6GHz * 2 * 2)
$ uname -a Linux vm-ubuntu64 3.19.0-25-generic #26~14.04.1-Ubuntu SMP Fri Jul 24 21:16:20 UTC 2015 x86_64 x86_64 x86_64 GNU/Linux $ lsb_release -a No LSB modules are available. Distributor ID: Ubuntu Description: Ubuntu 14.04.3 LTS Release: 14.04 Codename: trusty $ grep "processor\|model name" /proc/cpuinfo processor : 0 model name : Intel(R) Core(TM) i5-4200U CPU @ 1.60GHz processor : 1 model name : Intel(R) Core(TM) i5-4200U CPU @ 1.60GHz processor : 2 model name : Intel(R) Core(TM) i5-4200U CPU @ 1.60GHz processor : 3 model name : Intel(R) Core(TM) i5-4200U CPU @ 1.60GHz
素数の生成
素数とその積の生成は、OpenSSLを用いて行うこととする。 128bitのmodulusを生成する場合は次のようになる。
$ openssl genrsa 128 | openssl rsa -noout -text
Generating RSA private key, 128 bit long modulus
.+++++++++++++++++++++++++++
.+++++++++++++++++++++++++++
e is 65537 (0x10001)
Private-Key: (128 bit)
modulus:
00:bf:27:4e:d6:73:f6:ca:f4:f6:85:f0:b6:a0:27:
1a:31
publicExponent: 65537 (0x10001)
privateExponent:
76:67:8f:d3:10:e3:d7:14:d1:a8:40:83:b3:81:1b:
ad
prime1: 16144814945699257943 (0xe00de7a77ad25e57)
prime2: 15737972014275192503 (0xda688229e3e966b7)
exponent1: 9570317772029090007 (0x84d0994ed566d8d7)
exponent2: 9046986552112602077 (0x7d8d5a560e685bdd)
coefficient: 9750135051940033389 (0x874f702bdaf8736d)
試し割り
もっとも単純なアルゴリズムであり、整数nに対してsqrt(n)以下の整数で順に割っていくというもの。 計算量はO(n1/2)である。
# trial_division.py
import sys
def isqrt(n):
x = n
y = (x + 1) // 2
while y < x:
x = y
y = (x + n//x) // 2
return x
def trial_division(n):
k = isqrt(n)
if n % 2 == 0:
return 2
if n % 3 == 0:
return 3
i = 6
while i < k:
if n % (i-1) == 0:
return i-1
if n % (i+1) == 0:
return i+1
i += 6
if __name__ == '__main__':
n = long(sys.argv[1])
p = trial_division(n)
if p:
print "%d = %d * %d" % (n, p, n/p)
else:
print "%d is prime" % n
上のコードでは、2と3の倍数となるものについては飛ばすよう若干の効率化を加えている。 また、isqrt関数は整数の範囲で平方根を求める関数である。
次の実行結果の通り、64bitの素因数分解に5分程度かかる。
$ time python trial_division.py 12814570762777948741 12814570762777948741 = 3318288047 * 3861801803 real 4m47.971s user 4m47.296s sys 0m0.596s $ time python trial_division.py 18366865165381711817 18366865165381711817 is prime real 6m15.418s user 6m15.348s sys 0m0.000s
Pollard-Rhoアルゴリズム + Miller-Rabin素数判定法
Pollard-Rhoアルゴリズムは、乱択アルゴリズムにより高速に因数を発見するアルゴリズムである。 整数nに対しその因数をpとしたとき、計算量はO(p1/2) <= O(n1/4)となる。
乱択アルゴリズムという性質上、与えられた数が素数だった場合、停止するまでに非常に時間がかかる。 そのため、同じく確率的アルゴリズムであるMiller-Rabin素数判定法と組み合わせ、あらかじめ素数かどうかを判定しておくことが多い。
# pollard_rho.py
import sys
from random import randint
from fractions import gcd
def miller_rabin(n, k=20):
s, d = 0, n-1
while d % 2 == 0:
s += 1
d /= 2
for i in xrange(k):
a = randint(2, n-1)
x = pow(a, d, n)
if x == 1:
continue
for r in xrange(s):
if x == n-1:
break
x = (x*x) % n
else:
return False
return True
def pollard_rho(n):
x, y, d = 2, 2, 1
while d == 1:
x = (x*x + 1) % n
y = (y*y + 1) % n
y = (y*y + 1) % n
d = gcd(abs(x-y), n)
if d != n:
return d
if __name__ == '__main__':
n = long(sys.argv[1])
is_prime = miller_rabin(n)
if is_prime:
print "%d is prime" % n
else:
p = pollard_rho(n)
print "%d = %d * %d" % (n, p, n/p)
Miller-Rabin素数判定法において、素数判定を誤る確率は繰り返し回数kに対し4-kとなることが知られている。
上のコードではk = 20としている。
次の実行結果の通り、96bitの素因数分解に1分程度かかる。 なお、素数の場合についてはMiller-Rabin素数判定法により1秒未満で判定できる。
$ time python pollard_rho.py 60766145992321225002169406923 60766145992321225002169406923 = 250117558771727 * 242950340194949 real 1m12.204s user 1m12.192s sys 0m0.000s $ time python pollard_rho.py 71939287897297826407363026419 71939287897297826407363026419 is prime real 0m0.029s user 0m0.024s sys 0m0.004s
フェルマー法
フェルマー法はRSAのように素数同士の積からなる合成数において、二つの素数の差が小さい場合に有効なアルゴリズムである。
具体的には、合成数nの平方根周辺のある値をxとしたとき、(x-k)*(x+k) == nとなるkを小さいものから求める。
# fermat.py
import sys
def isqrt(n):
x = n
y = (x + 1) // 2
while y < x:
x = y
y = (x + n//x) // 2
return x
def is_square(n):
if not n % 48 in (0, 1, 4, 9, 16, 25, 33, 36):
return False
x = isqrt(n)
return x*x == n
def fermat(n):
a = isqrt(n)
b2 = a*a - n
while not is_square(b2):
a += 1
b2 = a*a - n
return a - isqrt(b2)
if __name__ == '__main__':
n = long(sys.argv[1])
p = fermat(n)
if p:
print "%d = %d * %d" % (n, p, n/p)
else:
print "%d is prime" % n
上のコードにおいて、is_square関数は値が平方数かどうかを判定する関数である。 平方数の判定については、平方剰余をもとに平方数となりえない値を先に判定することで効率化が図ることができ、上では法48について一度判定することで83.4%の剰余について計算を省いている。
次の実行結果の通り、64bitの素因数分解に30秒程度かかる。 素数の場合については大きな時間がかかるが、Pollard-Rhoアルゴリズムの場合と同様にMiller-Rabin素数判定法を組み合わせれば高速に判定することが可能である。
$ time python fermat.py 12814570762777948741 12814570762777948741 = 3318288047 * 3861801803 real 0m23.931s user 0m23.924s sys 0m0.000s
また、kが小さい場合についてはnの大きさによらず高速に素因数分解できる。 ただし、そのような素数の組が選ばれる確率は非常に小さい。
次の実行結果では、nが4096bit、kが15bitの場合について、1秒未満で素因数分解できている（参考）。
$ time python fermat.py 1034186074668005446865661651578208369514158803552244005581133348807218259696136975186768057532805396894757503164834709575083359161158060701923008831149233177215596305395749862300700772298386522636562203611668211601359078939366452933448401372783149865070405032946502692195960337152490393257718681775969357953766124434670452760182942622011519373645653657784274264774773270847089540432369489508909290697992927057404389218509541644042197448296417276686963995927414541514778710278820030649227976349765070908207584414253880096569270321293357683127769728901749940165024782079154839570419235490028108256539125090626571373643257455521820767805808214404175900736109428513289545954241551995058498737472560704736894014465695117165251699935296644992760434412206519036326801985540077086806228464411438175416317958133959586691011057808141234300883692932436022122663781160823472736298914304029803683697389880156222528045623045436226522780569221934366978973936408257173404970887883159180079205465195884955186605310665826226975454171035075016280417683011623124380312403051425940705920625660567584606796218038626946234225864706878239646554381165540673005572649955471689205027195933802181344480474070023069194083134722979170599617466579222550596280381011 1034186074668005446865661651578208369514158803552244005581133348807218259696136975186768057532805396894757503164834709575083359161158060701923008831149233177215596305395749862300700772298386522636562203611668211601359078939366452933448401372783149865070405032946502692195960337152490393257718681775969357953766124434670452760182942622011519373645653657784274264774773270847089540432369489508909290697992927057404389218509541644042197448296417276686963995927414541514778710278820030649227976349765070908207584414253880096569270321293357683127769728901749940165024782079154839570419235490028108256539125090626571373643257455521820767805808214404175900736109428513289545954241551995058498737472560704736894014465695117165251699935296644992760434412206519036326801985540077086806228464411438175416317958133959586691011057808141234300883692932436022122663781160823472736298914304029803683697389880156222528045623045436226522780569221934366978973936408257173404970887883159180079205465195884955186605310665826226975454171035075016280417683011623124380312403051425940705920625660567584606796218038626946234225864706878239646554381165540673005572649955471689205027195933802181344480474070023069194083134722979170599617466579222550596280381011 = 32158763574926282399690427421751598974822750157866002942864427634852437380540017586451493854661729909380518733649186624385206737336324813109500237603304026009112696565510846849987937423619262973393969175056759821652138869783215378169757835283584660846583208812725733839059137580944002686113912792569631796916069732431775599320458346937859589815497525828622622652165709271152246464728489927670696601016248559515951932154686633599100945314921834227324381958751184684979824241375253606863601895383658582486045363570755445629865194046700806542078378801136397577730247660070033517187439537339428288763342861366560261414507 * 32158763574926282399690427421751598974822750157866002942864427634852437380540017586451493854661729909380518733649186624385206737336324813109500237603304026009112696565510846849987937423619262973393969175056759821652138869783215378169757835283584660846583208812725733839059137580944002686113912792569631796916069732431775599320458346937859589815497525828622622652165709271152246464728489927670696601016248559515951932154686633599100945314921834227324381958751184684979824241375253606863601895383658582486045363570755445629865194046700806542078378801136397577730247660070033517187439537339428288763342861366560261446073 real 0m0.040s user 0m0.028s sys 0m0.012s
Common Factor Attacks
厳密には素因数分解アルゴリズムではないが、異なる二つの合成数が共通の素因数を持つ場合、合成数のbit数によらず高速に素因数を発見することができる。 具体的には、二つの合成数の最大公約数を求めるだけでよい。 最大公約数はユークリッドの互除法によりO(log n)で計算できる。
# common_factor.py
import sys
from itertools import combinations
from fractions import gcd
def common_factor(nums):
for x, y in combinations(nums, 2):
p = gcd(x, y)
if p > 1:
yield (x, y, p)
if __name__ == '__main__':
nums = [long(line) for line in sys.stdin]
for (x, y, p) in common_factor(nums):
print "%d = %d * %d" % (x, p, x/p)
print "%d = %d * %d" % (y, p, y/p)
次の実行結果の通り、共通の素因数を持つ二つの合成数において、4096bitの素因数分解が1秒未満でできている。 ただし、フェルマー法同様、そのような合成数の組が見つかる確率は非常に小さい。
$ cat nums.txt
888441243015294102172292601067855336326287332703328874021611868561802676986077663943061900147565825499428361944793850926519424097012115662952834330500310063991436323100617725235375735298880799081272752074568088718329068470304011938234831484558254069480648619902970082268415966404044075480346440255780228196130123810436581349922045986882478809572371947975635740413437554435429379257619351276688537073428555742828688693292492859848882427167863812003144627707272854531325119129193487170520323833544069156767215691882313314486711034734316444571037196897906716976294880217146223065562518732008190599420275772517007856746981048559878253715189698713288512904998165116853386358827890646734214221734467879782398263142454822004494788965625480385484067919453526130347620614241545718696612964390044265731882751353123553453763779825959232927147708421269420912056959539697543946480513329283728353598683158139025560947346471410872276377801166133810455502905521701767968469458987925204027876348478451236285710225645598004719391738052126832244276934547442298148805426263554624932225480952697394164120418952036473986289321176380982895752273796339856601757641830247264277014250785556784199680489778112198639175574756714183440252052420135756327849404813
882780483105414285207808409377995483701534823624777052872370512205663687323012030665008836246779641378802514783938670119510387540080290753721597154663983983212286475091452225736450381989144354075664524345788834071829669406861239991491714910029461119046332512912612670883516196208270591616798036275990449980681900277329633024676002455118771425041101638785071781244260968665572605189724603728273275682283698062764859089709731324988332263001477198449644990702495668478665515369529879885560385782218986658365704612364761778013317060972940381244288115042449405877429592775144911345025062002681580353283005136150916508175831243673046417028814562130687594141584703791313070099191199276692294970536889058503070974721545833923098157515201071570986923012040954935154602172279935762508063846354593464196649444856619523404920553516291837255896358713883164374689514763790547373493144894084431571218539012930850277919063121821877687816010348113929847662784908586555748228552075893306686302435532676222920629055063379540147073916377160422832146045388449909494920704515440022193230206490791891945947544516000231822642964432381556023595277757817268369138829082830952126257601345645720002969113016584602034934040121873671897644510285530925114503625051
$ time { cat nums.txt | python common_factor.py; }
888441243015294102172292601067855336326287332703328874021611868561802676986077663943061900147565825499428361944793850926519424097012115662952834330500310063991436323100617725235375735298880799081272752074568088718329068470304011938234831484558254069480648619902970082268415966404044075480346440255780228196130123810436581349922045986882478809572371947975635740413437554435429379257619351276688537073428555742828688693292492859848882427167863812003144627707272854531325119129193487170520323833544069156767215691882313314486711034734316444571037196897906716976294880217146223065562518732008190599420275772517007856746981048559878253715189698713288512904998165116853386358827890646734214221734467879782398263142454822004494788965625480385484067919453526130347620614241545718696612964390044265731882751353123553453763779825959232927147708421269420912056959539697543946480513329283728353598683158139025560947346471410872276377801166133810455502905521701767968469458987925204027876348478451236285710225645598004719391738052126832244276934547442298148805426263554624932225480952697394164120418952036473986289321176380982895752273796339856601757641830247264277014250785556784199680489778112198639175574756714183440252052420135756327849404813 = 30094504000188239020151081364836978642515567675718773053464099309910028925560991426570939310582088936607583861708569872351068260656205876930160648487885561172643428908022945411031069520952948240317537081178371267446884744948251284183621948670880897605924632767128833848853213483976135804220279591102249759334548097546952384012455150437958938593286685443258091084425290523032227270109892468085164321266154276823025471164370162127401053735334983910367679071627280601662913673154344586359423856921429216035850978379728954418530655343614680839905243246402183753209499080677780527109087925474741203711791946310699728133647 * 29521710775154724738536772527584123117558710613185797920601145487104166977711662880879179271048876405095288045905017016692076383141051745304194231171302588112692719379484816312610044960068902592215853578316022870490297200043945853963854204074128066849603964094184178301412655771436359730063212919768104154891433003679579996888043953631799144984305859356177785182529634481307702685821400386110992830487859752626831497462372631255954636435957631272568434117771136620170717336511538556050260061771074028879608580486351337906026632036083341002564850216617332234416072504420109009193206508593690307919068066529130502806179
882780483105414285207808409377995483701534823624777052872370512205663687323012030665008836246779641378802514783938670119510387540080290753721597154663983983212286475091452225736450381989144354075664524345788834071829669406861239991491714910029461119046332512912612670883516196208270591616798036275990449980681900277329633024676002455118771425041101638785071781244260968665572605189724603728273275682283698062764859089709731324988332263001477198449644990702495668478665515369529879885560385782218986658365704612364761778013317060972940381244288115042449405877429592775144911345025062002681580353283005136150916508175831243673046417028814562130687594141584703791313070099191199276692294970536889058503070974721545833923098157515201071570986923012040954935154602172279935762508063846354593464196649444856619523404920553516291837255896358713883164374689514763790547373493144894084431571218539012930850277919063121821877687816010348113929847662784908586555748228552075893306686302435532676222920629055063379540147073916377160422832146045388449909494920704515440022193230206490791891945947544516000231822642964432381556023595277757817268369138829082830952126257601345645720002969113016584602034934040121873671897644510285530925114503625051 = 30094504000188239020151081364836978642515567675718773053464099309910028925560991426570939310582088936607583861708569872351068260656205876930160648487885561172643428908022945411031069520952948240317537081178371267446884744948251284183621948670880897605924632767128833848853213483976135804220279591102249759334548097546952384012455150437958938593286685443258091084425290523032227270109892468085164321266154276823025471164370162127401053735334983910367679071627280601662913673154344586359423856921429216035850978379728954418530655343614680839905243246402183753209499080677780527109087925474741203711791946310699728133647 * 29333611316534492642114695914361875827353936602039611693387782601066021747108940770284205684743192567435597448844270242635706940719055360846091860177842785258548515263955446733702873336801213730696067625102999759204582858629776679167503390124589638053007027785917641140606694678963338089301468801125801144536374695406033850449781597771312803137398214729077301657987716032422612267398441299477065463230629907039871180083779210842889144946322439175733816114014531390999811837034265672904377154925113301930277530009488651707667660960383253017461910446471490714077271746064722611111519237800570978700971166429844348550133
real 0m0.029s
user 0m0.028s
sys 0m0.000s
関連リンク
Msieveを使って大きな数を素因数分解してみる
「単純な素因数分解アルゴリズムを実装してみる」では、実装の容易な素因数分解アルゴリズムをPythonで実装し、一般的なPCで解くことができるbit数などを調べた。 より大きな数を素因数分解できるアルゴリズムには、楕円曲線法（ECM）、複数多項式二次ふるい法（MPQS）、一般数体ふるい法（GNFS）が知られている。 ここでは、これらの高速な実装であるMsieveを用いて素因数分解を行い、一般的なPCで解くことができるbit数などを調べてみる。
環境
Ubuntu 14.04.3 LTS 64bit版、Intel Core i5-4200U (1.6GHz * 2 * 2)
$ uname -a Linux vm-ubuntu64 3.19.0-25-generic #26~14.04.1-Ubuntu SMP Fri Jul 24 21:16:20 UTC 2015 x86_64 x86_64 x86_64 GNU/Linux $ lsb_release -a No LSB modules are available. Distributor ID: Ubuntu Description: Ubuntu 14.04.3 LTS Release: 14.04 Codename: trusty $ grep "processor\|model name" /proc/cpuinfo processor : 0 model name : Intel(R) Core(TM) i5-4200U CPU @ 1.60GHz processor : 1 model name : Intel(R) Core(TM) i5-4200U CPU @ 1.60GHz processor : 2 model name : Intel(R) Core(TM) i5-4200U CPU @ 1.60GHz processor : 3 model name : Intel(R) Core(TM) i5-4200U CPU @ 1.60GHz
Msieveのインストール
Msieveは下のURLからWindows用バイナリまたはソースコードの形でダウンロードできる。
ここでは、Linux上でソースコードをコンパイルして動かすこととする。 まず、必要となるパッケージをインストールする。
$ sudo apt-get install build-essential libgmp3-dev zlib1g-dev libecm-dev
ソースコードをダウンロード・展開し、ECMも使える形でコンパイルする。
$ wget "http://downloads.sourceforge.net/project/msieve/msieve/Msieve%20v1.52/msieve152.tar.gz?r=&ts=1452107977&use_mirror=jaist" -O msieve152.tar.gz $ tar xvf msieve152.tar.gz $ cd msieve-1.52/ $ make all ECM=1
コンパイルが終わると、次のようにしてMsieveが使えるようになる。
$ ./msieve --help
Msieve v. 1.52 (SVN unknown)
usage: ./msieve [options] [one_number]
numbers starting with '0' are treated as octal,
numbers starting with '0x' are treated as hexadecimal
options:
-s <name> save intermediate results to <name>
instead of the default msieve.dat
-l <name> append log information to <name>
instead of the default msieve.log
-i <name> read one or more integers to factor from
<name> (default worktodo.ini) instead of
from the command line
-m manual mode: enter numbers via standard input
-q quiet: do not generate any log information,
only print any factors found
-d <min> deadline: if still sieving after <min>
minutes, shut down gracefully (default off)
-r <num> stop sieving after finding <num> relations
-p run at idle priority
-v verbose: write log information to screen
as well as to logfile
-t <num> use at most <num> threads
elliptic curve options:
-e perform 'deep' ECM, seek factors > 15 digits
quadratic sieve options:
-c client: only perform sieving
number field sieve options:
[nfs_phase] "arguments"
where the first part is one or more of:
-n use the number field sieve (80+ digits only;
performs all NFS tasks in order)
-nf <name> read from / write to NFS factor base file
<name> instead of the default msieve.fb
-np perform only NFS polynomial selection
-np1 perform stage 1 of NFS polynomial selection
-nps perform NFS polynomial size optimization
-npr perform NFS polynomial root optimization
-ns perform only NFS sieving
-nc perform only NFS combining (all phases)
-nc1 perform only NFS filtering
-nc2 perform only NFS linear algebra
-ncr perform only NFS linear algebra, restarting
from a previous checkpoint
-nc3 perform only NFS square root
the arguments are a space-delimited list of:
polynomial selection options:
polydegree=X select polynomials with degree X
min_coeff=X minimum leading coefficient to search
in stage 1
max_coeff=X maximum leading coefficient to search
in stage 1
stage1_norm=X the maximum norm value for stage 1
stage2_norm=X the maximum norm value for stage 2
min_evalue=X the minimum score of saved polyomials
poly_deadline=X stop searching after X seconds (0 means
search forever)
X,Y same as 'min_coeff=X max_coeff=Y'
line sieving options:
X,Y handle sieve lines X to Y inclusive
filtering options:
filter_mem_mb=X try to limit filtering memory use to
X megabytes
filter_maxrels=X limit the filtering to using the first
X relations in the data file
filter_lpbound=X have filtering start by only looking
at ideals of size X or larger
target_density=X attempt to produce a matrix with X
entries per column
X,Y same as 'filter_lpbound=X filter_maxrels=Y'
linear algebra options:
skip_matbuild=1 start the linear algebra but skip building
the matrix (assumes it is built already)
la_block=X use a block size of X (512<=X<=65536)
la_superblock=X use a superblock size of X
cado_filter=1 assume filtering used the CADO-NFS suite
square root options:
dep_first=X start with dependency X, 1<=X<=64
dep_last=Y end with dependency Y, 1<=Y<=64
X,Y same as 'dep_first=X dep_last=Y'
楕円曲線法 (ECM; Elliptic Curve factorization Method)
素因数分解アルゴリズムは、計算量が最小の因数pの大きさに依存するものと合成数nの大きさのみに依存するものに分けることができる。 楕円曲線法（ECM）は前者のタイプのアルゴリズムであり、合成数nが比較的小さな因数pを持つか調べる際に用いられる。 ECMはp-1法と呼ばれるアルゴリズムを楕円曲線上の演算で行うものとなっている。
Msieveでは、-eオプションをつけることにより、ECMを用いて10進数15桁以上の因数を探すようにできる。
64bitと448bitの素数を掛け合わせた512bitの合成数について、ECMを有効にした上でMsieveを実行すると次のようになる。
$ ./msieve -q -v -e 0xef00fd8401c4621fa6db62a4d261d73bd0206e632c8c2910570b451f5cab294329d7febb6970381db71c61836d022cbab86a2672da80aea52dc163db9755463b Msieve v. 1.52 (SVN unknown) Thu Jan 7 22:39:29 2016 random seeds: ee949200 621c9755 factoring 12517648286097310891153129740876162488110165257671936603803904646250033327157265537697495023427880408632733768429321903568951073109470815757579744387089979 (155 digits) searching for 15-digit factors searching for 20-digit factors searching for 25-digit factors 200 of 214 curves completed 214 ECM curves searching for 30-digit factors ECM stage 1 factor found prp20 factor: 17406450469679373617 prp135 factor: 719138477307710723633158312040684228975153898269380132563831433336795351605636565609003425835826789865509189242440717516890584669490987 elapsed time 00:00:50
上の結果から、ECMを用いて1分程度で素因数分解できていることがわかる。
複数多項式二次ふるい法 (MPQS; Multiple Polynomial Quadratic Sieve)
複数多項式二次ふるい法（MPQS）は、計算量が合成数nの大きさに依存するタイプのアルゴリズムである。 MPQSは二次ふるい法（QS）をベースとしたアルゴリズムであり、その概要は次のようなものである。
- sqrt(n)以上のxに対して
y = x^2 - nを計算し、それらの中からfactor baseと呼ばれる素数の集合のみを因数に持つ値を残す（ふるいにかける） - 残った値の中から、各因数の個数をもとに掛け合わせたとき平方数となる組み合わせを探す
- 見つかった組み合わせから
(a-b)*(a+b) = 0 (mod n)となるa、bが計算でき、左辺の各項をnで割った余りがnの素因数として得られる
m個の y = x^2 - n を掛け合わせた結果が平方数 b^2 となるとき (x1^2 - n) * (x2^2 - n) * ... * (xm^2 - n) = y1 * y2 * ... * ym = b^2 (mod n) => (x1*x2*...*xm)^2 = b^2 (mod n) a = (x1*x2*...*xm) とおくと a^2 = b^2 (mod n) => (a-b)*(a+b) = 0 (mod n)
MPQSは、1の手順におけるy = x^2 - nの部分をy = (A*x+B)^2 - nとし、複数のA、Bからなる多項式を用いてyを計算する。
それぞれの多項式におけるふるい分けは独立して行うことができるので、並列化に向いている。
下の実行結果の通り、128bitの素数を掛け合わせた256bitの合成数が4分ほどで素因数分解できた。
$ ./msieve -q -v 0x00b84b8a6ba2595165e095dd8226e513f5e42fa956339f0f5b4743b24c223c8603 Msieve v. 1.52 (SVN unknown) Thu Jan 7 22:56:23 2016 random seeds: b33d7dc5 fd4ec6ff factoring 83359033011986879296513856404319061797526856544841606538718515630707528467971 (77 digits) searching for 15-digit factors commencing quadratic sieve (77-digit input) using multiplier of 19 using generic 32kb sieve core sieve interval: 12 blocks of size 32768 processing polynomials in batches of 17 using a sieve bound of 923773 (36471 primes) using large prime bound of 92377300 (26 bits) using trial factoring cutoff of 26 bits polynomial 'A' values have 10 factors sieving in progress (press Ctrl-C to pause) 36600 relations (18916 full + 17684 combined from 196785 partial), need 36567 36600 relations (18916 full + 17684 combined from 196785 partial), need 36567 sieving complete, commencing postprocessing begin with 215701 relations reduce to 51990 relations in 2 passes attempting to read 51990 relations recovered 51990 relations recovered 42678 polynomials attempting to build 36600 cycles found 36600 cycles in 1 passes distribution of cycle lengths: length 1 : 18916 length 2 : 17684 largest cycle: 2 relations matrix is 36471 x 36600 (5.3 MB) with weight 1096217 (29.95/col) sparse part has weight 1096217 (29.95/col) filtering completed in 4 passes matrix is 26182 x 26246 (4.1 MB) with weight 869648 (33.13/col) sparse part has weight 869648 (33.13/col) saving the first 48 matrix rows for later matrix includes 64 packed rows matrix is 26134 x 26246 (2.6 MB) with weight 628988 (23.97/col) sparse part has weight 424260 (16.16/col) commencing Lanczos iteration memory use: 2.6 MB lanczos halted after 415 iterations (dim = 26134) recovered 18 nontrivial dependencies prp39 factor: 288180072604771133716023733756993741403 prp39 factor: 289260226283275563579648656678444936057 elapsed time 00:03:45
一般数体ふるい法 (GNFS; General Number Field Sieve)
一般数体ふるい法（GNFS）は、10進数で100桁（332bit）以上の素因数分解において現時点でもっとも高速とされているアルゴリズムである。
これは、二次ふるい法（QS）における多項式f(x) = x^2 + 1を一般化し、代数体における整数環の素イデアルもfactor baseとして用いてふるいにかけるものである。
Msieveでは、-nオプションをつけることにより、10進数80桁以上の数に対してGNFSを用いるようにできる。
GNFSを用いる場合、素因数分解したい数はworktodo.iniというファイルに書く必要がある。
0x00d1696798485ddf46936b55068e1f36ef922d13748fa4e6e7a98afd74c2ab4148eaaf
上の136bitの素数を掛け合わせた272bitの合成数について素因数分解を試みたが、次の実行結果のように5分程度経過したところで異常終了し、解を得ることはできなかった。
$ time ./msieve -q -v -n Msieve v. 1.52 (SVN unknown) Fri Jan 8 20:28:33 2016 random seeds: 5f4b2d9f a557e83a factoring 6207544969206898368897868555163969060703241410703468520600725959726120571872275119 (82 digits) searching for 15-digit factors commencing number field sieve (82-digit input) commencing number field sieve polynomial selection polynomial degree: 4 max stage 1 norm: 6.85e+13 max stage 2 norm: 2.56e+13 min E-value: 8.57e-08 poly select deadline: 307 time limit set to 0.09 CPU-hours expecting poly E from 1.58e-07 to > 1.81e-07 searching leading coefficients from 1 to 2182505 deadline: 5 CPU-seconds per coefficient coeff 12 specialq 1 - 1456 other 2346 - 5630 aprogs: 159 entries, 524 roots 12 1633159927 150811568598986049453 12 28365837577 150811568151956063017 12 992431769 150811568599326420050 12 267101113 150811568597234510506 12 8869155629 150811568502795701193 12 11044194347 150811568615862920377 12 11414099191 150811568240381948784 save 1.412830e-11 -4.5116 4452360.73 8.775657e-08 rroots 2 save 1.412830e-11 -4.5116 4452360.73 8.775657e-08 rroots 2 12 14236060951 150811568631345977424 12 14224094453 150811568464692667845 12 14558241583 150811568594085519693 hashtable: 512 entries, 0.01 MB coeff 24 specialq 1 - 1889 other 2151 - 5162 aprogs: 153 entries, 420 roots 24 13870676035 126816907473782790413 24 8655358865 126816907417857728947 24 23990116715 126816907560250609681 save 1.860294e-11 -4.9666 2745370.32 1.002695e-07 rroots 2 save 1.860294e-11 -4.9666 2745370.32 1.002695e-07 rroots 2 24 13014280909 126816907670959530931 24 10695397835 126816907415129233556 24 4921213553 126816907405855042120 24 4244822651 126816907430348394592 24 7987263461 126816907609450597592 24 19252292627 126816907491466150455 24 12552609739 126816907094190571465 hashtable: 512 entries, 0.01 MB (snip) polynomial selection complete R0: -87070511241854355656 R1: 5110975483 A0: -2362525215465747901869945 A1: -36257592644644446639 A2: -10330966943699 A3: 49742175 A4: 108 skew 620857.53, size 2.743e-11, alpha -5.586, combined = 1.363e-07 rroots = 2 generating factor base factor base: found 120001 rational and 119728 algebraic entries factor base complete: 138618 rational roots (max prime = 1851271) 138552 algebraic roots (max prime = 1851259) a range: [-4410256, 4410256] b range: [1, 4294967295] number of hash buckets: 29 sieve block size: 65536 maximum RFB prime: 1851271 RFB entries: 138618 medium RFB entries: 6542 resieved RFB entries: 6374 small RFB prime powers: 2 projective RFB roots: 2 RFB trial factoring cutoff: 54 or 81 bits single large prime RFB range: 21 - 26 bits double large prime RFB range: 42 - 50 bits triple large prime RFB range: 65 - 76 bits maximum AFB prime: 1851259 AFB entries: 138552 medium AFB entries: 6518 resieved AFB entries: 6360 small AFB prime powers: 133 projective AFB roots: 4 AFB trial factoring cutoff: 54 or 81 bits single large prime AFB range: 21 - 26 bits double large prime AFB range: 42 - 50 bits triple large prime AFB range: 65 - 76 bits multiplying 939255 primes from 1851259 to 16777216 *** Error in `./msieve': realloc(): invalid old size: 0x0000000001155780 *** Aborted (core dumped) real 5m24.154s user 5m9.844s sys 0m13.660s
