平方根の計算

平方根の計算には「ニュートン・ラプソン」の方法と言う有名な解法があります。
収束が早く数回ループするだけで高精度の解が求まります。
昔懐かしい汎用コンピュータ & ＭＳ－ＤＯＳ時代の面白そうなプログラムを拾い出してみました。

平方根を求める

x を平方根を求める元の値(２～９)とします。
x1 に初期値として x を設定して、次の計算を繰り返すだけで高精度で解が求まります。

xw= x/x1;
x1= (x1+xw)/2.0;

ちなみに学生時代に覚えた平方根の暗記は次のとおりです。

値	Sqrt(N)	暗記方法	文章表現
２	141421356	ひとよひとよにひとみごろ	一夜一夜に人見ごろ
３	17320508	ひとなみにおごれや	人並みのおごれや
４	２
５	22360679	ふじさんろくおおむなく	富士山麓オーム鳴く
６	2449489	によよくよく	似よ良く良く
７	264575	なにむしいない	菜に虫いない
８	2828	にやにや	にやにや
９	３
10	3162277	みいろにならぶや	三色に並ぶや

プロジェクトの作成

新規プロジェクトから、空の Console Application を作成して、次のファイルをプロジェクトに加えて下さい。
CLI のアプリケーションと通常のアプリケーションでは、作成されるプロジェクトが違います。
CLI Code と Native Code の違いは平方根の計算を参照して下さい。

２～１０の平方根を計算するソースコードです。

    /*★ ニュートンラプソンで平方根を計算する     前田  稔 ★*/
    #include <stdio.h>

    // Function Prototype
    double  Sqrt(double x);

    //★ MAIN PROGRAM
    int  main()
    {   double  x;

        for(x=2.0; x<10.0; x++) Sqrt(x);
        return(0);
    }

    // SQRT() 関数
    double  Sqrt(double x)
    {   double  x1,xw;
        int     i;

        printf("\nSQRT(x=%f)\n",x);
        x1= x;
        for(i=0; i<5; i++)
        {   xw= x/x1;
            printf("I=%d    x1=%f    x/x1=%f",i,x1,xw);
            x1= (x1+xw)/2.0;
            printf("    ANS(%f)\n",x1);
        }
        return(x1);
    }

末尾の桁が合わないのは、印字の際に丸められるからで計算は正確に行われています。

SQRT(x=2.000000)
I=0    x1=2.000000    x/x1=1.000000    ANS(1.500000)
I=1    x1=1.500000    x/x1=1.333333    ANS(1.416667)
I=2    x1=1.416667    x/x1=1.411765    ANS(1.414216)
I=3    x1=1.414216    x/x1=1.414211    ANS(1.414214)
I=4    x1=1.414214    x/x1=1.414214    ANS(1.414214)

SQRT(x=3.000000)
I=0    x1=3.000000    x/x1=1.000000    ANS(2.000000)
I=1    x1=2.000000    x/x1=1.500000    ANS(1.750000)
I=2    x1=1.750000    x/x1=1.714286    ANS(1.732143)
I=3    x1=1.732143    x/x1=1.731959    ANS(1.732051)
I=4    x1=1.732051    x/x1=1.732051    ANS(1.732051)

SQRT(x=4.000000)
I=0    x1=4.000000    x/x1=1.000000    ANS(2.500000)
I=1    x1=2.500000    x/x1=1.600000    ANS(2.050000)
I=2    x1=2.050000    x/x1=1.951220    ANS(2.000610)
I=3    x1=2.000610    x/x1=1.999390    ANS(2.000000)
I=4    x1=2.000000    x/x1=2.000000    ANS(2.000000)

SQRT(x=5.000000)
I=0    x1=5.000000    x/x1=1.000000    ANS(3.000000)
I=1    x1=3.000000    x/x1=1.666667    ANS(2.333333)
I=2    x1=2.333333    x/x1=2.142857    ANS(2.238095)
I=3    x1=2.238095    x/x1=2.234043    ANS(2.236069)
I=4    x1=2.236069    x/x1=2.236067    ANS(2.236068)

SQRT(x=6.000000)
I=0    x1=6.000000    x/x1=1.000000    ANS(3.500000)
I=1    x1=3.500000    x/x1=1.714286    ANS(2.607143)
I=2    x1=2.607143    x/x1=2.301370    ANS(2.454256)
I=3    x1=2.454256    x/x1=2.444732    ANS(2.449494)
I=4    x1=2.449494    x/x1=2.449485    ANS(2.449490)

SQRT(x=7.000000)
I=0    x1=7.000000    x/x1=1.000000    ANS(4.000000)
I=1    x1=4.000000    x/x1=1.750000    ANS(2.875000)
I=2    x1=2.875000    x/x1=2.434783    ANS(2.654891)
I=3    x1=2.654891    x/x1=2.636643    ANS(2.645767)
I=4    x1=2.645767    x/x1=2.645736    ANS(2.645751)

SQRT(x=8.000000)
I=0    x1=8.000000    x/x1=1.000000    ANS(4.500000)
I=1    x1=4.500000    x/x1=1.777778    ANS(3.138889)
I=2    x1=3.138889    x/x1=2.548673    ANS(2.843781)
I=3    x1=2.843781    x/x1=2.813156    ANS(2.828469)
I=4    x1=2.828469    x/x1=2.828386    ANS(2.828427)

SQRT(x=9.000000)
I=0    x1=9.000000    x/x1=1.000000    ANS(5.000000)
I=1    x1=5.000000    x/x1=1.800000    ANS(3.400000)
I=2    x1=3.400000    x/x1=2.647059    ANS(3.023529)
I=3    x1=3.023529    x/x1=2.976654    ANS(3.000092)
I=4    x1=3.000092    x/x1=2.999908    ANS(3.000000)

SQRT(x=10.000000)
I=0    x1=10.000000    x/x1=1.000000    ANS(5.500000)
I=1    x1=5.500000    x/x1=1.818182    ANS(3.659091)
I=2    x1=3.659091    x/x1=2.732919    ANS(3.196005)
I=3    x1=3.196005    x/x1=3.128906    ANS(3.162456)
I=4    x1=3.162456    x/x1=3.162100    ANS(3.162278)

超初心者のプログラム入門(C/C++)

※・