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
