[JavaScript] JavaScriptで円周率を計算したい！

2017/12/23

こんにちは。きんくまです。

今回のお題は「JavaScriptで円周率を計算したい！」です。

ネットで調べてみると円周率計算には様々な方法があることがわかりました。

例）
>> 円周率 – Wikipedia

ガウス＝ルジャンドルのアルゴリズム

それで、今回はいろいろと調べた結果、実装方法が簡単でパフォーマンスもそこそこ出る「ガウス＝ルジャンドルのアルゴリズム」を使ってみることにしました。

>> ガウス＝ルジャンドルのアルゴリズム – Wikipedia

見ていただくとわかると思うのですが、非常に計算式が簡単です。

ただ、これを普通に実装するだけだと何回ステップ数を重ねても14桁ぐらいにしかなりません！
はっきりと理由はわからないのですが、たぶんJavaScriptのNumberの精度の問題なんじゃないかと思います。

そこで計算ライブラリを使います。

decimal.js

使ったのはこのライブラリです。

>> GitHub – MikeMcl/decimal.js: An arbitrary-precision Decimal type for JavaScript

これを使って計算精度を出してあげることができました。ありがたや。

円周率計算クラスを書く

というわけで円周率計算クラスを書きました。

pi_calclator.js

//ガウス＝ルジャンドルのアルゴリズム
//https://ja.wikipedia.org/wiki/%E3%82%AC%E3%82%A6%E3%82%B9%EF%BC%9D%E3%83%AB%E3%82%B8%E3%83%A3%E3%83%B3%E3%83%89%E3%83%AB%E3%81%AE%E3%82%A2%E3%83%AB%E3%82%B4%E3%83%AA%E3%82%BA%E3%83%A0

//正解のPIのページ
//http://www.angio.net/pi/piquery

//decimal.js
// https://github.com/MikeMcl/decimal.js


/**
 * 計算用データ
 * @constructor
 */
function CalcData(){

    var a0 = new Decimal(1.0);
    var b0 = (new Decimal(1.0)).dividedBy((new Decimal(2.0)).sqrt());
    var t0 = (new Decimal(1.0)).dividedBy(new Decimal(4.0));
    var p0 = new Decimal(1.0);

    this.a = a0;
    this.b = b0;
    this.t = t0;
    this.p = p0;
    this.n = 0;
}

/**
 * 円周率計算クラス
 * @param precision 有効桁数
 * @param step ステップ数
 * 数が多いほど精度が高くなる。しかしその分計算時間がかかるようになる。
 * @constructor
 */
var PiCalculator = function(precision, step){
    if(!precision || !step){
        throw new Error("precision or step must be set!!");
    }
    this.precision = precision;
    this.step = step;
}

/**
 * 円周率を計算する（最大10万桁）
 * @param precision 有効桁数
 */
PiCalculator.prototype.start = function(){
    Decimal.set({ precision: this.precision, rounding: 1 });

    var data = new CalcData();
    //console.log(data);
    for(var i = 0, len = this.step; i < len; i++){
        data = this.stepCalc(data);
        console.log('step', i, 'ended');
    }
    var pi = this.calcPi(data);

    //最大10万までの精度で返す
    return pi.toPrecision(100000);
}

PiCalculator.prototype.stepCalc = function(data)
{
    //現在の状態を仮保存
    var tempA = data.a;
    var tempB = data.b;
    var tempT = data.t;
    var tempP = data.p;

    data.a = this.calcA(tempA, tempB);
    data.b = this.calcB(tempA, tempB);
    data.t = this.calcT(tempA, data.a, tempT, tempP);
    data.p = this.calcP(tempP);

    data.n += 1;

    return data;
}

PiCalculator.prototype.calcPi = function(data)
{
    var a_plus_b = data.a.plus(data.b);
    var temp = (new Decimal(4.0)).times(data.t);
    return a_plus_b.times(a_plus_b).dividedBy(temp);
}

PiCalculator.prototype.calcA = function(a_n, b_n)
{
    return (a_n.plus(b_n)).times(0.5);
}

PiCalculator.prototype.calcB = function(a_n, b_n)
{
    return (a_n.times(b_n)).sqrt();
}

PiCalculator.prototype.calcT = function(a_n, a_n1, t_n, p_n)
{
    var tempA = a_n.minus(a_n1);
    return t_n.minus(p_n.times(tempA).times(tempA));
}

PiCalculator.prototype.calcP = function(p_n)
{
    return (new Decimal(2.0)).times(p_n);
}

/**
 * 数値が書いてあるテキストを比べる
 * @param text1 数値テキスト1
 * @param text2 数値テキスト2
 * @return 合っている最大文字数
 */
PiCalculator.prototype.compareNumberText = function(text1, text2)
{
    var text1Arr = text1.split('');
    var maxEqualDigit = 0;

    var text2Arr = text2.split('');
    for(var i = 0, len = text2Arr.length; i < len ; i++){
        if(i > text1Arr.length - 1){
            maxEqualDigit = i + 1;
            break;
        }
        var char = text1Arr[i];
        if(text2Arr[i] != char){
            maxEqualDigit = i + 1;
            break;
        }
    }
    return maxEqualDigit;
}

PiCalculatorというクラスで、コンストラクタに２つの引数をとります

/**
 * 円周率計算クラス
 * @param precision 有効桁数
 * @param step ステップ数
 * 数が多いほど精度が高くなる。しかしその分計算時間がかかるようになる。
 * @constructor
 */
var PiCalculator = function(precision, step){

有効桁数はDecimalライブラリで使う数値の精度です。計算ステップ数をどんなにあげてもこれより大きな桁数の精度は求められません。

ステップ数はガウス＝ルジャンドルのアルゴリズムで使う計算ステップ数です。ステップ数が小さいときは精度は出ません。大きくすると精度は出ますが、その分計算時間が増えます。

使ってみる

html内のjsから使ってみます。
ここでは、decimal.jsじゃなくて軽量版のdecimal.js-lightを使っていますが、２つのライブラリには互換性があるので上のソースは変えなくてもできました。

pi_answer_100000.txtファイルというのは、ネット上に公開されている正解の円周率をテキストファイルにコピペしたものになります。
今回は10万桁のものを保存しました。なのでこの場合、計算結果が合っているかのチェックは最大で10万桁になります。
>> The Pi-Search Page

<!DOCTYPE html>
<html lang="ja">
<head>
    <meta charset="UTF-8">
    <title>円周率計算</title>
</head>
<body>
<script src="jquery-3.1.1.min.js"></script>
<!--<script src="decimal.min.js"></script>-->
<script src="decimal-light.js"></script>
<script src="pi_calclator.js"></script>

<script>

    $(function(){
        main();
    });

    function main(){

        var piCalc = new PiCalculator(10000, 12);

        var t0 = performance.now();

        //円周率を計算
        var resultText = piCalc.start(100);

        var t1 = performance.now();
        //計算時間を出力
        console.log('calculate time = ', t1 - t0, 'ms');

        //何桁まで合っているかチェックする
        loadAnswerText(function(answerText){
            var maxEqualDigit = piCalc.compareNumberText(answerText, resultText);
            console.log('max digit = ', maxEqualDigit);
        });
    }

    function loadAnswerText(callback){
        $.get('pi_answer_100000.txt', function(answerText){
            //console.log(answerText);
            if(callback){
                callback(answerText);
            }
        });
    }
</script>
</body>
</html>

このデモでは、1万桁を12ステップで計算してみました。
出力結果をみてみます。

step 0 ended
step 1 ended
step 2 ended
step 3 ended
step 4 ended
step 5 ended
step 6 ended
step 7 ended
step 8 ended
step 9 ended
step 10 ended
step 11 ended
calculate time =  8236.519999999999 ms
max digit =  10001

1万桁で8.2秒でした。
ただし、ここから試しに3万桁とか5万桁とかやると、とたんに計算時間が長くなります！3分とか5分とか、、、。oh、、。

円周率をCで計算するのはネットに載っていたのですが、高い精度をJavaScriptでやるのをどうやるのかわからなかったので、今回作ってみました。
JSはCと比べてパフォーマンスは出ないのですが、手軽にできるのが良いのではないかと思います。

ではでは。