UVA 10236 The Fibonacci Primes

UVA_10236

    这个题目AC得太辛苦了,也再一次加深了对浮点数误差的理解,库函数能不调还是最好别调……其实后来感觉就像UVA论坛上一个人说的,AC这个题有两个准则,第一个是要用long double,第二个是不要相信fmod、sqrt、pow等这些库函数,它们都会带来精度上的问题。

    下面回到题目,这个题目涉及的知识还是满多的,又让我学到了不少。首先,我们要知道fibonacci素数的特征,fibonacci数列的第3、4项的2、3为fibonacci素数,第5项及其之后的素数项均为fibonacci素数,因此我们就要知道第22000个素数大概在什么范围。百度了之后发现10^6之内可以产生7万多个素数,因此我们可以先把10^6内的素数筛出来存到数组里,以便后面使用。

    下面我们就要知道从什么时候开始fibonacci数会超过9位,我测了一下第44个fibonacci数是最后一个9位及9位以内的fibonacci数,对应的就是第14个fibonacci素数,因此我们先可以预处理出15以内的结果,剩下的再用函数去算。

    在算的时候,我们首先要找到第n个fibonacci素数的项数x,然后再求F(x)的前9位即可,在求的时候为了确保有9位精度,我们选择用long double形式的快速幂去计算,并且每次计算完之后,要把结果卡在某个范围以下,这样是为了避免后续的运算会超过long double的表示范围。

#include<stdio.h>
#include<string.h>
#include<math.h>
long long int fib[60];
int N, prime[1000010], p[100010];
long double a[60][4], b[4];
void prepare()
{
int i, j, k, n;
fib[0] = 0, fib[1] = 1;
for(i = 2; i < 50; i ++)
fib[i] = fib[i - 1] + fib[i - 2];
k = (int)sqrt((double)1000001);
memset(prime, -1, sizeof(prime));
n = 0;
for(i = 2; i <= k; i ++)
if(prime[i])
{
p[n ++] = i;
for(j = i * i; j < 1000001; j += i)
prime[j] = 0;
}
for(; i < 1000001; i ++)
if(prime[i])
p[n ++] = i;
a[0][0] = a[0][1] = a[0][2] = 1;
a[0][3] = 0;
}
void pow_mod(int n, int e)
{
if(n == 1)
{
a[e][0] = a[e][1] = a[e][2] = 1;
a[e][3] = 0;
return ;
}
pow_mod(n / 2, e + 1);
a[e][0] = a[e + 1][0] * a[e + 1][0] + a[e + 1][1] * a[e + 1][2];
a[e][1] = a[e + 1][0] * a[e + 1][1] + a[e + 1][1] * a[e + 1][3];
a[e][2] = a[e + 1][2] * a[e + 1][0] + a[e + 1][3] * a[e + 1][2];
a[e][3] = a[e + 1][2] * a[e + 1][1] + a[e + 1][3] * a[e + 1][3];
if(n % 2)
{
b[0] = a[e][0] * a[0][0] + a[e][1] * a[0][2];
b[1] = a[e][0] * a[0][1] + a[e][1] * a[0][3];
b[2] = a[e][2] * a[0][0] + a[e][3] * a[0][2];
b[3] = a[e][2] * a[0][1] + a[e][3] * a[0][3];
a[e][0] = b[0], a[e][1] = b[1], a[e][2] = b[2], a[e][3] = b[3];
}
while(a[e][0] > 10000)
a[e][0] /= 10, a[e][1] /= 10, a[e][2] /= 10, a[e][3] /= 10;
}
void solve()
{
int i, j, n;
long long int t;
long double x;
n = p[N - 1];
pow_mod(n - 1, 1);
x = a[1][0];
while(x < 100000000)
x *= 10;
t = (long long int)x;
printf("%lld\n", t);
}
int main()
{
prepare();
while(scanf("%d", &N) == 1)
{
if(N <= 2)
printf("%s\n", N == 1 ? "2" : "3");
else if(N <= 14)
printf("%lld\n", fib[p[N - 1]]);
else
solve();
}
return 0;
}


原文地址:https://www.cnblogs.com/staginner/p/2290388.html