Mercurial > hg > Game > Cerium
view example/Miller_Rabin/spe/Prime.cc @ 1244:cd50c48f45e7 draft real_matrix
fix
author | Daichi Toma <amothic@gmail.com> |
---|---|
date | Thu, 03 Nov 2011 20:40:17 +0900 |
parents | 9d37fa6bc1da |
children | 163207b736c5 |
line wrap: on
line source
#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; /* 判定する範囲 */ U64 index_range = range >> 1; /* 判定結果を収める配列を受け取る */ bool *output = (bool*)smanager->get_output(wbuf, 0); /* 初期化 */ for (U64 i = 0; i < index_range; i++){ output[i] = false; } for (U64 i = start + 1,index = 0; i < end ; i += 2, index++) { if (isPrime(i)) { output[index] = true; } } return 0; }