投稿者:しばっち
投稿日:2016年 2月28日(日)17時56分29秒
|
|
|
エラトステネスの篩をC++で書いてみました
64bit整数(long long)を使用しています
sqrt関数がdouble型のため、isqrt関数を定義しています
また、ライブラリーのatoll関数がunsignedに対応していないので
atoull関数を定義しています。これにより、2^64-1まで対応できます
bool型(1byte)を使用し、余計な処理はしていません)
ASSIGN文でdllファイルのパスを指定してください
OPTION ARITHMETIC DECIMAL_HIGH
LET S=1E+19+1 !'奇数 (1E+19 < 2^64 < 1E+20))
LET E=S+1000 !'奇数
LET S$=STR$(S)
LET E$=STR$(E)
PRINT ERATOS(S$,E$)
END
EXTERNAL FUNCTION ERATOS(S$,E$)
OPTION ARITHMETIC DECIMAL_HIGH
ASSIGN ".\DLL\eratos.dll","eratos",FPU
END FUNCTION
-----------------------------------------------------------------------
eratos.cpp
using namespace std;
unsigned long long isqrt(unsigned long long num)
{
unsigned long long n,n1;
if(num==0) return 0;
n=num/2+1;
n1=(n+(num/n))/2;
while (n1<n) {
n=n1;
n1=(n+(num/n))/2;
}
return n;
}
unsigned long long atoull(char *str)
{
unsigned long long result = 0;
while (*str>='0' && *str<='9'){
result=(result*10)+(*str++ - '0');
}
return result;
}
extern "C" __declspec(dllexport) double eratos(char *s,char *e)
{
unsigned long long n,nn,m,i,j,count=0;
bool *x;
n=atoull(s);
m=atoull(e);
if (n>=m) return -1.0;
x=new bool [m-n+2];
for(i=0;i<=m-n+1;i++) x[i]=true;
for(i=3;i<=isqrt(m);i+=2){
if (n % i==0) nn=n;
else
nn=((n/i)+1)*i;
if (n<=i && x[nn-n]==true) nn+=i;
for(j=nn;j<=m;j+=i) x[j-n]=false;
}
for(j=n;j<=m;j++)
if (j % 2==1 && x[j-n]==true) count++;
delete [] x;
return (double)count;
}
|
|
|
投稿者:たろさ
投稿日:2016年 3月 7日(月)17時56分20秒
|
|
|
> No.3985[元記事へ]
しばっちさんへのお返事です。
> ASSIGN文でdllファイルのパスを指定してください
パスが通らないのでルートで動作報告です。
Intel i7-4790K@4.7GHz Win7 32bt
100億 455052511
200億 882206716
300億 1300005926
400億 1711955433
500億 2119654578
600億 2524038155
700億 2925699539
800億 3325059246
900億 3722428991
1000億 4118054813
100億から1000億 1億刻みで901回 3596.68 秒で計算しました
素数定理
100億と1000億の素数の個数を確認しました。
十進BASICプログラムの計算時間です。
http://blogs.yahoo.co.jp/donald_stinger
|
|
|
投稿者:たろさ
投稿日:2016年 3月 8日(火)08時56分24秒
|
|
|
投稿者:たろさ
投稿日:2016年 3月 9日(水)11時41分43秒
|
|
|
> No.3985[元記事へ]
しばっちさんへのお返事です。
1億STEPと10億STEPの比較DATA
100000000000=1E+11から1000000000000=1E+12までの素数の個数
100000000000 : 4118054813個
1000000000000 : 37607912018個
37607912018個 確認しました。
【テスト機】 【範囲】 【計算時間】
Intel i7-4790K@4.7GHz Win7 32bt FOR i=1E+11 TO 1E+12 STEP 1E+8(11時間27分35秒)
Intel G3258@4.6GHz Win7 32bt FOR i=1E+11 TO 1E+12 STEP 1E+8(12時間21分56秒)
Intel G3258@4.6GHz Win7 32bt FOR i=1E+11 TO 1E+12 STEP 1E+9(13時間53分44秒)
nまで : π(x)
100000000000 : 4118054813個 【STEP 1E+9】(13時間53分44秒)
100000000000 : 4118054813 【STEP 1E+8】(11時間27分35秒)
200000000000 : 8007105059個
200000000000 : 8007105059
300000000000 : 11818439135個
300000000000 : 11818439135
400000000000 : 15581005657個
400000000000 : 15581005657
500000000000 : 19308136142個
500000000000 : 19308136142
600000000000 : 23007501786個
600000000000 : 23007501786
700000000000 : 26684074310個
700000000000 : 26684074310
800000000000 : 30341383527個
800000000000 : 30341383527
900000000000 : 33981987586個
900000000000 : 33981987586
1000000000000 : 37607912018個
1000000000000 : 37607912018
http://blogs.yahoo.co.jp/donald_stinger
|
|
|
投稿者:しばっち
投稿日:2016年 7月17日(日)16時12分34秒
|
|
|
> No.3985[元記事へ]
旧版 eratos.dllをOPEN MP(マルチスレッド化)にとりあえず成功しました。
#3986
VC++2015(x86)にて、ビルド(コンパイル)しています。
こちらの環境では、およそ1.4倍高速化しました。(8スレッドで実行)
但し、OPEN MP上の制約により、unsignedが使えず、2^63-1までになります。
実行にはOPEN MPランタイム(vcomp140.dll)が必要です。(basic.exeと同じフォルダに入れてください)
なお、vcomp140.dllはライセンス上再配布可能らしいので、ビルドしたdllファイルと同梱し
アップローダーにupしておきました。(eratos.zip) 下記 URLの有効期限は現時点より1ヶ月間となります。
ダウンロードパス:shibacchi
OPTION ARITHMETIC NATIVE
LET S=1E+13+1
LET E=S+1E+8
LET S$=STR$(S)
LET E$=STR$(E)
LET T=TIME
PRINT "シングルスレッド"
PRINT ERATOS(S$,E$);"個"
LET T1=TIME-T
PRINT T1;"秒"
LET T=TIME
PRINT "マルチスレッド"
PRINT ERATOS2(S$,E$,0);"個"
LET T2=TIME-T
PRINT T2;"秒"
PRINT T1/T2;"倍"
END
EXTERNAL FUNCTION ERATOS(S$,E$)
OPTION ARITHMETIC NATIVE
ASSIGN "eratos.dll","eratos",FPU
END FUNCTION
EXTERNAL FUNCTION ERATOS2(S$,E$,THREADS) !'THREADS=実行スレッド数(0の時 自動設定)
OPTION ARITHMETIC NATIVE
ASSIGN "eratos_parallel.dll","eratos",FPU
END FUNCTION
-------------------------------------------------------------------------------------------
eratos_parallel.cpp
#include <cstdlib>
#include <omp.h>
using namespace std;
long long isqrt(long long num)
{
long long n,n1;
if(num==0) return 0;
n=num/2+1;
n1=(n+(num/n))/2;
while (n1<n) {
n=n1;
n1=(n+(num/n))/2;
}
return n;
}
extern "C" __declspec(dllexport) double eratos(char *s,char *e,int threads)
{
long long n,nn,m,i,j,count=0;
bool *x;
n=atoll(s);
m=atoll(e);
if (n>=m || n<1) return -1.0;
if (n==1 || n==2) count++;
if (n==1) n+=2;
if (n % 2==0) n++;
x=new bool [m-n+2];
if (x==NULL) return -999.0;
if (threads>0) omp_set_num_threads(threads);
#pragma omp parallel for
for(i=0;i<=m-n+1;i++) x[i]=true;
#pragma omp parallel for private(j,i,nn) schedule(dynamic)
for(i=3;i<=isqrt(m);i+=2){
if (n % i==0) nn=n;
else
nn=((n/i)+1)*i;
if (n<=i && x[nn-n]==true) nn+=i;
for(j=nn;j<=m;j+=i) x[j-n]=false;
}
#pragma omp parallel for
for(j=n;j<=m;j+=2)
if (x[j-n]==true) {
#pragma omp atomic
count++;
}
delete [] x;
return (double)count;
}
|
|
|
戻る