ゆっちのBlog

Java

素数を求める vol.3

Java

素数を求めるプログラムをシンプルな方法で作ってきました。

総当たり方式のものをシングルスレッドとマルチスレッドのプログラムで。

３億までの素数を検出表示させるプログラムです。

実行結果はシングルスレッドのものが約４時間半、

マルチスレッドのものが約１時間半でした。

ここで最終兵器の登場です。

エラトステネスの篩

これはどういったものかと言うと素数の性質を利用して次々と篩い落としていく感じのものです。

例えば 20 までの素数を求める場合

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

最初の 1 は素数の定義から外れるので素数ではない。

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

次の 2 は当然素数の最小値である。

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

と言うことは 2 の倍数である 4 6 8 10 12 14 16 18 20 は素数ではない！

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

次の 3 は消されてないから（1 と自分自身の 3 しか割り切れないから）素数である。

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

同様に 3 の倍数は素数では無いので消していく。9 と 15 は素数ではない。

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

次に現れる 5 は３番目の素数となる。

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

同様に 5 の倍数を消していく。と言っても 5 の倍数は残ってない。

ここで 20 までと言う条件から 5 は 20 の平方根より大きいのでこれで終了になり、

残りの 7 11 13 17 19 は素数となる。

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

実にシンプルで解りやすいですね。

これをプログラムにすればいいだけです。

私が最終兵器と言ったのは処理速度が先の二つのプログラムのものより異常に速いからです。

結果から先に紹介すると３億までの素数を求めるのに１時間もかからないんです。

run:
2 is Prime Number. —> 1
3 is Prime Number. —> 2
5 is Prime Number. —> 3
7 is Prime Number. —> 4
11 is Prime Number. —> 5

　　　　　略

299999939 is Prime Number. —> 16252323
299999957 is Prime Number. —> 16252324
299999977 is Prime Number. —> 16252325
16252325個の素数を検出しました。
0時間37分10秒6889942880002309
構築成功 (合計時間: 37 分 11 秒)

今までのプログラムもこのプログラムも標準出力にすべての素数と何番目か、

それに総数と処理時間を出力していました。

ここで素数を調べるターゲットの数を非常に大きなものにし、標準出力へは

総数と処理時間を表示させるだけのものとしてパフォーマンスの違いを調べてみましょう。

では 2,147,483,640 までの素数を検出させてみます。

まず、前回の並行処理プログラムです。

jp\yucchi\primenumber\PrimeNumber.java

package jp.yucchi.primenumber;

import java.lang.management.ManagementFactory;
import java.util.ArrayList;
import java.util.List;
import java.util.concurrent.*;
import java.util.logging.Level;
import java.util.logging.Logger;

public class PrimeNumber implements Callable<List<Long>> {

    private static final long TARGET_NUMBER = 2_147_483_640L; // 8 の倍数を指定すること(；´Д｀)
    private final long from;
    private final long to;

    public PrimeNumber(final long from, final long to) {
        this.from = from;
        this.to = to;
    }

    public static void main(String[] args) {
        long startTime = System.nanoTime();
        long primeCounter = 0L;
        final int procs = ManagementFactory.getOperatingSystemMXBean().getAvailableProcessors();
        ExecutorService executor = Executors.newFixedThreadPool(procs);
        final long range = TARGET_NUMBER / procs;
        List<Future<List<Long>>> futures = new ArrayList<>();
        for (int i = 0; i < procs; i++) {
            final long from = i * range + 1;
            final long to = (i + 1) * range;
            futures.add(executor.submit(new PrimeNumber(from, to)));
        }
        executor.shutdown();
        List<Long> totalPrimes = new ArrayList<>();
        for (Future<List<Long>> future : futures) {
            try {
                List<Long> primes = future.get();
                totalPrimes.addAll(primes);
            } catch (InterruptedException | ExecutionException ex) {
                Logger.getLogger(PrimeNumber.class.getName()).log(Level.SEVERE, null, ex);
            }
        }
//        for (Long prime : totalPrimes) {
//            primeCounter++;
//            System.out.println(prime + " is Prime Number." + " ---> " + primeCounter);
//        }
//        System.out.println(primeCounter + "個の素数を検出しました。");
        System.out.println(totalPrimes.size());
        long time = System.nanoTime() - startTime;
        System.out.println((int) (time * 1e-9) / 3_600 + "時間"
                + (int) ((time * 1e-9) / 60) % 60 + "分"
                + (int) (time * 1e-9 % 60) + "秒"
                + Double.toString((time * 1e-9 % 60) % 1).substring(2));
    }

    @Override
    public List<Long> call() throws Exception {
        List<Long> primes = new ArrayList<>();
        for (long i = from; i < to + 1L; i++) {
            // 2 より大きな偶数は処理をスキップ
            if ((i & 0b1) == 0 && i > 0b10L) {
                continue;
            }
            long j;
            for (j = (long) Math.sqrt(i); i % j != 0; j--) {
                // i を割り切る j を見つけ出す処理
            }
            // i を割り切れず j が 1 の場合素数と判定、ただし i が 1 の場合は除く        
            if (j == 1L && i != 1L) {
                primes.add(i);
            }
        }
        return primes;
    }
}

package jp.yucchi.primenumber;

import java.lang.management.ManagementFactory;

import java.util.ArrayList;

import java.util.List;

import java.util.concurrent.*;

import java.util.logging.Level;

import java.util.logging.Logger;

public class PrimeNumber implements Callable<List<Long>> {

private static final long TARGET_NUMBER = 2_147_483_640L; // 8 の倍数を指定すること(；´Д｀)

private final long from;

private final long to;

public PrimeNumber(final long from, final long to) {

this.from = from;

this.to = to;

}

public static void main(String[] args) {

long startTime = System.nanoTime();

long primeCounter = 0L;

final int procs = ManagementFactory.getOperatingSystemMXBean().getAvailableProcessors();

ExecutorService executor = Executors.newFixedThreadPool(procs);

final long range = TARGET_NUMBER / procs;

List<Future<List<Long>>> futures = new ArrayList<>();

for (int i = 0; i < procs; i++) {

final long from = i * range + 1;

final long to = (i + 1) * range;

futures.add(executor.submit(new PrimeNumber(from, to)));

}

executor.shutdown();

List<Long> totalPrimes = new ArrayList<>();

for (Future<List<Long>> future : futures) {

try {

List<Long> primes = future.get();

totalPrimes.addAll(primes);

} catch (InterruptedException | ExecutionException ex) {

Logger.getLogger(PrimeNumber.class.getName()).log(Level.SEVERE, null, ex);

}

// for (Long prime : totalPrimes) {

// primeCounter++;

// System.out.println(prime + " is Prime Number." + " ---> " + primeCounter);

// }

// System.out.println(primeCounter + "個の素数を検出しました。");

System.out.println(totalPrimes.size());

long time = System.nanoTime() - startTime;

System.out.println((int) (time * 1e-9) / 3_600 + "時間"

+ (int) ((time * 1e-9) / 60) % 60 + "分"

+ (int) (time * 1e-9 % 60) + "秒"

+ Double.toString((time * 1e-9 % 60) % 1).substring(2));

}

@Override

public List<Long> call() throws Exception {

List<Long> primes = new ArrayList<>();

for (long i = from; i < to + 1L; i++) {

// 2 より大きな偶数は処理をスキップ

if ((i & 0b1) == 0 && i > 0b10L) {

continue;

}

long j;

for (j = (long) Math.sqrt(i); i % j != 0; j--) {

// i を割り切る j を見つけ出す処理

}

// i を割り切れず j が 1 の場合素数と判定、ただし i が 1 の場合は除く

if (j == 1L && i != 1L) {

primes.add(i);

}

return primes;

}

実行結果は次のようにかなりの時間を必要としました。( 15 時間半(；´Д｀) 電気代が・・・）

run:
105097564
15時間26分32秒6862315460020909
構築成功 (合計時間: 926 分 46 秒)

では、最終兵器のエラトステネスの篩ではどうでしょうか？

jp\yucchi\primenumber_eratosthenes\PrimeNumber_Eratosthenes.java

package jp.yucchi.primenumber_eratosthenes;

public class PrimeNumber_Eratosthenes {

    private static final int TARGET_NUMBER = 2_147_483_640;

    public static void main(String[] args) {
        long startTime = System.nanoTime();
        int primeCounter = 0;
        boolean primes[] = new boolean[TARGET_NUMBER];
        for (int i = 2; i < TARGET_NUMBER; i++) {
            // true は素数
            primes[i] = true;
        }
        for (int i = 2; i <= (int) Math.sqrt(TARGET_NUMBER); i++) {
            if (primes[i]) {
                // 素数の倍数を篩にかける
                for (int j = i * 2; 1 < j && j < TARGET_NUMBER; j += i) {
                    primes[j] = false;
                }
            }
        }
        for (int i = 2; i < primes.length; i++) {
            if (primes[i]) {
                primeCounter++;
//                System.out.println(i + " is Prime Number." + " ---> " + primeCounter);
            }
        }
        System.out.println(primeCounter + "個の素数を検出しました。");
        long time = System.nanoTime() - startTime;
        System.out.println((int) (time * 1e-9) / 3_600 + "時間"
                + (int) ((time * 1e-9) / 60) % 60 + "分"
                + (int) (time * 1e-9 % 60) + "秒"
                + Double.toString((time * 1e-9 % 60) % 1).substring(2));
    }
}

package jp.yucchi.primenumber_eratosthenes;

public class PrimeNumber_Eratosthenes {

private static final int TARGET_NUMBER = 2_147_483_640;

public static void main(String[] args) {

long startTime = System.nanoTime();

int primeCounter = 0;

boolean primes[] = new boolean[TARGET_NUMBER];

for (int i = 2; i < TARGET_NUMBER; i++) {

// true は素数

primes[i] = true;

}

for (int i = 2; i <= (int) Math.sqrt(TARGET_NUMBER); i++) {

if (primes[i]) {

// 素数の倍数を篩にかける

for (int j = i * 2; 1 < j && j < TARGET_NUMBER; j += i) {

primes[j] = false;

}

for (int i = 2; i < primes.length; i++) {

if (primes[i]) {

primeCounter++;

// System.out.println(i + " is Prime Number." + " ---> " + primeCounter);

}

System.out.println(primeCounter + "個の素数を検出しました。");

long time = System.nanoTime() - startTime;

System.out.println((int) (time * 1e-9) / 3_600 + "時間"

+ (int) ((time * 1e-9) / 60) % 60 + "分"

+ (int) (time * 1e-9 % 60) + "秒"

+ Double.toString((time * 1e-9 % 60) % 1).substring(2));

}

さて、実行結果はどうなったでしょうか。

run:
105097564個の素数を検出しました。
0時間0分32秒7029516360000017
構築成功 (合計時間: 32 秒)

速っ！＼(◎o◎)／！

この違いは何なんだ！

凄いぜ！　エラトステネスの篩

Hatena タグ: Java

TAGS: Java | 2012年7月5日2:45 PM | Comment : 0

素数を求める vol.2

Java

昨日、素数を求めるプログラムを作ったが 3 億までの素数を求めるのに約 4 時間半かかった。

これが速いのか遅いのかは置いといて、もう少し処理時間の短縮をはかりたい。

最近の PC はマルチコアプロセッサを搭載しているので並列処理させれば処理時間短縮が可能かもしれない。

幸いにも私の PC は 8 コアなので早速試してみました。

ちなみに下記プログラムは、java.util.concurrent パッケージについて調べていた時にあるサイトにあった

サンプルプログラムを参考にさせていただきました。

jp\yucchi\primenumber\PrimeNumber.java

package jp.yucchi.primenumber;

import java.lang.management.ManagementFactory;
import java.util.ArrayList;
import java.util.List;
import java.util.concurrent.*;
import java.util.logging.Level;
import java.util.logging.Logger;

public class PrimeNumber implements Callable<List<Long>> {

    private static final long TARGET_NUMBER = 300_000_000L; // 8 の倍数を指定すること(；´Д｀)
    private final long from;
    private final long to;

    public PrimeNumber(final long from, final long to) {
        this.from = from;
        this.to = to;
    }

    public static void main(String[] args) {
        long startTime = System.nanoTime();
        long primeCounter = 0L;
        final int procs = ManagementFactory.getOperatingSystemMXBean().getAvailableProcessors();
        ExecutorService executor = Executors.newFixedThreadPool(procs);
        final long range = TARGET_NUMBER / procs;
        List<Future<List<Long>>> futures = new ArrayList<>();
        for (int i = 0; i < procs; i++) {
            final long from = i * range + 1;
            final long to = (i + 1) * range;
            futures.add(executor.submit(new PrimeNumber(from, to)));
        }
        executor.shutdown();
        List<Long> totalPrimes = new ArrayList<>();
        for (Future<List<Long>> future : futures) {
            try {
                List<Long> primes = future.get();
                totalPrimes.addAll(primes);
            } catch (InterruptedException | ExecutionException ex) {
                Logger.getLogger(PrimeNumber.class.getName()).log(Level.SEVERE, null, ex);
            }
        }
        for (Long prime : totalPrimes) {
            primeCounter++;
            System.out.println(prime + " is Prime Number." + " ---> " + primeCounter);
        }
        System.out.println(primeCounter + "個の素数を検出しました。");
        long time = System.nanoTime() - startTime;
        System.out.println((int) (time * 1e-9) / 3_600 + "時間"
                + (int) ((time * 1e-9) / 60) % 60 + "分"
                + (int) (time * 1e-9 % 60) + "秒"
                + Double.toString((time * 1e-9 % 60) % 1).substring(2));
    }

    @Override
    public List<Long> call() throws Exception {
        List<Long> primes = new ArrayList<>();
        for (long i = from; i < to + 1L; i++) {
            // 2 より大きな偶数は処理をスキップ
            if ((i & 0b1) == 0 && i > 0b10L) {
                continue;
            }
            long j;
            for (j = (long) Math.sqrt(i); i % j != 0; j--) {
                // i を割り切る j を見つけ出す処理
            }
            // i を割り切れず j が 1 の場合素数と判定、ただし i が 1 の場合は除く        
            if (j == 1L && i != 1L) {
                primes.add(i);
            }
        }
        return primes;
    }
}

package jp.yucchi.primenumber;

import java.lang.management.ManagementFactory;

import java.util.ArrayList;

import java.util.List;

import java.util.concurrent.*;

import java.util.logging.Level;

import java.util.logging.Logger;

public class PrimeNumber implements Callable<List<Long>> {

private static final long TARGET_NUMBER = 300_000_000L; // 8 の倍数を指定すること(；´Д｀)

private final long from;

private final long to;

public PrimeNumber(final long from, final long to) {

this.from = from;

this.to = to;

}

public static void main(String[] args) {

long startTime = System.nanoTime();

long primeCounter = 0L;

final int procs = ManagementFactory.getOperatingSystemMXBean().getAvailableProcessors();

ExecutorService executor = Executors.newFixedThreadPool(procs);

final long range = TARGET_NUMBER / procs;

List<Future<List<Long>>> futures = new ArrayList<>();

for (int i = 0; i < procs; i++) {

final long from = i * range + 1;

final long to = (i + 1) * range;

futures.add(executor.submit(new PrimeNumber(from, to)));

}

executor.shutdown();

List<Long> totalPrimes = new ArrayList<>();

for (Future<List<Long>> future : futures) {

try {

List<Long> primes = future.get();

totalPrimes.addAll(primes);

} catch (InterruptedException | ExecutionException ex) {

Logger.getLogger(PrimeNumber.class.getName()).log(Level.SEVERE, null, ex);

}

for (Long prime : totalPrimes) {

primeCounter++;

System.out.println(prime + " is Prime Number." + " ---> " + primeCounter);

}

System.out.println(primeCounter + "個の素数を検出しました。");

long time = System.nanoTime() - startTime;

System.out.println((int) (time * 1e-9) / 3_600 + "時間"

+ (int) ((time * 1e-9) / 60) % 60 + "分"

+ (int) (time * 1e-9 % 60) + "秒"

+ Double.toString((time * 1e-9 % 60) % 1).substring(2));

}

@Override

public List<Long> call() throws Exception {

List<Long> primes = new ArrayList<>();

for (long i = from; i < to + 1L; i++) {

// 2 より大きな偶数は処理をスキップ

if ((i & 0b1) == 0 && i > 0b10L) {

continue;

}

long j;

for (j = (long) Math.sqrt(i); i % j != 0; j--) {

// i を割り切る j を見つけ出す処理

}

// i を割り切れず j が 1 の場合素数と判定、ただし i が 1 の場合は除く

if (j == 1L && i != 1L) {

primes.add(i);

}

return primes;

}

このプログラムの実行結果は次のようになった。

2 is Prime Number. —> 1
3 is Prime Number. —> 2
5 is Prime Number. —> 3
7 is Prime Number. —> 4
11 is Prime Number. —> 5
13 is Prime Number. —> 6
17 is Prime Number. —> 7
19 is Prime Number. —> 8
23 is Prime Number. —> 9
29 is Prime Number. —> 10

　　　　略

299999801 is Prime Number. —> 16252317
299999807 is Prime Number. —> 16252318
299999813 is Prime Number. —> 16252319
299999827 is Prime Number. —> 16252320
299999897 is Prime Number. —> 16252321
299999923 is Prime Number. —> 16252322
299999939 is Prime Number. —> 16252323
299999957 is Prime Number. —> 16252324
299999977 is Prime Number. —> 16252325
16252325個の素数を検出しました。
1時間26分23秒42106727300051716

なんと 3 時間も処理時間が短縮できました。(^_^)v

それでも 1 時間半近くかかります。

こんなものなのかなぁ・・・？

素数を求めるアルゴリズムにエラトステネスの篩ってのがあります。

それも試してみましょう。

それでは

おっとそろそろ時間が・・・

~~続きは Web で！~~

つづきはまた今度ね♪

Hatena タグ: Java

TAGS: Java | 2012年6月25日8:15 AM | Comment : 0

素数を求める

Java

素数を求めるプログラムを考えてみた。

まず素数とはどのように定義されているか確認した。

1 と自分自身以外に正の約数を持たない、1 でない自然数のことである。

ということは、1 より大きくて自分自身でないと割り切れない自然数であるってことかな。

完全に総当たり方式でいけば楽勝ですね。

でも、2 より大きな偶数は 2 で割り切れるから除外すればいいし、

ターゲット n が合成数（素数でない数）の場合、n=ab となる自然数 a , b が存在する。

a と b を掛けると、n になるのだから、a , b のうちどちらかは n の平方根以下であり、どちらかは n の平方根以上だ。

従って、素数でないかどうかを確認するには、nの平方根以下の自然数で割り切れるかどうか調べればいいはず。

これらの特徴をふまえた上で、 3 億までの素数を求めるプログラムをくんでみた。

jp\yucchi\primenumber_loop\PrimeNumber_Loop.java

package jp.yucchi.primenumber_loop;

public class PrimeNumber_Loop {

    private static final long TARGET_NUMBER = 300_000_000L;

    public static void main(String[] args) {
        long startTime = System.nanoTime();
        long n, i, primeCounter = 0L;
        for (n = 2L; n <= TARGET_NUMBER; n++) {
            // 2 より大きな偶数は処理をスキップ           
            if ((n & 0b1) == 0 && n > 2) {
                continue;
            }
            for (i = (long) Math.sqrt(n); n % i != 0; i--) {
                // n を割り切る i を見つけ出す処理
            }
            if (i == 1) { // n を割り切れず i が 1 の場合は素数と判定
                primeCounter++;
                System.out.println(n + " is Prime Number." + " ---> " + primeCounter);
            }
        }
        System.out.println(primeCounter + "個の素数を検出しました。");
        long time = System.nanoTime() - startTime;
        System.out.println((int) (time * 1e-9) / 3_600 + "時間"
                + (int) ((time * 1e-9) / 60) % 60 + "分"
                + (int) (time * 1e-9 % 60) + "秒"
                + Double.toString((time * 1e-9 % 60) % 1).substring(2));
    }
}

package jp.yucchi.primenumber_loop;

public class PrimeNumber_Loop {

private static final long TARGET_NUMBER = 300_000_000L;

public static void main(String[] args) {

long startTime = System.nanoTime();

long n, i, primeCounter = 0L;

for (n = 2L; n <= TARGET_NUMBER; n++) {

// 2 より大きな偶数は処理をスキップ

if ((n & 0b1) == 0 && n > 2) {

continue;

}

for (i = (long) Math.sqrt(n); n % i != 0; i--) {

// n を割り切る i を見つけ出す処理

}

if (i == 1) { // n を割り切れず i が 1 の場合は素数と判定

primeCounter++;

System.out.println(n + " is Prime Number." + " ---> " + primeCounter);

}

System.out.println(primeCounter + "個の素数を検出しました。");

long time = System.nanoTime() - startTime;

System.out.println((int) (time * 1e-9) / 3_600 + "時間"

+ (int) ((time * 1e-9) / 60) % 60 + "分"

+ (int) (time * 1e-9 % 60) + "秒"

+ Double.toString((time * 1e-9 % 60) % 1).substring(2));

}

このプログラムを実行結果は次のようになった。

299999771 is Prime Number. —> 16252316
299999801 is Prime Number. —> 16252317
299999807 is Prime Number. —> 16252318
299999813 is Prime Number. —> 16252319
299999827 is Prime Number. —> 16252320
299999897 is Prime Number. —> 16252321
299999923 is Prime Number. —> 16252322
299999939 is Prime Number. —> 16252323
299999957 is Prime Number. —> 16252324
299999977 is Prime Number. —> 16252325
16252325個の素数を検出しました。
4時間29分24秒10574628300128097

凄く時間がかかった（×_×）

ちなみに Dual Xeon 3GHz で総コア数 8 です。

ハイパースレッドは無い古い CPU です。

ちょっとでも時間短縮する為に次回は並列処理をさせてみたいと思います。

実は並列処理の方法は java.util.concurrent パッケージについて調べていた時にあるサイトに

サンプルプログラムとしてありました。

それをちょっとだけ変更して使用させていただきます。

と言いつつ睡魔に襲われてきたのでまた今度ってことにします。

それでは、つづきは ~~Web で！~~

おやすみ～～～～( ^.^)( -.-)( _ _)

Hatena タグ: Java

TAGS: Java | 2012年6月24日11:55 PM | Comment : 0

円周率をモンテカルロ法で求める vol.3

Java

先日、モンテカルロ法で円周率を求めるプログラムをくんだ。

ちょっと気になることがあって確かめてみたら驚いた。

そのプログラムを振り返ってみよう。

はじめに二重ループを使ったシンプルなプログラム。

jp\yucchi\loop_motecarlo_pi\Loop_Motecarlo_PI.java

package jp.yucchi.loop_motecarlo_pi;

public class Loop_Motecarlo_PI {

    private static final double TOTAL_THREAD = 16.0;
    private final static double NEEDLE = 1_000_000.0;

    public static void main(String[] args) {
        long stime = System.currentTimeMillis();
        double totalPi = 0.0;

        for (double j = 0.0; j < TOTAL_THREAD; j++) {
            double hitCount = 0.0;
            double x, y, pi;
            for (double i = 0.0; i < NEEDLE; i++) {
                x = Math.random();
                y = Math.random();
                if (Math.pow(x, 2.0) + Math.pow(y, 2.0) < 1.0) {
                    hitCount++;
                }
            }
            pi = 4.0 * hitCount / NEEDLE;
            System.out.println("Montecarlo PI :  " + pi);
            totalPi += pi;
        }
        double avgPi = totalPi / TOTAL_THREAD;
        System.out.println("Ｍａｔｈクラスによる円周率は " + Math.PI);
        System.out.println("モンテカルロ法による円周率は " + avgPi);
        stime = System.currentTimeMillis() - stime;
        System.out.println("Time : " + stime + "ミリ秒");
    }
}

package jp.yucchi.loop_motecarlo_pi;

public class Loop_Motecarlo_PI {

private static final double TOTAL_THREAD = 16.0;

private final static double NEEDLE = 1_000_000.0;

public static void main(String[] args) {

long stime = System.currentTimeMillis();

double totalPi = 0.0;

for (double j = 0.0; j < TOTAL_THREAD; j++) {

double hitCount = 0.0;

double x, y, pi;

for (double i = 0.0; i < NEEDLE; i++) {

x = Math.random();

y = Math.random();

if (Math.pow(x, 2.0) + Math.pow(y, 2.0) < 1.0) {

hitCount++;

}

pi = 4.0 * hitCount / NEEDLE;

System.out.println("Montecarlo PI : " + pi);

totalPi += pi;

}

double avgPi = totalPi / TOTAL_THREAD;

System.out.println("Ｍａｔｈクラスによる円周率は " + Math.PI);

System.out.println("モンテカルロ法による円周率は " + avgPi);

stime = System.currentTimeMillis() - stime;

System.out.println("Time : " + stime + "ミリ秒");

}

その実行結果は次のようになりました。

run:
Montecarlo PI : 3.144016
Montecarlo PI : 3.143836
Montecarlo PI : 3.140812
Montecarlo PI : 3.142616
Montecarlo PI : 3.141592
Montecarlo PI : 3.139224
Montecarlo PI : 3.143332
Montecarlo PI : 3.140608
Montecarlo PI : 3.141352
Montecarlo PI : 3.139116
Montecarlo PI : 3.140224
Montecarlo PI : 3.141696
Montecarlo PI : 3.144996
Montecarlo PI : 3.142304
Montecarlo PI : 3.14044
Montecarlo PI : 3.14482
Ｍａｔｈクラスによる円周率は 3.141592653589793
モンテカルロ法による円周率は 3.141936500000001
Time : 1516ミリ秒
構築成功 (合計時間: 1 秒)

次に今時の PC はマルチコアが当たり前になってきてるので

並行処理をさせるように変更しました。

jp\yucchi\montecarlo_pi\Montecarlo_PI.java

package jp.yucchi.montecarlo_pi;

import java.lang.management.ManagementFactory;
import java.util.ArrayList;
import java.util.List;
import java.util.concurrent.*;
import java.util.logging.Level;
import java.util.logging.Logger;

public class Montecarlo_PI implements Callable<List<Double>> {

    private static final double TOTAL_THREAD = 16.0;
    private final double NEEDLE = 1_000_000.0;

    public static void main(String[] args) {
        long stime = System.currentTimeMillis();
        int procs = ManagementFactory.getOperatingSystemMXBean().getAvailableProcessors();
        ExecutorService executor = Executors.newFixedThreadPool(procs);
        List<Future<List<Double>>> futures = new ArrayList<>();
        for (double i = 0; i < TOTAL_THREAD; i++) {
            futures.add(executor.submit(new Montecarlo_PI()));
        }
        executor.shutdown();
        List<Double> mPi = new ArrayList<>();
        for (Future<List<Double>> future : futures) {
            try {
                List<Double> pi = future.get();
                mPi.addAll(pi);
            } catch (InterruptedException | ExecutionException ex) {
                Logger.getLogger(Montecarlo_PI.class.getName()).log(Level.SEVERE, null, ex);
            }
        }
        double totalPi = 0.0;
        for (Double pi : mPi) {
            totalPi += pi;
        }
        double avgPi = totalPi / TOTAL_THREAD;
        System.out.println("Ｍａｔｈクラスによる円周率は " + Math.PI);
        System.out.println("モンテカルロ法による円周率は " + avgPi);
        stime = System.currentTimeMillis() - stime;
        System.out.println("Time : " + stime + "ミリ秒");
    }

    @Override
    public List<Double> call() throws Exception {
        List<Double> pi = new ArrayList<>();
        double hitCount = 0.0;
        for (double i = 0; i < NEEDLE; i++) {
            double x = Math.random();
            double y = Math.random();
            if (Math.pow(x, 2.0) + Math.pow(y, 2.0) < 1.0) {
                hitCount++;
            }
        }
        double _pi = 4.0 * hitCount / NEEDLE;
        System.out.println("Montecarlo PI : " + _pi);
        pi.add(_pi);
        return pi;
    }
}

package jp.yucchi.montecarlo_pi;

import java.lang.management.ManagementFactory;

import java.util.ArrayList;

import java.util.List;

import java.util.concurrent.*;

import java.util.logging.Level;

import java.util.logging.Logger;

public class Montecarlo_PI implements Callable<List<Double>> {

private static final double TOTAL_THREAD = 16.0;

private final double NEEDLE = 1_000_000.0;

public static void main(String[] args) {

long stime = System.currentTimeMillis();

int procs = ManagementFactory.getOperatingSystemMXBean().getAvailableProcessors();

ExecutorService executor = Executors.newFixedThreadPool(procs);

List<Future<List<Double>>> futures = new ArrayList<>();

for (double i = 0; i < TOTAL_THREAD; i++) {

futures.add(executor.submit(new Montecarlo_PI()));

}

executor.shutdown();

List<Double> mPi = new ArrayList<>();

for (Future<List<Double>> future : futures) {

try {

List<Double> pi = future.get();

mPi.addAll(pi);

} catch (InterruptedException | ExecutionException ex) {

Logger.getLogger(Montecarlo_PI.class.getName()).log(Level.SEVERE, null, ex);

}

double totalPi = 0.0;

for (Double pi : mPi) {

totalPi += pi;

}

double avgPi = totalPi / TOTAL_THREAD;

System.out.println("Ｍａｔｈクラスによる円周率は " + Math.PI);

System.out.println("モンテカルロ法による円周率は " + avgPi);

stime = System.currentTimeMillis() - stime;

System.out.println("Time : " + stime + "ミリ秒");

}

@Override

public List<Double> call() throws Exception {

List<Double> pi = new ArrayList<>();

double hitCount = 0.0;

for (double i = 0; i < NEEDLE; i++) {

double x = Math.random();

double y = Math.random();

if (Math.pow(x, 2.0) + Math.pow(y, 2.0) < 1.0) {

hitCount++;

}

double _pi = 4.0 * hitCount / NEEDLE;

System.out.println("Montecarlo PI : " + _pi);

pi.add(_pi);

return pi;

}

実行結果は先ほどのプログラムを遙かに上回る結果となるはずだったんだけど・・・

run:
Montecarlo PI : 3.143472
Montecarlo PI : 3.141672
Montecarlo PI : 3.1408
Montecarlo PI : 3.140352
Montecarlo PI : 3.141668
Montecarlo PI : 3.13966
Montecarlo PI : 3.139112
Montecarlo PI : 3.142532
Montecarlo PI : 3.13982
Montecarlo PI : 3.143956
Montecarlo PI : 3.141428
Montecarlo PI : 3.138096
Montecarlo PI : 3.141932
Montecarlo PI : 3.143548
Montecarlo PI : 3.14172
Montecarlo PI : 3.141612
Ｍａｔｈクラスによる円周率は 3.141592653589793
モンテカルロ法による円周率は 3.14133625
Time : 22172ミリ秒
構築成功 (合計時間: 22 秒)

(゜◇゜)ガーン

な、なんと　遙かに遅い！

ExecutorService でスレッドプール作って処理させるほうがコストが高くつく！

この程度の計算ならシングルスレッドでもいいみたい。

そこで時間のかかる処理を無意味に追加した場合どうなるか確認した。

次のコードをそれぞれのプログラムに追加した。

for (int k = 0; k < 10_000; ++k) {
    for (int n = 0; n < 10_000; ++n) {
        double p = k * n;
                p = Math.sqrt(p);
    }
}

その結果は次のようになった。

二重ループのプログラムでは

run:
Montecarlo PI : 3.138036
Montecarlo PI : 3.144924
Montecarlo PI : 3.139508
Montecarlo PI : 3.14336
Montecarlo PI : 3.140848
Montecarlo PI : 3.14178
Montecarlo PI : 3.141744
Montecarlo PI : 3.143868
Montecarlo PI : 3.14006
Montecarlo PI : 3.143568
Montecarlo PI : 3.137936
Montecarlo PI : 3.141968
Montecarlo PI : 3.143396
Montecarlo PI : 3.142736
Montecarlo PI : 3.137876
Montecarlo PI : 3.141624
Ｍａｔｈクラスによる円周率は 3.141592653589793
モンテカルロ法による円周率は 3.141452
Time : 250579ミリ秒
構築成功 (合計時間: 4 分 10 秒)

時間が結構かかってますね！

平行処理プログラムでは次のようになりました。

run:
Montecarlo PI : 3.142172
Montecarlo PI : 3.140212
Montecarlo PI : 3.142544
Montecarlo PI : 3.14186
Montecarlo PI : 3.14264
Montecarlo PI : 3.143252
Montecarlo PI : 3.140788
Montecarlo PI : 3.140544
Montecarlo PI : 3.14346
Montecarlo PI : 3.144144
Montecarlo PI : 3.141128
Montecarlo PI : 3.140564
Montecarlo PI : 3.140828
Montecarlo PI : 3.13794
Montecarlo PI : 3.140336
Montecarlo PI : 3.141248
Ｍａｔｈクラスによる円周率は 3.141592653589793
モンテカルロ法による円周率は 3.1414787499999997
Time : 32468ミリ秒
構築成功 (合計時間: 32 秒)

(^_^)v　速い!　

浮動小数点演算なんかをさせるとかなり違いがでるのかな？

なんでもかんでも並行処理させればいいってことじゃないのかなぁ・・・？

リソースモニターですべてのコアが 100% で動いてるのは圧巻だったけど

処理が遅くては本末転倒ですね。(^_^;)

Hatena タグ: Java

TAGS: Java | 2012年6月17日9:02 PM | Comment : 0

円周率をモンテカルロ法で求める vol.2

Java

昨日、モンテカルロ法にて円周率を求めるプログラムをつくった。

乱数シュミレーションとしては、まぁまぁの数値を導き出していると個人的には思った。

さらに精度を上げようとして何回か同じ処理をして平均値を取ればいいかな？と思い

次のようにプログラムを変更してみました。

はじめに二重ループを使い平均値を求めるプログラムを組んでしまって

非常に効率が悪いことに気が付いたのは秘密です。(；´Д｀)

相変わらず汚いコードですが素人のやることですからプロフェッショナルな人は

優しく解りやすくご指摘くださればうれしいです。

jp\yucchi\montecarlo_pi\Montecarlo_PI.java

package jp.yucchi.montecarlo_pi;

import java.lang.management.ManagementFactory;
import java.util.ArrayList;
import java.util.List;
import java.util.concurrent.*;
import java.util.logging.Level;
import java.util.logging.Logger;

public class Montecarlo_PI implements Callable<List<Double>> {

    public static void main(String[] args) {
        double TOTAL_THREAD = 1_000_000.0;
        int procs = ManagementFactory.getOperatingSystemMXBean().getAvailableProcessors();
        ExecutorService executor = Executors.newFixedThreadPool(procs);
        List<Future<List<Double>>> futures = new ArrayList<>();
        for (double i = 0; i < TOTAL_THREAD; i++) {
            futures.add(executor.submit(new Montecarlo_PI()));
        }
        executor.shutdown();
        List<Double> mPi = new ArrayList<>();
        for (Future<List<Double>> future : futures) {
            try {
                List<Double> pi = future.get();
                mPi.addAll(pi);
            } catch (InterruptedException | ExecutionException ex) {
                Logger.getLogger(Montecarlo_PI.class.getName()).log(Level.SEVERE, null, ex);
            }
        }
        double totalPi = 0.0;
        for (Double pi : mPi) {
            totalPi += pi;
        }
        double avgPi = totalPi / mPi.size();
        System.out.println("Ｍａｔｈクラスによる円周率は " + Math.PI);
        System.out.println("モンテカルロ法による円周率は " + avgPi);
    }

    @Override
    public List<Double> call() throws Exception {
        List<Double> pi = new ArrayList<>();
        double hitCount = 0.0;
        double NEEDLE = 100_000_000.0;
        for (double i = 0; i < NEEDLE; i++) {
            double x = Math.random();
            double y = Math.random();
            if (Math.pow(x, 2.0) + Math.pow(y, 2.0) < 1.0) {
                hitCount++;
            }
        }
        double _pi = 4.0 * hitCount / NEEDLE;
        System.out.println("Montecarlo PI : " + _pi);
        pi.add(_pi);
        return pi;
    }
}

package jp.yucchi.montecarlo_pi;

import java.lang.management.ManagementFactory;

import java.util.ArrayList;

import java.util.List;

import java.util.concurrent.*;

import java.util.logging.Level;

import java.util.logging.Logger;

public class Montecarlo_PI implements Callable<List<Double>> {

public static void main(String[] args) {

double TOTAL_THREAD = 1_000_000.0;

int procs = ManagementFactory.getOperatingSystemMXBean().getAvailableProcessors();

ExecutorService executor = Executors.newFixedThreadPool(procs);

List<Future<List<Double>>> futures = new ArrayList<>();

for (double i = 0; i < TOTAL_THREAD; i++) {

futures.add(executor.submit(new Montecarlo_PI()));

}

executor.shutdown();

List<Double> mPi = new ArrayList<>();

for (Future<List<Double>> future : futures) {

try {

List<Double> pi = future.get();

mPi.addAll(pi);

} catch (InterruptedException | ExecutionException ex) {

Logger.getLogger(Montecarlo_PI.class.getName()).log(Level.SEVERE, null, ex);

}

double totalPi = 0.0;

for (Double pi : mPi) {

totalPi += pi;

}

double avgPi = totalPi / mPi.size();

System.out.println("Ｍａｔｈクラスによる円周率は " + Math.PI);

System.out.println("モンテカルロ法による円周率は " + avgPi);

}

@Override

public List<Double> call() throws Exception {

List<Double> pi = new ArrayList<>();

double hitCount = 0.0;

double NEEDLE = 100_000_000.0;

for (double i = 0; i < NEEDLE; i++) {

double x = Math.random();

double y = Math.random();

if (Math.pow(x, 2.0) + Math.pow(y, 2.0) < 1.0) {

hitCount++;

}

double _pi = 4.0 * hitCount / NEEDLE;

System.out.println("Montecarlo PI : " + _pi);

pi.add(_pi);

return pi;

}

一億個の針を落として計算する処理を百万回して平均値をとってみました。

って・・・計算処理が終わってない(￣。￣;)

しかたないから一億個の針を落とす計算処理を千回に変更してみました。

上記のソースコードの NEEDLE の値を 1000.0 に変更して実行してみました。

これでも終わりそうにないくらい時間がかかります。

不本意ですが計算処理回数を３２回に変更しました。(..;)

実行結果は次のようになりほんの気持ち程度精度は上がりました。

run:
Montecarlo PI : 3.14146332
Montecarlo PI : 3.1416326
Montecarlo PI : 3.1415482
Montecarlo PI : 3.141632
Montecarlo PI : 3.14127944
Montecarlo PI : 3.14178148
Montecarlo PI : 3.1416116
Montecarlo PI : 3.14141588
Montecarlo PI : 3.1416632
Montecarlo PI : 3.14117084
Montecarlo PI : 3.14161116
Montecarlo PI : 3.14168484
Montecarlo PI : 3.14166664
Montecarlo PI : 3.14151984
Montecarlo PI : 3.14131844
Montecarlo PI : 3.14161456
Montecarlo PI : 3.14160796
Montecarlo PI : 3.14157636
Montecarlo PI : 3.14174076
Montecarlo PI : 3.14137984
Montecarlo PI : 3.14161576
Montecarlo PI : 3.14160444
Montecarlo PI : 3.14165568
Montecarlo PI : 3.14140136
Montecarlo PI : 3.1414314
Montecarlo PI : 3.14134932
Montecarlo PI : 3.14164312
Montecarlo PI : 3.14159468
Montecarlo PI : 3.141492
Montecarlo PI : 3.14139636
Montecarlo PI : 3.14149788
Montecarlo PI : 3.1414862
Ｍａｔｈクラスによる円周率は 3.141592653589793
モンテカルロ法による円周率は 3.14153397375
構築成功 (合計時間: 22 分 47 秒)

時間に余裕のあるときに計算処理回数を増やしてみて結果を確認したいと思います。

ちなみにこの計算は Intel Core i7-2600K (3.4GHz) の CPU を積んだハードウェアでしました。

OS は、Windows 7 64bit 版です。

Dual Xeon 4 Core (3.0GHz) のシステムより２倍くらい速いです！

ハイエンドマシンをお持ちの方は是非試してみてくださいね。(ё_ё)

実行結果の合計時間: 22 分 47 秒ってのは使用した NetBeans 7.2 開発版での実行時の値です。

Hatena タグ: Java

TAGS: Java | 2012年6月13日12:02 AM | Comment : 0

« 古い記事新しい記事 »

ゆっちのBlog

素数を求める vol.3

素数を求める vol.2

素数を求める

円周率をモンテカルロ法で求める vol.3

円周率をモンテカルロ法で求める vol.2

Recent Comments

Java

NetBeans

カテゴリー

アーカイブ

RSS FEED

ねこび～ん

メタ情報

ゆっちのカレンダー

2024年10月
月	火	水	木	金	土	日
« 1月
	1	2	3	4	5	6
7	8	9	10	11	12	13
14	15	16	17	18	19	20
21	22	23	24	25	26	27
28	29	30	31

ゆっちのBlog

素数を求める vol.3

素数を求める vol.2

素数を求める

円周率をモンテカルロ法で求める vol.3

円周率をモンテカルロ法で求める vol.2

Popular post ranking

Recent Comments

Java

NetBeans

カテゴリー

アーカイブ

RSS FEED

ねこび～ん

メタ情報

ゆっちのカレンダー