Codeforces D. The Sum of the k-th Powers（拉格朗日插值）

题目描述：

The Sum of the k-th Powers

time limit per test

2 seconds

memory limit per test

256 megabytes

input

standard input

output

standard output

There are well-known formulas: img, img, img. Also mathematicians found similar formulas for higher degrees.

Find the value of the sum img modulo 109 + 7 (so you should find the remainder after dividing the answer by the value 109 + 7).

Input

The only line contains two integers n, k (1 ≤ n ≤ 109, 0 ≤ k ≤ 106).

Output

Print the only integer a — the remainder after dividing the value of the sum by the value 109 + 7.

Examples

Input

Copy

4 1

Output

Copy

10

Input

Copy

4 2

Output

Copy

30

Input

Copy

4 3

Output

Copy

100

Input

Copy

4 0

Output

Copy

4

思路：

我从来没见过这样的题，今天是开了眼界了。听说是数学分析上的拉格朗日插值多项式（当时学的时候迷迷糊糊的）现在竟然就直接跳过理论来实际应用了。还是写代码，虽然知道是数值分析上的，本来就和计算机挂钩，但还是第一次直观感受到了数学在计算机的应用。以前学的时候老是疑惑我学了这些能干嘛？都说是应用广，但对我而言，应用只是一句空洞的话，现在看到了现实，还是有点小兴奋嗷~

在网上找了找资料，先摆一摆公式，以后弄明白了可能会回来补（可能吧）

拉格朗日插值公式

在数值分析中，拉格朗日插值法是以法国十八世纪数学家约瑟夫·拉格朗日命名的一种多项式插值方法。许多实际问题中都用函数来表示某种内在联系或规律，而不少函数都只能通过实验和观测来了解。如对实践中的某个物理量进行观测，在若干个不同的地方得到相应的观测值，拉格朗日插值法可以找到一个多项式，其恰好在各个观测的点取到观测到的值。这样的多项式称为拉格朗日（插值）多项式。数学上来说，拉格朗日插值法可以给出一个恰好穿过二维平面上若干个已知点的多项式函数。——摘自百度百科

$$f(x)=\sum_{i=0}^n y_i\prod_{j!=i}\frac{x-x_j}{x_i-x_j}$$.

简单的解释就是带入一个\(x_i\)，对第\(i\)项来说后面的连乘是1，而对于非\(i\)项来说，后面的连乘是零。所以有\(f(x_i)=y_i\).

一种特殊的情况：\(x_i\)是连续的。

令\(P_i=\prod_{j=0}^i(x-x_j)\)，\(S_i=\prod_{j=i}^n(x-x_j)\).

而由于\(x_i\)是连续的正整数，连称号中的分母就是第一项\((1-2)(1-3)(1-4)...(1-n)\)和第二项\((2-1)(2-3)(2-4)...(2-n)等等\)，可表示为\((-1^{n-i})(n-i)!(i-1)!\)

这是后面就变成了$$\prod_{i!=j} \frac{P_{i-1}S_{i+1}}{(-1){(n-i)}(n-i)!(i-1)!}=\frac{\prod_{j=0}{(x-j)!}}{(-1)^{n-i}(x-i)(n-i)!(i-1)!}$$.

计算时用费马小定理求出\((x-i),(n-i)!,(i-1)!\)的逆元累加即可。

求逆元的复杂度是\(O(log_2n)\)，累加是\(O(n)\)，总的复杂度是\(O(nlong_2n)\).

求幂和

幂和可表示为$$S_k(n)=\dfrac{(n+1){k+1}-\sum\limits_{j=0}\binom{k+1}{j}S_j(n)-1}{k+1}$$(具体推导见参考文章1)

可见为K+1次多项式。我们知道k+1次需要有k+2个点来确定，因此我们可以暴力计算出前k+2项的幂和，如果n小于等于k+2，我们通过这样的预处理可以直接得到答案。如果n>k+2，则使用上面提到的拉格朗日插值求得n的k次幂和。

注意最后得答案时加一些模在模一下保证是个正数。

代码：

#include <iostream>
#define max_n 1000006
const long long mod = 1000000007;
using namespace std;
int k;
int n;
long long fac[max_n];
long long y[max_n];
long long pi = 1;
long long q_mod(long long a,long long b,long long mod)
{
    long long res = 1;
    while(b)
    {
        if(b&1)
        {
            res = (res*a)%mod;
        }
        a = (a*a)%mod;
        b>>=1;
    }
    return res;
}

int main()
{
    cin >> n >> k;
    y[0] = 0;
    for(int i = 1;i<=k+2;i++)
    {
        y[i] = (y[i-1]+q_mod(i,k,mod))%mod;
    }
    if(n<=k+2)
    {
        cout << y[n] << endl;
        return 0;
    }
    fac[0] = 1;
    for(int i = 1;i<=max_n;i++)
    {
        fac[i] = (fac[i-1]*i)%mod;
    }
    for(int i = 1;i<=k+2;i++)
    {
        pi = (pi*(n-i))%mod;
    }
    long long ans = 0;
    for(int i = 1;i<=k+2;i++)
    {
        long long f = ((k+2-i)%2)?-1:1;
        long long inv = q_mod(fac[k+2-i]%mod*fac[i-1]%mod,mod-2,mod);
        long long inv2 = q_mod(n-i,mod-2,mod);
        ans = (ans+f*y[i]*pi%mod*inv%mod*inv2%mod)%mod;

    }
    cout << (ans+mod)%mod << endl;
    return 0;
}

参考文章：

Wolfycz，自然数幂和，https://www.cnblogs.com/Wolfycz/p/10622577.html，浅谈算法——拉格朗日插值，https://www.cnblogs.com/Wolfycz/p/10624678.html

qscqesze，Educational Codeforces Round 7 F. The Sum of the k-th Powers 拉格朗日插值法，https://www.cnblogs.com/qscqesze/p/5207132.html

sdfzchy，【CF622F】The Sum of the k-th Powers （拉格朗日插值法），https://blog.csdn.net/sdfzchy/article/details/78307534

另外还有一篇比较不错的差分求正整数的k次幂和的文章：

张永强.差分的应用及正整数的k 次方幂求和[J].唐山师范学院学报，2007,29(5) :140-141
posted @ 2019-08-13 21:23  小张人  阅读(366)  评论(0)  编辑  收藏  举报
分享到：