求解自然数幂和的若干种方法

问题的引入

给定\(n,k\)求\[\sum_{i=1}^ni^k\]

1. 循环

四年级应该会循环了。

能做到\(O(nk)\)的优秀时间复杂度。

2. 快速幂

五年级学了快速幂之后就能做到\(O(nlog_2k)\)

请不要小看这个算法。有时候在特定的情况下(例如\(n\)很小,或\(1\rightarrow n\)的距离变得很小时),这个复杂度真的很优秀。

3. 差分法

六年级应该知道差分和二项式定理了。那么:\[(a+1)^k-a^k=\sum_{i=0}^{k-1}C_k^ia^i\]

于是:
\[ (n+1)^k-1 =\sum_{i=1}^n (i+1)^k-i^k \\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ =\sum_{i=1}^n\sum_{j=0}^{k-1}C_k^ji^j\\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ =\sum_{i=0}^{k-1}C_k^i\sum_{j=1}^n j^i\\ \ \ \ \ \ \ \ \ \ \ \ \ \ =\sum_{i=0}^{k-1}C_k^i S(i)
\]\[\therefore (n+1)^{k+1}-1=\sum_{i=0}^kC_{k+1}^iS(i)\]

把\(i=k\)时移项,可以得到\[(k+1)S(k)=(n+1)^{k+1}-1-\sum_{i=0}^{k-1}C_{k+1}^iS(i)\]

所以\[S(k)=\frac{(n+1)^{k+1}-1-\sum_{i=0}^{k-1}C_{k+1}^iS(i)}{k+1}\]

同时仔细观察这个式子,我们发现,\(k\)次方和的求和公式是\(k+1\)次的。归纳证明即可。

3. 倍增

初一应该会倍增了,所以我们令\(f_{n,k}=\sum_{i=1}^ni^k\)。

当\(n\)是奇数的时候直接由\(f_{n-1,k}+n^k\)转移过来。偶数的时候拆开来,运用简单的二项式定理,一波式子推得:\[f(n,k)=f(\frac{n}{2},k)+\sum_{j=0}^kC_k^j*f(\frac{n}{2},j)*\frac{n}{2}^{k-j}\]

每一层的\(f_{n,k}\)我们计算的时间复杂度都是\(O(k^2)\)的,\(log_n\)层,时间复杂度\(O(k^2log_n)\).

4. 高斯消元

初一应该会高斯消元了。这是个大脑洞。虽然时间复杂度比上一个还劣一些。

根据\(k\)次方和的求和公式是\(k+1\)次的,所以列出\(k+2\)条式子就可以唯一确定这个多项式。

时间复杂度\(O(k^3)\).

5. 第一类斯特林数

初二来学习一下斯特林数。

第一类斯特林数我们一般清楚的是它的组合意义,即把\(n\)个元素分成\(k\)个圆排列的方案。根据组合意义,我们不难推出它的式子是\[S_u(n,m)=S_u(n-1,m-1)+(n-1)S_u(n-1,m)\]

但事实上,我们求解自然数幂和需要用到的是它的原始定义
\[x^{n\downarrow}=x\cdot (x-1)\cdot (x-2)\cdots (x-n+1)=\sum_{k=0}^ns_s(n,k)\cdot x^k\]\[x^{n\uparrow}=x\cdot (x+1)\cdot (x+2)\cdots(x+n-1)=\sum_{k=0}^ns_u(n,k)\cdot x^k\]

这里需要注意,第一类斯特林数根据定义分成了有符号\(S_s\)和无符号\(S_u\)两种。事实上,我们可以很轻松的从这个原始定义推出它的组合意义

因为\[\sum_{k=0}^nS_u(n,k)·x^k=x^{n\uparrow}=x^{(n-1)\uparrow}·(x+n-1)\]\[=\sum_{k=0}^{n-1}S_u(n-1,k)x^{k+1}+(n-1)\sum_{k=0}^{n-1}S_u(n-1,k)x^k\]

对比两边\(x^m\)的系数,可以得到\[S_u(n,m)=S_u(n-1,m-1)+(n-1)S_u(n-1,m)\]继续推有符号的,可以得到\[S_s(n,m)=S_s(n-1,m-1)-(n-1)S_s(n-1,m)\]事实上,我们可以完全不用记第一类斯特林数的组合意义,通过公式直接推出来即可。当然记了更好,还可以验证。

所以,根据第一类斯特林数的的定义,得到:\[\prod_{x=0}^{k-1}(n-x)=\sum_{k=0}^nS_s(n,k)x^k\]

于是我们可以得到一个显然的式子是:\[n^m=n^{m\downarrow}-\sum_{k=0}^{m-1}S_s(m,k)·n^k\]

我们继续推,发现下降幂的和是可以写成一个组合数的形式的,比方说\[\sum_{i=m}^ni^{m\downarrow}=\sum_{i=m}^n\frac{i!m!}{(i-m)!m!}=m!\sum_{i=m}^n\binom{i}{m}=m!\binom{n+1}{m+1}\]

而后面那一坨式子也是可以化简的,比方说\[\sum_{i=0}^n\sum_{k=0}^{m-1}S_s(m,k)·i^k=\sum_{k=0}^{m-1}S_s(m,k)\sum_{i=0}^ni^k\]

发现后面那条式子\(\sum_{i=0}^ni^k\)的\(k\)是降了阶的,所以可以边处理边记录一下,就不用重新算了,时间复杂度就变成了\(O(k^2)\)。而处理\(S_s(m,k)\)也是\(O(k^2)\)级别。

事实上如果当你升入初三,\(S_s(m,k)\)就可以运用分治\(NTT\)做到\(O(klog^2k)\)了,虽然然并卵。

但请注意,这种方法虽然时间复杂度是\(O(k^2)\)级别的,但是它并非没有什么用,因为它——不用做除法

6. 第二类斯特林数

第二类斯特林数的组合意义就是\(n\)个元素分成\(m\)个集合,且集合非空的方案数。

基本性质是\[\begin{Bmatrix}n\\m\end{Bmatrix}=\begin{Bmatrix}n-1\\m-1\end{Bmatrix}+m\cdot \begin{Bmatrix}n-1\\m\end{Bmatrix}\]

考虑它的通项公式,可以先把所有集合标号,最后除以集合的阶乘即可,那么考虑容斥,枚举非空集合个数\(i\),可以得到\[\begin{Bmatrix}n\\m\end{Bmatrix}=\frac 1 {m!}\sum_{i=0}^m(-1)^i\binom mi(m-i)^n\]

接下来继续推导自然数幂和。

显然!!\[i^k=\sum_{j=0}^k\left\{\begin{array}{c}{k}\\{j}\end{array}\right\}i^{j\downarrow}\]

继续推导\[\sum_{i=1}^ni^k=\sum_{i=1}^n\sum_{j=0}^k\left\{\begin{array}{c}{k}\\{j}\end{array}\right\}i^{j\downarrow}\]\[=\sum_{i=1}^n\sum_{j=0}^k\left\{\begin{array}{c}{k}\\{j}\end{array}\right\} j!\left(\begin{array}{c}{i}\\{j}\end{array}\right)\]\[=\sum_{j=0}^k\left\{\begin{array}{c}{k}\\{j}\end{array}\right\} j!\sum_{i=j}^n\left(\begin{array}{c}{i}\\{j}\end{array}\right)\]\[=\sum_{j=0}^k\left\{\begin{array}{c}{k}\\{j}\end{array}\right\} j!\binom{n+1}{j+1}\]

很明显,除去预处理第二类斯特林数的复杂度,后面是一样不用做除法的,可以做到\(O(k)\).

那么时间复杂度决定于预处理第二类斯特林数的复杂度。显然可以用\(O(k^2)\)递推。

然而事实上,我们来看看斯特林数的通项公式:\[\begin{Bmatrix}n\\m\end{Bmatrix}=\frac 1 {m!}\sum_{i=0}^m(-1)^i\binom mi(m-i)^n\]

一拼凑,咦~\[\begin{Bmatrix}n\\m\end{Bmatrix}=\sum_{i=0}^m\frac{(-1)^i}{i!}\cdot\frac{(m-i)^n}{(m-i)!}\]

这原来可以写成形如\(\sum_{i=0}^mf(i)*g(m-i)\)的卷积形式。于是第一个多项式的第\(i\)项系数是\(\frac {(-1)^i}{i!}\),另一个多项式的第\(i\)项系数是\(\frac {i^n}{i!}\),卷积后第\(i\)项的系数就是\(\begin{Bmatrix}n\\i\end{Bmatrix}\).

于是愉快的将时间变成了\(O(KlogK)\)

7. 差分表

初三来学习一下差分表吧。

对于任何一个序列\(a_0, a_1, ... , a_n, ...\)我们都可以定义它的差分序列\(\Delta a_0, \Delta a_1, ... ,\Delta a_n, ...\),其中\(\Delta a_i=a_{ i+1 }-a_i\)

类似的,我们可以构造序列\(\{\Delta a_n\}\)的二阶、三阶…\(k\)阶差分序列。 不妨记为\(\{\Delta^2 a_n\},...,\{\Delta^k a_n\}\)

令序列是一个\(p\)次多项式,那么差分表一个很重要且很显然的性质是:\[\forall n>=0, \Delta^{p+1}a_n = 0\]这是由于每次差分都必然会把最高次项消去!

另外一个很重要的性质就是差分表的线性性,即如果\(f_n=k_1g_n+k_2h_n\),那么一定有\[\forall p, n, \Delta^p f_n=k_1\Delta^p g_n + k_2\Delta^p h_n\]

而其最重要的一个性质就是,任何一个\(p\)阶多项式,都必定可以由其差分表的第一条对角线确定。为了证明这个结论,不妨先考虑最简单的情况:

差分表的一条对角线为\(0, ..., 0, 1, 0, ...\),即第一条对角线上只有第\(p\)个位置为\(1\),其他都为\(0\),那么可以写出这个序列的通项公式\[f_n=c(n)(n-1)(n-2)...(n-p+1)\]

代入\(n=p,f_p=1\),得到\[c=\frac{1}{p!}\]所以可以得到\[f_n=\frac{n!}{p!(n-p)!}=\tbinom{n}{p}\]

那么根据差分表的线性性,我们就可以得知\[f_n=\sum_{i=0}^p c_i\tbinom{n}{i}\]
由于\[\sum_{k = 0}^n f(k) = \sum_{k = 0}^p c_k {n + 1 \choose k}\]所以利用差分表,我们可以在\(O(p^2)\)的时间复杂度求解类似于\[\sum_{i=0}^n f_i\]的式子。

回到自然数幂和的问题上,我们把\(f_i=i^k\)代入计算前\(p\)项的值,通过\(p^2\)的时间复杂度处理出差分表的第一条对角线,设这个对角线为\(c(p, 0), c(p, 1), c(p, 2), \dots, c(p, p)\),那么答案就是\[\sum_{k = 0}^p c(p, k){n + 1 \choose k}\]

8. 伯努利数

初三再来学一学伯努利数吧。

根据伯努利数的生成函数定义,可知\[\frac{x}{e^{x}-1} = \sum_{i\geq 0} B_{i}*\frac{x^{i}}{i!}\]

由于\[\frac{x}{e^x-1}·(e^x-1)=x\]

且\[[x^n]\frac{x}{e^x-1}·(e^x-1)=\sum_{i=0}^{n-1}\frac{B_i}{i!}·\frac{1}{(n-i)!}=[n=1]\]

两边都乘上\(n!\)可以得到\[\sum_{i=0}^{n-1}\binom{n}{i}B_i=[n=1]\]

这是伯努利数的一个基本性质。后面会用到。

我们再定义一个多项式\(B_n(t)\)表示\[B_{n}(t) = \sum_{k=0}^{n-1} B_{k} * t^{n-k}*C_{n}^{k}\]

然后我们发现
\[B_n(t+1)-B_n(t)=\sum_{k=0}^{n-1} B_{k}*\lgroup(t+1)^{n-k} - t^{n-k}\rgroup C_{n}^{k}\]\[=\sum_{k=0}^{n-1} B_{k}*(\sum_{i=0}^{n-k-1} C_{n-k}^{i} * t^{i})* C_{n}^{k}\]\[=\sum_{k=0}^{n-1} B_{k}*(\sum_{i=0}^{n-k-1} C_{n-k}^{i} * t^{i}*C_{n}^{k})\]\[=\sum_{i=0}^{n-1}B_k*\sum_{i=0}^{n-k-1}\frac{n!t^i}{i!k!(n-k-i)!}\]\[=\sum_{k=0}^{n-1} B_{k}*(\sum_{i=0}^{n-k-1} C_{n-i}^{k} * t^{i}*C_{n}^{i})\]\[=\sum_{i=0}^{n-1} C_{n}^{i}*t^{i}*\sum_{k=0}^{n-1-i} B_{k}*C_{n-i}^{k}\]

注意到后面只有当\(n-i-1=0\)时值为\(1\),于是\[B_n(t+1)-B_n(t)=n*t^{n-1}\]

然后我们考虑差分,就有\[\sum_{t=0}^{n-1}B_k(t+1)-B_k(t)=k·\sum_{i=0}^{n-1}i^{k-1}\]即\[B_{k+1}(n+1)=(k+1)·\sum_{i=0}^ni^k\]

可得自然数幂和\[\sum_{i=0}^ni^k=\frac{1}{k+1}·\sum_{i=0}^kB_in^{k+1-i}\binom{k+1}{i}\]

问题转化成了求\(B_i\),注意到它的生成函数定义,事实上我们只需要求\[\sum_{i\ge 0}\frac{x^i}{(i+1)!}\]在模\(x^{k+1}\)的逆元即可。

时间复杂度\(O(KLogK)\),当然如果递推的话,也是可以轻松做到\(O(k^2)\)的。

9. 拉格朗日插值法

两大作用:

  1. 快速根据点值逼近原函数.
  2. 取点值对大于\(n\)时唯一确定\(n\)次多项式.
Example

例如对于\[\sum\limits_{i=1}^ni=\frac{n(n+1)}{2}\]

我们知道它的通项公式是二次的,所以我们只需要三个点值对就可以唯一确定这个多项式:\[(1,1),(2,3),(3,6)\]

General method

对于已知的\(n+1\)个点对\((x_0,y_0),(x_1,y_1)...(x_n,y_n)\),求\(n+1\)个函数\(f_i\),使得该函数在\(x_i\)处取得对应的\(y_i\)值,其余\(x_j\)处为\(0\),最后把这\(n+1\)个函数线性结合即可。

\[f_i(x)=\frac{\prod\limits_{j\neq i}(x-x_j)}{\prod\limits_{j\neq i}(x_i-x_j)}*y_i\]\[g(x)=\sum_{i=0}^nf_i(x)\]

Practice

例如我们要求自然数幂和.

各种方法可以证明\(i^k\)的和是\(k+1\)次的, 所以我们只需要给出\(k+2\)个点值表达,就可以求得通项公式.

记\(S(n)=\sum_{i=1}^n i^k\), 则
\[S(n)=\sum_{i=1}^{k+2}y_i\prod_{j=1,i\neq j}^{k+2}\frac {n-x_j}{x_i-x_j}=\sum_{i=1}^{k+2}y_i\frac {\prod_{j=1,i\neq j}^{k+2}(n-j)}{\prod_{j=1,i\neq j}^{k+2}(i-j)}\]

那么时间复杂度就在预处理\(y_i\)上面了, 利用线筛,可以做到\(O(k)\)级别.

后面的那一部分可以预处理,具体的说,就有:\[S(n)=\sum_{i=1}^{k+2}y_ipre[i-1]suf[i+1][(-1)^{k+2-i}(i-1)!(k+2-i)!]^{-1}\]

\(Code\)

时间复杂度:\(O(\frac{k}{ln_k}log_2k)=O(k)\).

另外,注意逆元要预处理,实测:\(k=1e7\)时\(0.9s\),可以说是非常优秀了

#include <bits/stdc++.h> 

#define F(i,a,b) for (int i = a; i <= b; i ++)
#define G(i,a,b) for (int i = a; i >= b; i --)

const int Mo = 998244353, M = 1e6 + 10;

using namespace std;

int l, r, k, m, y[M], z[M], jc[M], suf[M], pre[M];
bool bz[M];

int ksm(int x, int y) {
    int ans = 1;
    for (; y; y >>= 1, x = (1ll * x * x) % Mo)
        if (y & 1)
            ans = (1ll * ans * x) % Mo;
    return ans;
}

void Init() {
    scanf("%d%d%d", &l,&r,&k), y[1] = 1, m = k + 2;
    F(i, 2, m) {
        if (!bz[i])
            z[++ z[0]] = i, y[i] = ksm(i, k);
        F(j, 1, z[0]) {
            if (z[j] * i > m) break;
            bz[z[j] * i] = 1;
            y[z[j] * i] = (1ll * y[z[j]] * y[i]) % Mo;
            if (i % z[j] == 0) break;
        }
    }
    F(i, 2, m)
        y[i] = (y[i - 1] + y[i]) % Mo;
    jc[0] = 1;
    F(i, 1, m)
        jc[i] = 1ll * jc[i - 1] * i % Mo;
    jc[m] = ksm(jc[m], Mo - 2);
    G(i, m - 1, 1)
        jc[i] = 1ll * jc[i + 1] * (i + 1) % Mo;
}

int Solve(int n) {
    pre[0] = suf[m + 1] = 1;
    F(i, 1, m)
        pre[i] = 1ll * pre[i - 1] * (n - i) % Mo;
    G(i, m, 1)
        suf[i] = 1ll * suf[i + 1] * (n - i) % Mo;

    int Ans = 0;
    F(i, 1, m)
        Ans = (Ans + 1ll * y[i] * pre[i - 1] % Mo * suf[i + 1] % Mo * (((k-i+2)&1) ? (-1) : 1) * jc[i - 1] % Mo * jc[k + 2 - i] % Mo) % Mo;
    return Ans;
}

int main() {
    Init();

    printf("%d\n", (Solve(r) - Solve(l - 1) + Mo) % Mo);
}

原文地址:https://www.cnblogs.com/Pro-king/p/10664390.html

时间: 2024-10-11 11:16:48

求解自然数幂和的若干种方法的相关文章

Java 计算中英文长度的若干种方法

在项目开发中经常碰到到输入字符的校验,特别是中英文混合在一起的校验.而为了满足校验的需求,有时需要计算出中英文的长度. 本文将通过几种常用的方法实现长度的计算: <span style="font-size:18px;">import java.io.UnsupportedEncodingException; /** * 中英文校验的处理 * @author a123demi * */ public class EnChValidate { public static vo

footer固定在页面底部的若干种方法

1 <div class="header"><div class="main"></div></div> 2 <div class="container"><div class="main"></div></div> 3 <div class="footer"><div class="

离散数学中使用生成函数求解递推式的一种方法

感觉书(Rosen的离散数学,机械工业的)上的做法有些逆向思维了,没有说明为什么要那样构造,以致大多数同学是背板子上的考场.然而其实用同样的思路我们完全可以使用一种让人可以理解的求解生成函数的方法. 听同学说期末考了两道,我就搞了搞,然鹅缓考时老师换题了一道都没考Orz……我个人的这种做法也就没能施展.想了想还是顺手记在这里吧,万一有人好奇搜了搜呢,万一本篇能帮上点忙呢-可能过几个月我自己都忘了. 先举一例: 总结起来大概套路是: ------就酱------ 原文地址:https://www.

康复计划#3 简单常用的几种计算自然数幂和的方法

本篇口胡写给我自己这样的东西都忘光的残废选手 以及暂时还不会自然数幂和的人- 这里大概给出最简单的几种方法:扰动法(化为递推式),斯特林数(离散微积分),高阶差分(牛顿级数),伯努利数(指数生成函数)- 不同方法的思维难度.普适程度.实现难度.时间复杂度上面都有差异-同时自然数幂和是探究各种求和方法的经典例子,了解多一点它的做法对于处理各种求和问题是有所帮助的- 问题:求$\sum_{k=0}^{n} k^t$,其中$t \in \mathbb{N}$是一个常数.要求求解的时间复杂度与$n$无关

各类求自然数幂和方法

高斯消元 我们知道: \[\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,我们可以拿这个猜

三种方法求解约瑟夫环问题

约瑟夫环是一个数学的应用问题:已知n个人(以编号1,2,3...n分别表示)围坐在一张圆桌周围.从编号为k的人开始报数,数到m的那个人出列:他的下一个人又从1开始报数,数到m的那个人又出列:依此规律重复下去,直到圆桌周围的人全部出列. 方法1:使用stl::list模拟环形链表,参考剑指offer 代码: #include <iostream> #include <list> using namespace std; int lastNumber(unsigned int n,un

OC中动态创建可变数组的问题.有一个数组,数组中有13个元素,先将该数组进行分组,每3个元素为一组,分为若干组,最后用一个数组统一管理这些分组.(要动态创建数组).两种方法

<span style="font-size:24px;">//////第一种方法 // NSMutableArray *arr = [NSMutableArray array]; // for (int i = 0; i < 13; i ++) { // [arr addObject:[NSString stringWithFormat:@"lanou%d",i + 1]]; // } // NSLog(@"%@",arr);

Matlab求解线性方程组Ax=b的几种常见方法Matlab求解线性方程组Ax=b的几种常见方法

原文:http://blog.csdn.net/stzh_bk/article/details/70983856 例如方程组: 法1:左除法 >> A=[3 1 -1;1 2 4;-1 4 5];b=[3.6;2.1;-1.4]; >> x=A\b x = 1.4818 -0.4606 0.3848 法2:求逆法 >> A=[3 1 -1;1 2 4;-1 4 5];b=[3.6;2.1;-1.4]; >> x=inv(A)*b x = 1.4818 -0.

OC动态创建的问题变量数组.有数组,在阵列13要素,第一个数据包阵列,每3元素为一组,分成若干组,这些数据包的统一管理。最后,一个数组.(要动态地创建一个数组).两种方法

<span style="font-size:24px;">//////第一种方法 //        NSMutableArray *arr = [NSMutableArray array]; //        for (int i = 0; i < 13; i ++) { //            [arr addObject:[NSString stringWithFormat:@"lanou%d",i + 1]]; //