【数学】円周率の計算(多角形の円周の長さ+半角の公式)

多角形の円周の長さと半角の公式を使った円周率の計算を思いついたので記載。

本方法で一番難しいのが半角の公式というレベルなので、めちゃめちゃ分かりやすい部類の円周率の求め方なのでは??と思ったり。(収束速度とかは悪いのだろうと思います。。)

まぁおそらく誰かが既に考えているとは思いますが。。

■1. 準備

■1.1. 半角の公式

まず、半角の公式がどんなものだったかというと、教科書やネットなどで簡単に見つかるかと思いますが、下記ですね。

\displaystyle{
\begin{aligned}
\sin \dfrac{\theta }{2}&=\pm \sqrt{\dfrac{1-\cos \theta }{2}}\\
\cos \dfrac{\theta }{2}&=\pm \sqrt{\dfrac{1+\cos \theta }{2}}\\
\tan \dfrac{\theta }{2}&=\dfrac{\sin \dfrac{\theta }{2}}{\cos \dfrac{\theta }{2}}
\end{aligned}
}

これは、ある角度の時のsin,cos,tanが分かっていれば、その半分の角度の時のsin,cos,tanが分かるというものですね。

■1.2. 2003年の東大の入試問題

伝説の入試問題ともされている2003年の東大の入試問題で、円周率>3.05であることの証明があり、有名になっている。

これは、下記のように円の中に八角形を内接させて、八角形の周の長さと円周の長さを比較すればよい。(もちろん下記以外の証明方法は多数ありますが、本証明方法が一番簡単で一番有名なのではと思います。)

  

半径r=1とすれば、余弦定理を使って、詳細な計算をすっ飛ばして、下記となり、 円周率が3.05より大きいことが分かるのです。

\displaystyle{
\begin{aligned}8 \sqrt{2-\sqrt{2}} &>2 \pi \\
4 \sqrt{2-\sqrt{2}} &>\pi \\
3.061... &>\pi 
 \end{aligned}
}

本件においては、多角形から円周率が近似できるということだけが分かればいいので、詳細は別な方の説明に全部おまかせ。

ちなみに六角形で考えた場合では、円周の長さが3となり、3.05より大きいことを示せないのです。

■2. 導出

■2.1. 概要

「1.2. 2003年の東大の入試問題」の考えを使って、角数を多くして、円周率に近似していくことを考えていきます。

半角の公式を利用して、四角形 →(2倍)→ 八角 →(2倍)→ 十六角形 →(2倍)→ 三十二角形・・・と角数を倍々に大きくしていきます。

余弦定理を使って八角形から始めてもいいのですが、半角の公式だけでもお腹いっぱいなので、余弦定理なんてものを使わないようにします。そのため、四角形から開始します。)

また、せっかくなので、内接だけでなく、外接についても考えていきます。

■2.2. 変数定義

下図のように、変数を定義します。

半径は1とします。

 <内接の場合>           <外接の場合>    

\displaystyle{
\begin{aligned}
 n&:角数(4の倍数) ※n=4なら四角形、n=8なら八角形な感じ \\
 \theta _{n}&:角度 \\
 s_n&:辺の半分の長さ(内接の場合) \\
 c_n&:\cosを取ったもの(図にはないが、計算過程で使用する) \\
 t_n&:辺の半分の長さ(外接の場合) \\
\end{aligned}
}



■2.3. $s_n, c_n, t_n$の導出

まず、$n=4$の四角形の場合(上図の場合)の$s_n, c_n, t_n$は、

\displaystyle{
\begin{aligned}
 \theta_4&=45°  \\
 s_4 &= \sin \theta _4 = \sin 45° = \dfrac{\sqrt{2}}{2} \\
 c_4 &= \cos \theta _4 = \cos 45° = \dfrac{\sqrt{2}}{2} \\
 t_4 &= \tan \theta _4 = 1
\end{aligned}
}


次に、$2n$を$n$から作成することを考えます。

(すなわち、半分の角度$\theta _{2n} = \frac{\theta_n}{2}$の場合を考えます。)

※これを考えることができれば、$n=4$が分かっているので、$n=4$から$n=8$がわかり、$n=8$から$n=4$がわかり、・・・と、ドミノ倒しのように4の倍数の$n$が全てわかります。

半角の公式より、下記となる。($0<\theta<90°$のため、$\pm$のマイナスは考慮不要。)

\displaystyle{
\begin{aligned}
 s_{2n} &= \sin \theta _{2n} = \sin \dfrac{\theta_n}{2} =\sqrt{\dfrac{1-\cos \theta_n }{2}} = \sqrt{\dfrac{1-c_n }{2}} \\
 c_{2n} &= \cos \theta _{2n} = \cos \dfrac{\theta_n}{2} = \sqrt{\dfrac{1+\cos \theta_n }{2}} = \sqrt{\dfrac{1+c_n }{2}} \\
 t_{2n} &= \tan \theta _{2n} = \dfrac{\sin \dfrac{\theta_{2n} }{2}}{\cos \dfrac{\theta_{2n} }{2}} = \dfrac{s_{2n}}{c_{2n}} \\
\end{aligned}
}


以上より、$s_{2n}, c_{2n}, t_{2n}$が$s_{n}, c_{n}, t_{n}$から計算できる。

(計算は、和差積商とルートのみ。)

■2.4. 円周率の導出

$s_{n}, t_{n}$から多角形の円周の長さを計算します。

<内接の場合>

$s_{n}$に2を掛けると、1辺の長さ。

さらにnを掛けると、周の長さ。

よって、

\displaystyle{
\begin{aligned}
 2\pi &= 2ns_{n} \\
 \pi &= ns_{n}
\end{aligned}
}


<外接の場合>

$t_{n}$に2を掛けると、1辺の長さ。

さらにnを掛けると、周の長さ。

よって、

\displaystyle{
\begin{aligned}
 2\pi &= 2nt_{n} \\
 \pi &= nt_{n}
\end{aligned}
}


となって、円周率を求めることができるのです。

■3. プログラム(C++言語)

下記プログラムとその結果を記載。

結果は、内接の場合、外接の場合の平均としました。

桁数をたくさん表示したいので、Boost Multiprecision Libraryを使用しました。

<プログラム>

#include <iostream>
#include <cmath>

using namespace std;

// Boost Multiprecision Library を include
#include <boost/multiprecision/cpp_dec_float.hpp>

int main(void)
{
    int num;
    int i;
    boost::multiprecision::cpp_dec_float_100 n;
    boost::multiprecision::cpp_dec_float_100 s, c, t;
    boost::multiprecision::cpp_dec_float_100 pi1, pi2, pi3;

    num = 25;
    n = 4;
    s = sqrt(2) / 2;
    c = sqrt(2) / 2;
    t = 1;

    pi1 = n * s;
    pi2 = n * t;
    pi3 = (n * s + n * t) / 2;

    cout << std::setprecision(numeric_limits<decltype(n)>::digits10 + 1) << "角数:" << n << endl;
    cout << std::setprecision(numeric_limits<decltype(pi3)>::digits10 - 30) << "円周率:" << pi3 << endl;

    for (i = 1; i <= num; i++) {
        s = sqrt((1 - c) / 2);
        c = sqrt((1 + c) / 2);
        t = s / c;

        n *= 2;

        pi1 = n * s;
        pi2 = n * t;
        pi3 = (n * s + n * t) / 2;

        cout << endl;
        cout << std::setprecision(numeric_limits<decltype(n)>::digits10 + 1) << "角数:" << n << endl;
        cout << std::setprecision(numeric_limits<decltype(pi3)>::digits10 - 30) << "円周率:" << pi3 << endl;
    }

    return 0;
}

<実行結果>

角数:4
円周率:3.4142135623730951454746218587388284504413604736328125

角数:8
円周率:3.187587978952738995638360618576982603094104740481711631531989631548283

角数:16
円周率:3.152021515166289921864346688035169105826944550560999613192868755739024

角数:32
円周率:3.144136698987597407042045831508800876173960653974967244515617753424677

角数:64
円周率:3.142224771100328313931922566098387118466988501380916749287786963045404

角数:128
円周率:3.141750440437614583250554257543883938743795574657618525333868560016284

角数:256
円周率:3.14163208515663349427522123514159754702370759932784160717791443339265

角数:512
円周率:3.141602510535136321492337553922590596832276571284802741590945720708669

角数:1024
円周率:3.141595117766984050263282232691365008267211045324845330385440364502857

角数:2048
円周率:3.141593269630394238743024839635736552438027984179288978679497379265493

角数:4096
円周率:3.141592807599712253132530870679134771313837411690570090949953093800364

角数:8192
円周率:3.141592692092258347672573904294432563312215322918591853298209298105053

角数:16384
円周率:3.14159266321540840822960009059583505299860181897968836888224516653819

角数:32768
円周率:3.141592655996196769426296644873828321465550314442478187888364688059035

角数:65536
円周率:3.141592654191393912604057878354481940682812738171009093605341892412497

角数:131072
円周率:3.141592653740193201703409834757059946870990834262754510725626436420853

角数:262144
円周率:3.141592653627393024184804801150362674619278502901085041952755476935821

角数:524288
円周率:3.141592653599192979818063353818395625635439222968016111203990559449286

角数:1048576
円周率:3.141592653592142968727184855177087382237419816587431798387218101468274

角数:2097152
円周率:3.141592653590380465954515659466237835292950892038030372810835291495926

角数:4194304
円周率:3.141592653589939840261351512347867750894289653351486057170900240128743

角数:8388608
円周率:3.141592653589829683838060672556359123030065206452174981235876434121603

角数:16777216
円周率:3.141592653589802144732237974920237209380901490263844907331666936793383

角数:33554432
円周率:3.141592653589795259955782301280691433675755044400500319306175072814399

角数:67108864
円周率:3.141592653589793538761668382918897783668662442953377587658784576629569

角数:134217728
円周率:3.141592653589793108463139903331455170786838878839949771251389658904192