プロジェクトオイラー問187

http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%20187
Problem 187 「半素数」 †
合成数とは2つ以上の素因数を含む整数のことである. 例えば15 = 3 × 5, 9 = 3 × 3, 12 = 2 × 2 × 3が合成数である.

30以下には丁度2つの素因数を含む合成数 (異なる素因数でなくてもよい) が, 10個存在する. 4, 6, 9, 10, 14, 15, 21, 22, 25, 26がそうである.

合成数n < 10^8について, 丁度2つの素因数を含む合成数 (異なる素因数でなくてもよい) はいくつあるか.


素数の個数を全部求めて足し算してみました。
2.4ギガヘルツPCで計算時間5秒。
遅いのか早いのか微妙?

1億くらいならbool型をうまく使えば100メガだからふるいにかけたほうが速いです。
今回の私の手法は100億とかで篩よりもメモリが節約できる利点があります。

私の手法の場合
100億でもメモリ使用量100メガバイト未満、CPU使用率20%で計算時間10分。
n<10^10の場合の答え
ans=1493776443
となります。


100億の場合をためし割りで行うと素数が10億個として一つ試すのに10000回かかるとします。
10億*10000=10兆回の計算となり一秒1億回の計算としても10万秒計算時間がかかることになります。
計算に1日以上かかります。

それでは、篩の場合私の手法より早いですが、bool型を使っても10ギガバイトのメモリを食います。

なので私の手法は省メモリでそれなりの速度ということです。



#include<stdio.h>
#include<map>
#include<vector>
#include<set>
#include<iostream> 

const __int64 L=100000;
const __int64 LIMIT=L*L;

bool isP(__int64 n){
	if(n<2)return false;
	for(__int64 i=2;i*i<=n;i++){
		if(n%i==0)return false;
	}
	return true;
} 

std::set<__int64> ps2;
std::vector<__int64> ps;

void f(__int64 p,__int64 n){
	if(n>1){
		ps2.insert(n);
}
	for(__int64 i=p;i<ps.size();i++){
		__int64 p2=ps[i];
		if(L<p2*n)break;
		f(i,n*p2);
	}
}
 
__int64 main(){
	std::map<__int64,__int64> memo;
	std::map<__int64,__int64>::iterator it,it2;
	
	for(__int64 i=2;i<=L;i++){
		if(isP(i))ps.push_back(i);
	}
	memo[LIMIT]=1;
	for(__int64 i=0;i<ps.size();i++){
		__int64 p=ps[i];
		it=memo.end();
		it--;
		while(1){
			__int64 f=(*it).first;
			__int64 s=(*it).second;
			if(f<p)break;
			memo[f/p]+=s;
			it=memo.lower_bound(f);
		it--;
		}
	}
	
	
 	for(it=memo.begin(),it2=it,it2++;(*it).first<L;it++,it2++){
	__int64 f1=(*it).first;
		__int64 f2=(*it2).first;
		memo[f1]=LIMIT/f1-LIMIT/f2-memo[f1];
	}
	it=memo.begin();
	it2=it;
	it2++;
	while((*it).first<L){
		__int64 f1=(*it).first;
		__int64 f2=(*it2).first;
		if(f2-f1>2){
			printf("error");
		}
		if(f2-f1==2){
			__int64 d=LIMIT/(f1+1)-LIMIT/(f1+2);
			memo[f1+1]=d;
			memo[f1]-=d;
			it=memo.lower_bound(f1+1);
			it2=it;
			it2++;
		}else{
 			it++;
			it2++;
		}
	}
	f(0,1);
	it=memo.lower_bound(L);
	
	memo.erase(it,memo.end());
	
	std::map<__int64,__int64>::reverse_iterator rIt;
	for(rIt=memo.rbegin();rIt!=memo.rend();rIt++){
		std::set<__int64>::iterator sIt;
		__int64 f=(*rIt).first;
		__int64 s=(*rIt).second;
 		for(sIt=ps2.begin();sIt!=ps2.end();sIt++){
			__int64 p2=(*sIt);
			if(f<p2)break;
			memo[f/p2]-=s;
		}	
	}
	__int64 sum=0;
 	for(rIt=memo.rbegin();rIt!=memo.rend();rIt++){
		sum+=(*rIt).second;
		memo[(*rIt).first]=sum;
	}
	__int64 ans=0;
	for(__int64 i=0;i<ps.size();i++){
		ans+=ps.size()-i;
		it=memo.lower_bound(ps[i]);
		if(memo.end()==it)continue;
		ans+=(*it).second;
	}
	std::cout<<"ans="<<ans<<"\n";
}
最終更新:2016年01月21日 05:45