CODE:

#!/usr/bin/perl

use strict;
use Math::BigInt;
use Math::BigInt::Random qw/ random_bigint /;

my($result);
if(@ARGV==0 || @ARGV>2) {
	die("Usage: miraarabinntesuto.pl  [loop count]\n");
} elsif(@ARGV>=2) {
	$result=&miraarabinntesuto($ARGV[0],$ARGV[1]);
} else {
	$result=&miraarabinntesuto($ARGV[0]);
}
print $result?"Prime!\n":"Not prime...\n";
exit(0);

# http://d.hatena.ne.jp/pashango_p/20090704/1246692091

# if the number is prime, return 1. otherwise, return 0.
sub miraarabinntesuto {
	my($q,$loopnum)=@_;
	my($s,$d,$a,$r,$i,$j,$ok);
	my($powresult,$pownow,$powcurrent);
	if(!$loopnum){$loopnum=50;}
	$q=Math::BigInt->new($q);
	if($qnew("2");
		$a=random_bigint( min => "1", max => $q-1);

		# 「pow(a,d,q) != 1」がTrueだった場合
		$powresult=Math::BigInt->new("1");
		$pownow=$a;
		$powcurrent=$d;
		while($powcurrent>0) {
			if($powcurrent & 1){$powresult=($powresult*$pownow)%$q;}
			$pownow=($pownow*$pownow)%$q;
			$powcurrent>>=1;
		}
		if($powresult==1) {next;}

		# 0～s-1までのカウンタをrとし、すべての「pow(a, 2**(r*d), d)!=-1」がTrue
		# もしかして：pow(a, d*(2**r), q)!=-1
		$ok=1;
		# powresultは前の結果を利用
		for($r=0;$r [ループ回数]