圆周率

计算圆周率，最简单的是莱布尼茨公式：

\[\begin{align} \arcsin x &= x-\frac{x^3}{3}+\frac{x^5}{5}-\cdot \cdot \cdot \\ 代入x=1得：\frac{\pi}{4} &=\sum_{i=0}^{\infty}{\frac{(-1)^{i}}{2i+1}} \end{align} \]

但这个公式很慢，一秒内只能计算20位左右，于是我们要更强的公式：Chudnovsky

\[\frac{1}{\pi}=\frac{1}{53360\sqrt{640320}}\sum_{i=0}^{\infty}(-1)^i\frac{(6i)!}{(i!)^3(3i)!}\frac{13591409+545140134i}{640320^{3i}} \]

#include <bits/stdc++.h>
#include <quadmath.h>
//-std=gnu++11 -lquadmath
using namespace std;
__float128 pi;
__float128 xi;
__float128 ge;
__float128 jie(int x){
	__float128 r=1.;
	for(int i=1;i<=x;i++) r=r*i;
	return r;
}
char s[1000];
int n;
int main(){
	__float128 l=800.199,r=800.200,m;
        // 二分
	for(int i=1;i<=1e7;i++){
		m=(l+r)/2.;
		if(m*m>640320.) r=m;
		else l=m; 
	}
	ge=m;
	xi=1./(53360.*ge);
	pi=0;
	cin >> n;
	for(int i=0;i<n;i++){
		__float128 v=1.;
		if(i%2) v=-1.;
		__float128 res;
		res=jie(6*i)/(jie(i)*jie(i)*jie(i)*jie(3*i));
		res=res*(13591409.+545140134.*i);
		for(int j=1;j<=3*i;j++){
			res=res/640320.;
		}
		pi=pi+res*v;
	}
	quadmath_snprintf(s,sizeof(s),"%.100Qf",1./(xi*pi));
	puts(s);
        // 正确的圆周率
	puts("3.1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679");
	return 0;
}