各类求自然数幂和方法

Posted 2021-02-06 zjlcnblogs

高斯消元

我们知道：
[sum_{i=1}^{n}i=frac{n(n+1)}{2}]
以及：
[sum_{i=1}^{n}i^2=frac{n(n+1)(2n+1)}{6}]
以及：
[sum_{i=1}^{n}i^3=(sum_{i=1}^{n}i)^2=(frac{n(n+1)}{2})^2]
那我们可以猜想，自然数的(k)次幂和对应的公式是一个次数为(k+1)的没有常数项的多项式（实际上也是的）。
证明吗，暂时不会。。。

However，我们可以拿这个猜想做题。

设这个(k+1)次的多项式(f(x)=sum_{i=1}^{k+1}a_ix^i)
利用待定系数法，我们只需要知道(k+2)对((x,f(x)))，列出方程组就能解出所有的(a_i)，从而就能代入更大的(x)求出(f(x))。

由于解方程组需要用到高斯消元算法，时间复杂度是(O(k^3))，在(kleq 100)的范围内还是能无压力解决的。

总结

时间复杂度：(O(k^3))
空间复杂度：(O(k^2))

由于高斯消元时要在模意义下做除法，对于模数不是质数的情况无法适应，而且时间复杂度难以接受，不是一种较常用的方法。

第二类斯特林数

分析

定义(S(n,m))表示(n)个有差别的球放入m个无差别的盒子中的方案数，要求盒子不能为空。

容易得到下面的递推式：
[S(n,m)=S(n-1,m-1)+mS(n-1,m)]
考虑新加入的球，要么放在新的盒子里，要么放在之前的盒子里。因为球是有差别的，所以放在任意一个盒子里的方案都是不一样的，因此(S(n-1,m))要乘上一个(m)。

要用它解决自然数幂和问题，还是要用到第二类斯特林数的一个性质：
[a^k=sum_{i=0}^{k}S(k,i)i!C_a^i]
这个性质还是很好解释的，我们可以把(a^k)当做(k)个有差别的球，放入(a)个有差别的盒子的方案数，盒子可以为空。
那么我们就枚举(i)个盒子被放满了，(S(k,i))只保证了球有差别，乘以(i!)相当于给盒子编号，令盒子也有差别，最后乘上一个(C_a^{i})表示在(a)个盒子中选(i)个的方案数。

那么就可以开始化自然数幂求和的式子：
(sum_{a=1}^{n}a^k)
(=sum_{a=1}^{n}sum_{i=0}^{k}S(k,i)i!C_a^i)
两个sigma没有关联，我们可以交换枚举顺序：
(=sum_{i=0}^{k}S(k,i)i!sum_{a=1}^{n}C_a^i)
由于(a<i)时(C_a^{i}=0)，又可以化成：
(=sum_{i=0}^{k}S(k,i)i!sum_{a=i}^{n}C_a^i)

继续化简需要用到一个性质：
(sum_{a=i}^{n}C_a^i=C_{n+1}^{i+1})
证明考虑运用组合数递推公式即：
(C_i^j=C_{i-1}^j+C_{i-1}^{j-1})
(C_{n+1}^{i+1})
(=C_n^i+C_n^{i+1})
(=C_n^i+C_{n-1}^i+C_{n-1}^{i+1})
(=C_n^i+C_{n-1}^i+C_{n-2}^i+C_{n-2}^{i+1})
继续化下去就会得到：
(=sum_{a=i}^{n}C_a^i)

性质就得证了，上面的式子就化简为：
(=sum_{i=0}^{k}S(k,i)i!C_{n+1}^{i+1})
组合数有点麻烦，我们展开为阶乘形式：
(=sum_{i=0}^{k}S(k,i)i!frac{(n+1)!}{(i+1)!(n-i)!})
拆开((i+1)!=i!*(i+1))：
(=sum_{i=0}^{k}S(k,i)frac{(n+1)!}{(i+1)(n-i)!})
发现(frac{(n+1)!}{(n-i)!})其实是(i+1)个连续整数相乘，其中必有一个是(i+1)的倍数，因此式子一定取整数，就不用考虑模数的问题了。

那么直接枚举这(i+1)个连续的整数，得到了时间复杂度为(O(k^2))的算法。计算斯特林数和求自然数幂的复杂度都是(O(k^2))，总复杂度就是(O(k^2))。

Code

附带分解乘法黑科技

#include <cstdio>
#include <cstring>

typedef long long ll;
const int N = 2007;

ll k, n, p, s[N][N];

ll multi(ll a, ll b)
{
    ll x1 = a / 1000000, x2 = a % 1000000, y1 = b / 1000000, y2 = b % 1000000;
    return (x1 * 1000000 % p * y1 % p * 1000000 % p + x1 * 1000000 % p * y2 % p + y1 * 1000000 % p * x2 % p + x2 * y2 % p) % p;
}

ll solve(ll n)
{
    if (n == 0) return 0;
    ll ret = 0;
    for (int i = 1; i <= k && i <= n; i++)
    {
        ll sum = s[k][i];
        for (ll j = n - i + 1; j <= n + 1; j++)
            if (j % (i + 1) == 0) sum = multi(sum, j / (i + 1));
            else sum = multi(sum, j);
        ret = (ret + sum) % p;
    }
    return ret;
}

int main()
{
    scanf("%lld%lld%lld%lld", &k, &n, &p);
    s[0][0] = 1;
    for (int i = 1; i <= k; i++)
        for (int j = 1; j <= i; j++)
            s[i][j] = (s[i - 1][j - 1] + multi(j, s[i - 1][j])) % p;
    printf("%lld
", solve(n) % p);
    return 0;
}

总结

时间复杂度：(O(k^2))
空间复杂度：(O(k^2))

这种做法由于不用除法而适用于模数为任意数的情况，但是求斯特林数复杂度是(O(k^2))的，当(k)较大时不再适用。

拉格朗日插值法

大坑待填。。。。