|
> No.4326[元記事へ]
しばっちさんへのお返事です。
> WAIT EVENTとSIGNAL文を外すと何故か動きました。
> 私もエラトステネスの篩を素数リストで篩ってみました。
> VC++2017 + OPEN MPより速くなってしまいました。驚きです。
素晴らしいプログラムを有難うございます。毎度参考にしています。
実は、エラトステネスの篩を素数リストで篩って成功したのは、初めてです。
昨年から失敗続きでした。素数カウントの精度については未確認です。
下記のプログラムで1~1000億まで確認です。1~1兆は、何故か? 誤差が出ます。なお計算時間は5時間です。
最小値の設定です。
LET S=1
LET E=1E11
LET TOTAL=0 2は排除されますが、1をカウントしてます。
------------------------------------
LET S=99E10+1 !pi(1E12),37607912018
LET E=1E12
これは、うまくいく。
なので、まだ暫定版です。私のPARACT PART1~4 実験では、10兆間を計算するのは難してです。
EXTERNAL SUB sqprime ですが、1億step用に作成しましたので、
LET B(c)=i+1
1~1E15 までは確認しました。
並列プログラムについては、
DECLARE STRUCTURE struct3: 1 OF NUMERIC(1 TO 100001)
ここの配列をPARACT PART1 から設定できると、6n+k篩と同じにできるのですが、残念です。
VC++2017 + OPEN MP 500兆までまもなく計算完了します。スレッド自動割り当ては便利です。
8コアCPUが無いので、8スレッドのプログラムを4コアCPUで、これから試してみます。15.17秒
確かに、驚きます。
失礼しました。
Intel 4790K Devil's Canyon S-spec SR219 CPU Overclocking Report
ここ見ると、#of Threads 8
! Sieve of Eratosthenes 1兆 4/28
DECLARE STRUCTURE struct1: 1 OF NUMERIC
DECLARE STRUCTURE struct2: 1 OF NUMERIC(78505)
DECLARE STRUCTURE struct3: 1 OF NUMERIC(1 TO 100001)
DECLARE STRUCTURE struct4: 4 OF NUMERIC
DECLARE SHARED sha OF struct2
DECLARE SHARED shb OF struct3
DECLARE MESSAGE mes2 OF struct1
DECLARE MESSAGE mes3 OF struct1
DECLARE MESSAGE mes4 OF struct1
DECLARE MESSAGE mes5 OF struct4
DECLARE MESSAGE mes6 OF struct4
DECLARE MESSAGE mes7 OF struct4
PARACT PART1
OPTION ARITHMETIC NATIVE
LET k=1000117
LET k3=78505
DECLARE EXTERNAL SUB prime
CALL prime(k,k3)
WAIT EVENT Ok5
LET S=99E10+1 !pi(1E12),37607912018
LET E=1E12 !pi(1E11),4118054813
LET ST=1E7
LET cc=(S-1)/ST+1
DECLARE EXTERNAL SUB sqprime
CALL sqprime(1,100001,k3)
DECLARE NUMERIC B(1 TO 100001)
WAIT EVENT Ok6
GET FROM shb TO B
START Part2
START PART3
START PART4
SEND TO mes5 FROM S,E,ST,CC
SEND TO mes6 FROM S,E,ST,CC
SEND TO mes7 FROM S,E,ST,CC
LET TOTAL=37245934597!1e7=664579,1e9=50847534
DECLARE EXTERNAL FUNCTION ERATOS
!DECLARE NUMERIC X,Y,Z
FOR I=S TO E STEP ST
LET k2=B(cc)
LET L=ERATOS(I,I+ST/4-1,k2)
RECEIVE FROM mes2 TO X
RECEIVE FROM mes3 TO Y
RECEIVE FROM mes4 TO Z
LET L=L+X+Y+Z
LET TOTAL=TOTAL+L
PRINT I-1;I+ST-1;TOTAL;L
!PRINT TOTAL
LET cc=cc+1
NEXT I
END PARACT
PARACT PART2
OPTION ARITHMETIC NATIVE
DECLARE NUMERIC B(1 TO 100001)
GET FROM shb TO B
RECEIVE FROM mes5 TO X,Y,Z,G
LET S=X
LET E=Y
LET ST=Z
DECLARE EXTERNAL FUNCTION ERATOS
DECLARE NUMERIC L
LET cc=G
FOR I=S TO E STEP ST
LET k2=B(cc)
LET cc=cc+1
LET L=ERATOS(I+ST/4,I+ST/2-1,k2)
SEND TO mes2 FROM L
NEXT I
END PARACT
PARACT PART3
OPTION ARITHMETIC NATIVE
DECLARE NUMERIC B(1 TO 100001)
GET FROM shb TO B
RECEIVE FROM mes6 TO X,Y,Z,G
LET S=X
LET E=Y
LET ST=Z
DECLARE EXTERNAL FUNCTION ERATOS
DECLARE NUMERIC L
LET cc=G
FOR I=S TO E STEP ST
LET k2=B(cc)
LET cc=cc+1
LET L=ERATOS(I+ST/2,I+3*ST/4-1,k2)
SEND TO mes3 FROM L
NEXT I
END PARACT
PARACT PART4
OPTION ARITHMETIC NATIVE
DECLARE NUMERIC B(1 TO 100001)
GET FROM shb TO B
RECEIVE FROM mes7 TO X,Y,Z,G
LET S=X
LET E=Y
LET ST=Z
DECLARE EXTERNAL FUNCTION ERATOS
DECLARE NUMERIC L
LET cc=G
FOR I=S TO E STEP ST
LET k2=B(cc)
LET cc=cc+1
LET L=ERATOS(I+3*ST/4,I+ST-1,k2)
SEND TO mes4 FROM L
NEXT I
END PARACT
EXTERNAL FUNCTION ERATOS(N,M,k2)
OPTION ARITHMETIC NATIVE
DECLARE NUMERIC G(78505) !素数
GET FROM sha TO G
DIM X(0 TO M-N+1)
LET COUNT=0
FOR I=2 TO k2 STEP 1
LET P=G(I)
IF MOD(N,P)=0 THEN
LET NN=N
ELSE
LET NN=INT(N/P+1)*P
END IF
IF N<=I THEN LET NN=NN+P
FOR J=NN TO M STEP P
LET X(J-N)=1
!PRINT i;p;j;n;j-n
NEXT J
NEXT I
FOR J=N TO M STEP 2
IF X(J-N)=0 THEN
LET COUNT=COUNT+1
END IF
NEXT J
LET ERATOS=COUNT
SIGNAL Ok7
END FUNCTION
EXTERNAL SUB prime(k,k2)
OPTION ARITHMETIC NATIVE
DECLARE NUMERIC G(78505) !素数
!エラトステネスの篩
LET Fu=1013 !5633
LET Fm=170 !739
DIM P(Fu)
DIM A(Fm)
MAT P=ZER
MAT A=ZER
LET A(1)=2
LET H1=1
FOR I=3 TO SQR(Fu) STEP 2
IF P(I)=0 THEN
FOR J=I*I TO Fu STEP I
LET P(J)=1
NEXT J
END IF
NEXT I
FOR I=3 TO Fu STEP 2
IF P(I)=0 THEN
LET H1=H1+1
LET A(H1)=I
END IF
NEXT I
LET Q=6
LET k7=k !篩の計算範囲
LET k5=IP(k7/Q)+1
DIM Au(k5),Av(k5)
MAT Au = ZER !(6*n-1)
MAT Av = ZER !(6*n+1)
FOR n=3 TO Fm
LET Pu=A(n)
IF MOD(Pu+1,Q)=0 THEN !(6*n-1)
LET ru=(Pu+1)/Q
FOR i=1 TO k5
IF Pu*i+ru>k5 THEN EXIT FOR
LET Au(Pu*i+ru)=1
NEXT i
END IF
IF MOD(Pu-1,Q)=0 THEN
LET ru=(Pu-1)/Q
FOR i=1 TO k5
IF Pu*i-ru>k5 THEN EXIT FOR
LET Au(Pu*i-ru)=1
NEXT i
END IF
IF MOD(Pu+1,Q)=0 THEN !(6*n+1)
LET ru=(Pu+1)/Q
FOR i=1 TO k5
IF Pu*i-ru>k5 THEN EXIT FOR
LET Av(Pu*i-ru)=1
NEXT i
END IF
IF MOD(Pu-1,Q)=0 THEN
LET ru=(Pu-1)/Q
FOR i=1 TO k5
IF Pu*i+ru>k5 THEN EXIT FOR
LET Av(Pu*i+ru)=1
NEXT i
END IF
NEXT n
LET G(1)=2
LET G(2)=3
LET cz=2
FOR n=1 TO k5
IF 6*n-1>k7 THEN GOTO 100
IF Au(n)=0 THEN
LET cz=cz+1
LET G(cz)=6*n-1
END IF
100 IF 6*n+1>k7 THEN EXIT FOR
IF Av(n)=0 THEN
LET cz=cz+1
LET G(cz)=6*n+1
END IF
NEXT n
PUT TO sha FROM G
SIGNAL Ok5
END SUB
EXTERNAL SUB sqprime(SA,SB,k2)
OPTION ARITHMETIC NATIVE
DECLARE NUMERIC G(78505) !素数
GET FROM sha TO G
DECLARE NUMERIC B(1 TO 100001)
LET x=1E7
LET x1=x
FOR i=1 TO k2-1
LET v=G(i)^2
LET vi=G(i+1)^2
LET vv=vi-v
IF x =< v THEN
CALL sqrp(1)
IF vv>=x1 THEN
LET vt=INT(vv/x1)
CALL sqrp(vt)
END IF
END IF
NEXT i
SUB sqrp(vo)
FOR ns=1 TO vo
LET C=C+1
LET X=X+x1
IF C>=SA AND SB>=C THEN
LET B(c)=i+1
END IF
NEXT ns
END SUB
PUT TO shb FROM B
SIGNAL Ok6
END SUB
http://blogs.yahoo.co.jp/donald_stinger
|
|