HDU6217 BBP Formula

简略题意：问π的十六进制表示的小数点后第n位是多少。

关于BBP公式详见WIKI百科

这个公式的作用就是计算十六进制下的小数点后第n位，而不必计算前n-1项。
$\pi = \sum_{k=0}^{\infty}[\frac{1}{16^k}(\frac{4}{8k+1})-(\frac{2}{8k+4})-(\frac{1}{8k+5})-(\frac{1}{8k+6})]$。
抽取第一个求和部分。
$\sum_{k=0}^{\infty}\frac{1}{16^k}(\frac{4}{8k+1})$
= $\sum_{k=0}^{\infty}(\frac{4}{16^k(8k+1)})$
= 4$\sum_{k=0}^{\infty}(\frac{1}{16^k(8k+1)})$

$\sum_{k=0}^{\infty}(\frac{1}{16^k(8k+1)})$，我们想要获取第n位，把原式子拆成两部分。
$\sum_{k=0}^{\infty}(\frac{1}{16^k(8k+1)})$
=$\sum_{k=0}^{n}(\frac{1}{16^k(8k+1)}) + \sum_{k=n + 1}^{\infty}(\frac{1}{16^k(8k+1)})$
将等式乘上$16^n$,使得小数点恰好落在第n位。
$\sum_{k=0}^{n}(\frac{1}{16^k(8k+1)}) + \sum_{k=n + 1}^{\infty}(\frac{1}{16^k(8k+1)})$ = $\sum_{k=0}^{n}(\frac{16^{n-k}}{(8k+1)}) + \sum_{k=n + 1}^{\infty}(\frac{16^{n-k}}{(8k+1)})$
我们只关心上式的小数部分，为了避免高精度计算，$\sum_{k=0}^{n}(\frac{16^{n-k}}{(8k+1)})$ = $\sum_{k=0}^{n}(\frac{16^{n-k} mod (8k+1)}{(8k+1)})$,两者的小数部分是相同的，并且对于后半部分，$\sum_{k=n + 1}^{\infty}(\frac{16^{n-k}}{(8k+1)})$，这部分对精度要求并不高，我们向后处理一个足够多的位数就够了。
令$\sum1 = \sum_{k=0}^{n}(\frac{16^{n-k}}{(8k+1)}) + \sum_{k=n + 1}^{\infty}(\frac{16^{n-k}}{(8k+1)})$， 那么答案就是$4\sum1 - 2\sum2 - \sum3 - \sum4$的小数部分，乘以16之后，得到的整数部分转化成16进制即可。
复杂度为$O(nlogn)$。

#define others
#ifdef poj
#include <iostream>
#include <cstring>
#include <cmath>
#include <cstdio>
#include <algorithm>
#include <vector>
#endif // poj
#ifdef others
#include <bits/stdc++.h>
#endif // others
//#define file
#define all(x) x.begin(), x.end()
using namespace std;
const double eps = 1e-8;
int dcmp(double x) { if(fabs(x)<=eps) return 0; return (x>0)?1:-1;};
typedef long long LL;
typedef unsigned long long ULL;

namespace solver {
    char out(int x) {
        if(0 <= x && x <= 9) return x + '0';
        else if(x == 10) return 'A';
        else if(x == 11) return 'B';
        else if(x == 12) return 'C';
        else if(x == 13) return 'D';
        else if(x == 14) return 'E';
        else if(x == 15) return 'F';
    }
    LL Pow(LL a, LL b, LL mod) {
        LL res = 1;
        while(b) {
            if(b & 1) res *= a, res %= mod;
            b >>= 1;
            a *= a, a %= mod;
        }
        return res;
    }
    double BBP(int n, LL k, LL b) {
        double val = 0;
        for(int i = 0; i <= n; i++)  val += (Pow(16, n-i, 8*i+b)*1.0/(8*i+b));
        for(int i = n + 1; i <= n + 1 + 1000; i++) val += (powf(16, n - i)/(8*i+b));
        return k*val;
    }
    void solve() {
        int t;
        scanf("%d", &t);
        int f = 0;
        while(t--) {
            double v = 0;
            int n;
            scanf("%d", &n);
            n--;
            v = BBP(n, 4, 1) + BBP(n, -2, 4) + BBP(n, -1, 5) + BBP(n, -1, 6);
            v = v - (int)v;
            if(v < 0) v += 1;
            v *= 16;
            printf("Case #%d: %d %c\n", ++f, n + 1, out((int)v));
        }
    }
}

int main() {
    solver::solve();
    return 0;
}