HDU 多校训练赛第一场1005 Fibonacci Sum

标签：斐波那契数   rac   namespace   c代码   init   数字   --   表示   names   

题意：

? 给定三个数字 \(n,c,k\) ，求以下式子

? \(\sum_{i=0}^nF(ic)^k\%(10^9+9)\)

? 其中\(F(x)\)为斐波那契数列第\(x\)项。

? \(1\leq n,c\leq10^{18},1\leq k \leq10^5\)

分析：

? 在比赛的时候我搜索了一波斐波那契数列的性质，意外发现了类似题目，ZOJ3774，只不过那个题里\(c=1\)，这里解释一下过程。

? 斐波那契有一个公式：

? \(F(n)=\frac 1 {\sqrt5}((\frac{1+\sqrt5}{2})^n-(\frac{1-\sqrt5}{2})^n)\)

? 我们先看\(c=1\)的情况，也就是每一项都取。即：

? \(sum=\sum_{i=0}^nF(i)^k\%(10^9+9)\)

? \(F(n)^k=(\frac 1 {\sqrt5}((\frac{1+\sqrt5}{2})^n-(\frac{1-\sqrt5}{2})^n))^k\)

? \(F(n)^k=\frac 1 {(\sqrt5)^k}((\frac{1+\sqrt5}{2})^n-(\frac{1-\sqrt5}{2})^n)^k\)

? 我们假设\(a=(\frac{1+\sqrt5}{2}),b=(\frac{1-\sqrt5}{2})\)

? 那么有：

? \(F(n)^k=\frac 1 {(\sqrt5)^k}(a^n-b^n)^k\)

? 而由二项展开式可以得到：

\((a-b)^k=C_k^0(a^n)^0(-b^n)^k+C^{k-1}_1(a^n)^1(-b^n)^{k-1}+...+C_k^k(a^n)^k(-b^n)^0\)

? 那么神奇的来了：

? \(F(n)^k+F(n-1)^k=\frac 1 {(\sqrt5)^k}((a^n-b^n)^k+(a^{n-1}-b^{n-1})^k)\)

? \(((a^n-b^n)^k+(a^{n-1}-b^{n-1})^k)\)

\(=C_k^0(a^n)^0(-b^n)^k+C_k^0(a^{n-1})^0(-b^{n-1})^k+...\)

\(=C_k^0(a^0)^n((-b)^k)^n+C_k^0(a^0)^{n-1}((-b)^k)^{n-1}+...\)

\(=(-1)^kC_k^0[(a^0b^k)^n+(a^0b^k)^{n-1}]+...\)

? 如果继续加上\(F[n-2]+F[n-3]+...\)

? 会发现构成一个等比数列求和的问题。

\((-1)^kC_k^0[(a^0b^k)^n+(a^0b^k)^{n-1}+(a^0b^k)^{n-2}+...++(a^0b^k)^{1}]\)

? 可以把 \(a^0b^k\) 看成公比，那么由求和公式：

? \(S_n=a_1\times \frac{q^n-1}{q-1}\)

? 令\(a^0b^k=x\) 得：

? \(\sum_{i=1}^nF(i)^k=\frac{1}{\sqrt5^k}\sum_{i=0}^k(-1)^kx\frac{x^n-1}{x-1}(x=a^ib^{k-i})\)

? 大功告成，如果看懂了这个式子，也就成功了一半了,现在我们加上\(c\) 的限制。

? \(F(cn)^k=\frac 1 {(\sqrt5)^k}(a^{cn}-b^{cn})^k\)

? 接下来的步骤是一样的，会发现结果是：

? \(\sum_{i=1}^nF(ci)^k=\frac{1}{\sqrt5^k}\sum_{i=0}^k(-1)^kx\frac{x^n-1}{x-1}(x=(a^ib^{k-i})^c)\)

? 现在需要处理一下根号的问题，因为是分数取模。

? 首先得到\(i^2mod(1e9+9)==5\)时\(i=383008016\)

? 然后根据这个值得出上述\(a,b\) 的值：

#include <bits/stdc++.h>
using namespace std;

long long mod=1e9+9;

long long fastpow(long long a, long long b)
{
    long long res = 1;
    a%=mod;
    while(b){
        if(b & 1)
            res = res * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return res;
}

int main(){
    printf("%lld\n",(1+383008016)*fastpow(2,mod-2)%mod);
    printf("%lld\n",((1-383008016)%mod+mod)%mod*fastpow(2,mod-2)%mod);
}

? 得到\(a=691504013,b=308495997\)

? 接下来就是代码的事情了，先放出来比赛时疯狂TLE的代码来帮助理解：

#include <bits/stdc++.h>
using namespace std;

const int MAXN = 1e5 + 10;
const long long mod = 1e9 + 9;
long long fac[MAXN], A[MAXN], B[MAXN];
//fac数组表示阶乘，之后求C用
//A[i]表示a的i次方，B[i]同理

void Init() {
    fac[0] = 1;
    for (int i = 1; i < MAXN; i++)
        fac[i] = fac[i - 1] * i % mod;
    A[0] = B[0] = 1;
    for (int i = 1; i < MAXN; i++) {
        A[i] = A[i - 1] * 691504013 % mod;
        //691504013表示的意义上文提到了
        B[i] = B[i - 1] * 308495997 % mod;
        //同上
    }
}

long long fastpow(long long a, long long b) {//快速幂
    long long ans = 1;
    a %= mod;
    while (b) {
        if (b & 1)ans = ans * a % mod;
        b >>= 1;
        a = a * a % mod;
    }
    return ans;
}

long long Solve(long long n, long long c, long long k) {
    long long ans = 0;
    for (int i = 0; i <= k; i++) {
        long long x = fastpow(A[k - i] * B[i] % mod, c);
        //首先求出公比
        long long C = fac[k] * fastpow(fac[k - i] * fac[i] % mod, mod - 2) % mod;
        //求出组合数
        long long tmp = x * (fastpow(x, n) - 1) % mod * fastpow(x - 1, mod - 2) % mod;
        //求出等比数列之和
        if (x == 1) tmp = n % mod;
        //如果公比为1，则是一个恒等数列，特判
        tmp = tmp * C % mod;
        //求出最后的值
        if (i & 1) ans -= tmp;
        //判断是加是减
        else ans += tmp;
        ans %= mod;
    }
    long long num = fastpow(383008016, mod - 2);
    //算出根号5在分母的值
    ans = ans * fastpow(num, k) % mod;
    //算出k次方
    ans = (ans % mod + mod) % mod;
    //使结果为正
    return ans;
}

int main() {
    int t;
    long long n, k, c;
    Init();
    scanf("%d", &t);
    while (t--) {
        scanf("%lld%lld%lld", &n, &c, &k);
        printf("%lld\n", Solve(n, c, k));
    }
    return 0;
}

? 然后这个代码TLE了十几次。

优化

? 我感觉是卡了常数，优化的话首先有一个等式：

? \(a^b\%x=a^{b\%(x-1)}\%x\)

? 第二点，上面的代码每次都会对\(x\)取快速幂一次，这个也是很关键的一点，因为每次循坏都取。观察发现，其实只需要第一次求出

? \((a^c)^0(b^c)^k\)的值，之后把这个值乘上\(\frac{a^c}{b^c}\)即可得到下一项\((a^c)^1(b^c)^{k-1}\)

? 那我们可以省去一个常数。

? AC代码：

#include <bits/stdc++.h>
using namespace std;

const int MAXN = 1e5 + 10;
const long long mod = 1e9 + 9;
long long fac[MAXN], inv[MAXN];

long long fastpow(long long a, long long b) {
    long long ans = 1;
    a %= mod;
    while (b) {
        if (b & 1)ans = ans * a % mod;
        b >>= 1;
        a = a * a % mod;
    }
    return ans;
}

void Init() {
    fac[0] = inv[0] = 1;
    for (int i = 1; i < MAXN; i++) {
        fac[i] = fac[i - 1] * i % mod;
        inv[i] = fastpow(fac[i], mod - 2);
    }
}

long long Solve(long long n, long long c, long long k) {
    long long ans = 0;
    long long A = fastpow(691504013, c % (mod - 1)), B = fastpow(308495997, c % (mod - 1));
    long long a = 1, b = fastpow(B, k), ib = fastpow(B, mod - 2);
    for (int i = 0; i <= k; i++) {
        long long x = a * b % mod;
        long long C = fac[k] * inv[i] % mod * inv[k - i] % mod;
        long long sum = x * (fastpow(x, n % (mod - 1)) - 1 + mod) % mod * fastpow(x - 1, mod - 2) % mod;
        if (x == 1) sum = n % mod;
        if ((k - i) & 1) ans -= sum * C % mod;
        else ans += sum * C % mod;
        ans %= mod;
        a = a * A % mod;
        b = b * ib % mod;
    }
    long long num = fastpow(383008016ll, mod - 2);
    ans = ans * fastpow(num, k) % mod;
    ans = (ans % mod + mod) % mod;
    return ans;
}

int main() {
    int t;
    long long n, k, c;
    Init();
    scanf("%d", &t);
    while (t--) {
        scanf("%lld%lld%lld", &n, &c, &k);
        printf("%lld\n", Solve(n, c, k));
    }
    return 0;
}

HDU 多校训练赛第一场1005 Fibonacci Sum

标签：斐波那契数   rac   namespace   c代码   init   数字   --   表示   names   

原文地址：https://www.cnblogs.com/jianjinxiani/p/13358249.html