diff example/Miller_Rabin/ppe/Prime.cc @ 1242:9d37fa6bc1da draft

add Miller_Rabin
author Daichi Toma <amothic@gmail.com>
date Tue, 01 Nov 2011 19:01:13 +0900
parents
children cd50c48f45e7
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/example/Miller_Rabin/ppe/Prime.cc	Tue Nov 01 19:01:13 2011 +0900
@@ -0,0 +1,80 @@
+#include <stdio.h>
+#include "SchedTask.h"
+#include "Prime.h"
+#include "Func.h"
+
+typedef unsigned long long U64;
+
+static U64 base[] = { 2, 3, 5, 7, 11 };
+#define MRMAX 5
+
+SchedDefineTask1(Prime, prime);
+
+static U64
+powMod(U64 base, U64 power, U64 mod)
+{
+	U64 result = 1;
+
+	while (power > 0) {
+		if ( power & 1 == 1 ) {
+			result = (result * base) % mod;
+		}
+		base = ( base * base ) % mod;
+		power >>= 1;
+	}
+	return result;
+}
+
+static bool
+millarRabin(U64 a, U64 s, U64 d, U64 n) {
+	U64 x = powMod(a, d, n);
+	if ( x == 1) return true;
+	for ( U64 r = 0; r < s; r++) {
+		if ( x == n-1 ) return true;
+		x = (x * x) % n;
+	}
+	return false;
+}
+
+static bool
+isPrime(U64 n) {
+
+	if ( n <= 1 ) return false;
+
+	U64 d = n - 1;
+	U64 s = 0;
+
+	while( (d & 1) == 0 ) {
+		d >>= 1;
+		s++;
+	}
+
+	for ( int i = 0; i < MRMAX; i++ ) {
+		if ( n == base[i] ) return true;
+		if (!millarRabin(base[i], s, d, n)) return false;
+	}
+	return true;
+}
+
+	static int
+prime(SchedTask *smanager, void *rbuf, void *wbuf)
+{
+	U64 start = (U64)smanager->get_param(0);	/* 素数判定の開始地点 */
+	U64 end   = (U64)smanager->get_param(1);	/* 素数判定の終了地点 */
+	U64 range = end - start;					/* 判定する範囲		  */
+
+	/* 判定結果を収める配列を受け取る */
+	bool *output = (bool*)smanager->get_output(wbuf, 0);
+
+	/* 初期化 */
+	for (U64 i = 0; i < range; i++){
+		output[i] = true;
+	}
+
+	for (U64 i = start + 1,index = 0; index < range; i += 2, index++) {
+		if (!isPrime(i)) {
+			output[index]  =  false;
+		}
+	}
+	return 0;
+}