【CF622F】The Sum of the k-th Powers (拉格朗日插值法)

用的dls的板子，因为看不懂调了好久...果然用别人的板子就是这么蛋疼- -||

num数组0~k+1储存了k+2个值，且这k+2个值是自然数i的k次方而不是次方和，dls的板子自己帮你算和的...搞得我弄了好久

  1 #include <iostream>
  2 #include <string.h>
  3 #include <cstdio>
  4 #include <vector>
  5 #include <queue>
  6 #include <assert.h>
  7 #include <math.h>
  8 #include <string>
  9 #include <algorithm>
 10 #include <functional>
 11 
 12 #define SIGMA_SIZE 26
 13 #define lson rt<<1
 14 #define rson rt<<1|1
 15 #define lowbit(x) (x&-x)
 16 #define foe(i, a, b) for(int i=a; i<=b; i++)
 17 #define fo(i, a, b) for(int i=a; i<b; i++)
 18 #define pii pair<int,int>
 19 #pragma warning ( disable : 4996 )
 20 
 21 using namespace std;
 22 typedef long long LL;
 23 inline LL LMax(LL a, LL b) { return a>b ? a : b; }
 24 inline LL LMin(LL a, LL b) { return a>b ? b : a; }
 25 inline LL lgcd(LL a, LL b) { return b == 0 ? a : lgcd(b, a%b); }
 26 inline LL llcm(LL a, LL b) { return a / lgcd(a, b)*b; }  //a*b = gcd*lcm
 27 inline int Max(int a, int b) { return a>b ? a : b; }
 28 inline int Min(int a, int b) { return a>b ? b : a; }
 29 inline int gcd(int a, int b) { return b == 0 ? a : gcd(b, a%b); }
 30 inline int lcm(int a, int b) { return a / gcd(a, b)*b; }  //a*b = gcd*lcm
 31 const LL INF = 0x3f3f3f3f3f3f3f3f;
 32 const LL mod = 1e9+7;
 33 const double eps = 1e-8;
 34 const int inf = 0x3f3f3f3f;
 35 const int maxk = 3e6 + 5;
 36 const int maxn = 1e6+10;
 37 
 38 /// 注意mod，使用前须调用一次 polysum::init(int M);
 39 namespace polysum {
 40 #define rep(i,a,n) for (int i=a;i<n;i++)
 41 #define per(i,a,n) for (int i=n-1;i>=a;i--)
 42     typedef long long ll;
 43     const ll mod = 1e9 + 7; /// 取模值
 44     ll powmod(ll a, ll b) { ll res = 1; a %= mod; assert(b >= 0); for (; b; b >>= 1) { if (b & 1)res = res*a%mod; a = a*a%mod; }return res; }
 45 
 46     const int D = 1010000; /// 最高次限制
 47     ll a[D], f[D], g[D], p[D], p1[D], p2[D], b[D], h[D][2], C[D];
 48     ll calcn(int d, ll *a, ll n) {
 49         if (n <= d) return a[n];
 50         p1[0] = p2[0] = 1;
 51         rep(i, 0, d + 1) {
 52             ll t = (n - i + mod) % mod;
 53             p1[i + 1] = p1[i] * t%mod;
 54         }
 55         rep(i, 0, d + 1) {
 56             ll t = (n - d + i + mod) % mod;
 57             p2[i + 1] = p2[i] * t%mod;
 58         }
 59         ll ans = 0;
 60         rep(i, 0, d + 1) {
 61             ll t = g[i] * g[d - i] % mod*p1[i] % mod*p2[d - i] % mod*a[i] % mod;
 62             if ((d - i) & 1) ans = (ans - t + mod) % mod;
 63             else ans = (ans + t) % mod;
 64         }
 65         return ans;
 66     }
 67     void init(int M) { /// M：最高次
 68         f[0] = f[1] = g[0] = g[1] = 1;
 69         rep(i, 2, M + 5) f[i] = f[i - 1] * i%mod;
 70         g[M + 4] = powmod(f[M + 4], mod - 2);
 71         per(i, 1, M + 4) g[i] = g[i + 1] * (i + 1) % mod;
 72     }
 73     ll polysum(ll n, ll *arr, ll m) { /// a[0].. a[m] \sum_{i=0}^{n-1} a[i]
 74         for (int i = 0; i <= m; i++) 
 75             a[i] = arr[i];
 76         a[m + 1] = calcn(m, a, m + 1);
 77         rep(i, 1, m + 2) a[i] = (a[i - 1] + a[i]) % mod;
 78         return calcn(m + 1, a, n - 1);
 79     }
 80     ll qpolysum(ll R, ll n, ll *a, ll m) { /// a[0].. a[m] \sum_{i=0}^{n-1} a[i]*R^i
 81         if (R == 1) return polysum(n, a, m);
 82         a[m + 1] = calcn(m, a, m + 1);
 83         ll r = powmod(R, mod - 2), p3 = 0, p4 = 0, c, ans;
 84         h[0][0] = 0; h[0][1] = 1;
 85         rep(i, 1, m + 2) {
 86             h[i][0] = (h[i - 1][0] + a[i - 1])*r%mod;
 87             h[i][1] = h[i - 1][1] * r%mod;
 88         }
 89         rep(i, 0, m + 2) {
 90             ll t = g[i] * g[m + 1 - i] % mod;
 91             if (i & 1) p3 = ((p3 - h[i][0] * t) % mod + mod) % mod, p4 = ((p4 - h[i][1] * t) % mod + mod) % mod;
 92             else p3 = (p3 + h[i][0] * t) % mod, p4 = (p4 + h[i][1] * t) % mod;
 93         }
 94         c = powmod(p4, mod - 2)*(mod - p3) % mod;
 95         rep(i, 0, m + 2) h[i][0] = (h[i][0] + h[i][1] * c) % mod;
 96         rep(i, 0, m + 2) C[i] = h[i][0];
 97         ans = (calcn(m, C, n)*powmod(R, n) - c) % mod;
 98         if (ans<0) ans += mod;
 99         return ans;
100     }
101 }
102 
103 LL num[maxn];
104 LL n, k;
105 
106 void init()
107 {
108     for( int i = 0; i <= k+1; i++ )
109         num[i] = polysum::powmod((LL)i+1, k); 
110 }
111 
112 int main(){
113     cin >> n >> k;
114     
115     polysum::init(k);
116     init();
117 
118     LL ans = polysum::polysum(n, num, k+1)%mod;
119     printf("%lld\n", ans);
120     return 0;
121 }
View Code

 

posted @ 2018-08-12 20:59  LBNOQYX  阅读(730)  评论(0)  编辑  收藏  举报