Top / ほぼ高速な素因数分解のプログラム
HTML convert time: 0.023 sec.

ほぼ高速な素因数分解のプログラム

Last-modified: 2015-03-21 (土) 15:41:23

ほぼ高速な素因数分解のプログラム

これまでに総当たりの方法とエラトステネスの篩を利用した方法を紹介しました。総当たりの方法では、小さい素因数を持っている数に対しては絶大な威力を見せつけました。しかし、この方法の弱点としては大きな数であり、大きな素因数を持つ数に対しては多くの時間を費やしてしまうという点があります。
そこで考えた方法は、ある程度の数までは総当たりで行い、それ以上の素因数を探すのにはエラトステネスの篩を利用する方法です。

しかし、タイトルに「ほぼ」という言葉が入っているのは、エラトステネスの篩で紹介したように用意できる配列に限界があるからです。そのため次の場合は正確な素因数分解ができません。

 ・202,500,000,000,000,000(約20京)を超える素数を含む
 ・450,000,000(4.5億)以上の素数を2つ含む場合

今、この弱点を克服した改良版を考え中です…
知識のある方からアドバイスをいただけたら嬉しいです。

<ソースコード>

/*                         */
/*    factorization.exe    */
/*                         */

#include <iostream>
#include <cmath>
#include <ctime>
using namespace std;

#define N  450000000      
#define PN 24000000      
int arr[N + 1];     // 篩の際に使う箱
int pn[PN];         // 素数を入れる箱

int main()
{
	__int64 X;
	cout << "\n素因数分解" << "\n* MAX:9223372036854775807 *" << endl;
	cout << "    X=";
	cin >> X;
	cout << "\n      " << X << " = ";

	clock_t start = clock();

	__int64 n = X;

	// paturn 1  ( X <= 9*10^14 (15桁))
	if (sqrt(X) <= 30000000){
		for (int i = 2; i <= sqrt(X); i++){
			int ruijo = 0;
			if (n%i == 0){
				cout << i;
				while (n%i == 0){        // 割り切れるまで行う
					ruijo++;
					n = n / i;
				}
				if (ruijo == 1){
					cout << " ";
				}
				if (ruijo >= 2){
					cout << "^" << ruijo << " ";
				}
				if (n == 1){             // 割り切れたら終了
					break;
				}
			}
		}

		if (n > sqrt(X)){
			if (n == X){
				cout << "Prime Number";
			}
			else{
				cout << n;
			}
		}
		cout << endl;
	}


	// paturn 2  ( X >10^16 )
	else{

		// stageⅠ  ( 30000000 の数までで割れるところまで割る )
		for (int i = 2; i <= 30000000; i++){
			int ruijo = 0;
			if (n%i == 0){
				cout << i;
				while (n%i == 0){      // 割り切れるまで行う
					ruijo++;
					n = n / i;
				}
				if (ruijo == 1){
					cout << " ";
				}
				if (ruijo >= 2){
					cout << "^" << ruijo << " ";
				}
				if (n == 1){
					break;              // 割り切れたら終了
				}
			}
		}

		// stageⅡ  ( stageⅠ で割ったところから篩を使って因数分解 )
		if (n != 1){
			__int64 quo = n;
			int pncount;

			// stageⅡ-1  ( sqrt(n) が 450000000 より小さいとき )
			if (sqrt(n) <= N){

				for (int i = 3; i <= sqrt(n); i = i + 2){
					arr[i] = i;
				}

				pncount = 1;
				pn[1] = 2;
				pncount++;

				for (int i = 3; i <= sqrt(n); i = i + 2){
					if (arr[i]){
						pn[pncount] = i;
						pncount++;
						for (int j = 3; i*j <= sqrt(n); j = j + 2){
							arr[i*j] = 0;
						}
					}
				}

				pncount--;

				for (int i = 1; i <= pncount; i++){
					int sub = pn[i];
					if (quo%sub == 0){
						int ruijo = 0;
						cout << sub;
						while (quo%sub == 0){
							quo = quo / sub;
							ruijo++;
						}
						if (ruijo == 1){
							cout << " ";
						}
						if (ruijo >= 2){
							cout << "^" << ruijo << " ";
						}
						if (quo == 1){
							break;
						}		
					}
				}
			}

			// stageⅡ-2  ( sqrt(n) が 450000000 より大きいとき )
			if (sqrt(n) > N){

				// 450000000 までの素数を入れる
				for (int i = 3; i < N; i = i + 2){
					arr[i] = i;
				}

				pncount = 1;
				pn[1] = 2;
				pncount++;

				for (int i = 3; i < N; i = i + 2){
					if (arr[i]){
						pn[pncount] = arr[i];
						pncount++;
						for (int j = 3; i*j < N; j = j + 2){
							arr[i*j] = 0;
						}
					}
				}

				pncount--;

				// 450000000 までの素数で割る
				for (int i = 1; i <= pncount; i++){
					int sub = pn[i];
					if (quo%sub == 0){
						int ruijo = 0;
						cout << sub;
						while (quo%sub == 0){
							quo = quo / sub;
							ruijo++;
						}
						if (ruijo == 1){
							cout << " ";
						}
						if (ruijo >= 2){
							cout << "^" << ruijo << " ";
						}
						if (quo == 1){
							break;
						}				
					}
				} 
			}

			// stageⅢ ( 割り切れなかった場合 )
			if (quo != 1){
				if (sqrt(quo) > 450000000){
					if (quo != X){
						cout << quo << endl;
						cout << " ※ただし " << quo << " は素数 もしくは 450000000以上の素数2つから成る数";
					}
					else{
						cout << "素数 もしくは 450000000以上の素数2つから成る数";
					}
				}
				else{
					if (quo == X){
						cout << "Prime Number";
					}
					else{
						cout << quo;
					}
				}
			}
			cout << endl;
		}
	}

	clock_t end = clock();
	cout << "\ntime=" << (double)(end - start) / CLOCKS_PER_SEC << "s" << endl;
	return 0;
}

<計算結果>

 X : time (総当たりでかかった時間)
 98765432123456789 : 0.378s (7.173s)
 (98765432123456789 = 449*494927*444444443)

 1234567890987654321 : 0.384s (25.313s)
 (1234567890987654321 = 3^2*7*19*928163*1111211111)

アドバイスやご意見をいただけたらなと思っております。
どんなささいなことでもいいのでよろしくお願いします。