|
|
試し割り法 + ρ法で素因数分解をします。
ISPRIME(X)はミラー・ラビン法などの素数判定ルーチンです。ルーチンは省略しています。
ポラード・ロー法(ρ法)をC++で書いています。
また、boostライブラリーを使い、任意精度整数型(cpp_int型)を使用しています。
ビルドしたdllファイルは大き過ぎるので、今回も掲示板への投稿は控えさせて頂きます。
C++開発環境を持たない方のため、VC++2015(x86)にてビルドしたファイルをアップローダーへ
upしておきました。(factor.zip) 下記URLは現時点より１ヶ月間有効です。
ダウンロードパス:shibacchi
試しに2^N-1 として Nを100～120(2^100-1 ～ 2^120-1)までを素因数分解してみました。
OPTION ARITHMETIC RATIONAL
LET L=1000000
DIM S(L)
FOR I=2 TO L !'エラトステネスの篩
IF S(I)=0 THEN
FOR J=I*I TO L STEP I
LET S(J)=1
NEXT J
END IF
NEXT I
FOR K=100 TO 120
LET N=2^K-1
PRINT "2 ^";K;"- 1"
PRINT N;"-->";
IF ISPRIME(N)=0 THEN !'ミラー・ラビン法など。合成数なら試し割り
LET NN=N
FOR I=2 TO L
IF S(I)=0 THEN
DO WHILE MOD(N,I)=0
LET N=N/I
PRINT I;
LOOP
END IF
NEXT I
IF NN<>N AND N>1 THEN PRINT N ELSE PRINT
CALL FACTOR(N)
ELSE
PRINT " 素数"
END IF
PRINT
NEXT K
END
EXTERNAL SUB FACTOR(N)
OPTION ARITHMETIC RATIONAL
IF N>1 AND ISPRIME(N)=0 THEN !'ミラー・ラビン法などで確認。合成数なら
PRINT N;"-->";
LET D=POLLARDRHO(N)
PRINT D;N/D
CALL FACTOR(D) !'再帰呼び出し
CALL FACTOR(N/D)
END IF
END SUB
EXTERNAL FUNCTION POLLARDRHO(X) !'ポラード・ロー素因数分解法(ρ法)
OPTION ARITHMETIC RATIONAL
LET N$=STR$(X)
CALL POLLARD(N$)
FOR I=1 TO LEN(N$) !' POLLARDRHO=VAL(N$) とできない
LET P$=N$(I:I)
IF ORD(P$)>=48 AND ORD(P$)<=57 THEN LET B$=B$ & P$
NEXT I
LET POLLARDRHO=VAL(B$)
END FUNCTION
EXTERNAL SUB POLLARD(X$) !'素数を入れないこと。無限ループのおそれ。
OPTION ARITHMETIC RATIONAL
ASSIGN "pollard.dll","pollard"
END SUB
----------------------------------------------------------------------------
pollard.cpp
#include <boost/multiprecision/cpp_int.hpp>
using namespace boost::multiprecision;
using namespace std;
cpp_int atocpp(char str[])
{
cpp_int result = 0;
int i=0;
while (str[i]>='0' && str[i]<='9'){
result=result*10+(str[i++] - '0');
}
return result;
}
cpp_int mulmod(cpp_int a,cpp_int b,cpp_int m)
{
cpp_int result=0;
a %= m;
while (b>0) {
if (b % 2 == 1) result = (result+a) % m;
a = (a*2) % m;
b = b / 2;
}
return result;
}
cpp_int pollardrho(cpp_int num,unsigned long c){
cpp_int x=1,y,p=1,n=2,i=1;
while (p == 1){
++i;
if (i == n){
y = x;
n *= 2;
}
x = (mulmod(x,x,num) + c) % num;
if (y == x) return pollardrho(num,c+1);
p = gcd(num,abs(x-y));
}
return p;
}
extern "C" __declspec(dllexport) void pollard(char x[])
{
cpp_int num,n;
int i=0,tmp,len;
num=atocpp(x);
len=strlen(x);
for(int j=0;j<len;j++) x[j]='\0';
n=pollardrho(num,1u);
while (n != 0) {
x[i] = (char)(n % 10) + '0';
n /= 10;
i++;
}
x[i] = '\0';
i--;
for(int j = 0;j <= i/2;j++) {
tmp = x[j];
x[j] = x[i-j];
x[i-j] = tmp;
}
}
実行結果
2 ^ 100 - 1
1267650600228229401496703205375 --> 3 5 5 5 11 31 41 101 251 601 1801 4051 8101 268501
2 ^ 101 - 1
2535301200456458802993406410751 -->
2535301200456458802993406410751 --> 7432339208719 341117531003194129
2 ^ 102 - 1
5070602400912917605986812821503 --> 3 3 7 103 307 2143 2857 6529 11119 43691 131071
2 ^ 103 - 1
10141204801825835211973625643007 -->
10141204801825835211973625643007 --> 2550183799 3976656429941438590393
2 ^ 104 - 1
20282409603651670423947251286015 --> 3 5 17 53 157 1613 2731 8191 858001 308761441
2 ^ 105 - 1
40564819207303340847894502572031 --> 7 7 31 71 127 151 337 29191 106681 122921 152041
2 ^ 106 - 1
81129638414606681695789005144063 --> 3 107 6361 69431 572263032673174337633
572263032673174337633 --> 20394401 28059810762433
2 ^ 107 - 1
162259276829213363391578010288127 --> 素数
2 ^ 108 - 1
324518553658426726783156020576255 --> 3 3 3 3 5 7 13 19 37 73 109 87211 246241 262657 279073
2 ^ 109 - 1
649037107316853453566312041152511 -->
649037107316853453566312041152511 --> 745988807 870035986098720987332873
2 ^ 110 - 1
1298074214633706907132624082305023 --> 3 11 11 23 31 89 683 881 2971 3191 201961 48912491
2 ^ 111 - 1
2596148429267413814265248164610047 --> 7 223 321679 5170159074941697420729913
5170159074941697420729913 --> 26295457 196617958567584409
196617958567584409 --> 616318177 319020217
2 ^ 112 - 1
5192296858534827628530496329220095 --> 3 5 17 29 43 113 127 257 5153 859166727965929937
859166727965929937 --> 15790321 54410972897
2 ^ 113 - 1
10384593717069655257060992658440191 --> 3391 23279 65993 1993423291715412685783
1993423291715412685783 --> 1868569 1066818132868207
2 ^ 114 - 1
20769187434139310514121985316880383 --> 3 3 7 571 32377 174763 524287 194620086937183
194620086937183 --> 1212847 160465489
2 ^ 115 - 1
41538374868278621028243970633760767 --> 31 47 14951 178481 10683848415441845139401
10683848415441845139401 --> 4036961 2646507710984041
2 ^ 116 - 1
83076749736557242056487941267521535 --> 3 5 59 233 1103 2089 174850288360342272373981
174850288360342272373981 --> 3033169 57646075230342349
57646075230342349 --> 107367629 536903681
2 ^ 117 - 1
166153499473114484112975882535043071 --> 7 73 79 937 6553 8191 86113 121369 7830118297
2 ^ 118 - 1
332306998946228968225951765070086143 --> 3 2833 37171 179951 5845385390147915655817
5845385390147915655817 --> 1824726041 3203431780337
2 ^ 119 - 1
664613997892457936451903530140172287 --> 127 239 20231 131071 8257410955834335790279
8257410955834335790279 --> 131105292137 62983048367
2 ^ 120 - 1
1329227995784915872903807060280344575 --> 3 3 5 5 7 11 13 17 31 41 61 151 241 331 1321 61681 4562284561
|
|