UVA766 Sum of powers

Luogu 链接
UVA 链接
Virtual Judge 链接

题意

给定 \(1\) 个非负整数 \(k\le 20\)，求 \(\displaystyle S_k(n)=\sum_{i=1}^ni^k\) 的求和公式。

多测，每组数据输出之间有 \(1\) 个空行（不得有多余的空格和空行）。

对于每组数据，输出 \(k+2\) 个整数

\[M,a_{k+1},a_k,\dots,a_1,a_0 \]

满足

\[\displaystyle S_k(n)=\frac{1}{M}\sum_{i=0}^{k+1}a_in^i \]

并且

\[\gcd(a_{k+1},a_k,\dots,a_1,a_0)=1 \]

思路

让我们来推导下如何用 \(S_0(n),S_1(n),\dots,S_{k-1}(n)\) 来表示 \(S_k(n)\)。

\(\displaystyle(i+1)^{k+1}-i^{k+1}=\sum_{j=0}^kC_{k+1}^ji^j\)

\(\displaystyle(n+1)^{k+1}-1=\sum_{i=1}^n[(i+1)^{k+1}-i^{k+1}]\)

\(\begin{aligned}\therefore(n+1)^{k+1}-1&=\sum_{i=1}^n\sum_{j=0}^kC_{k+1}^ji^j\\&=\sum_{j=0}^k(C_{k+1}^j\sum_{i=1}^ni^j)\\&=\sum_{i=0}^kC_{k+1}^iS_i(n)\\&=\sum_{i=0}^{k-1}C_{k+1}^iS_i(n)+(k+1)S_k(n)\end{aligned}\)

\(\displaystyle\therefore S_k(n)=\frac{\displaystyle(n+1)^{k+1}-1-\sum_{i=0}^{k-1}C_{k+1}^iS_i(n)}{k+1}\)

实现时，需要实现分数类和多项式类。

程序

AC 记录

#include<cstdlib>
#include<algorithm>
#include<cstdio>
#include<cstring>
#include<queue>
#include<cstdio>
#include<iostream>
#include<vector>
#include<map>
#include<cmath>
#include<iomanip>
#include<string>
#include<stack>
#define re register
#define ll long long
#define ull unsigned long long
#define vl __int128
#define ld long double
#define LL 2e18
#define INT 1e9
#define INF 0x3f3f3f3f
#define lb(x) (x&(-x))
#ifdef __linux__
#define gc getchar_unlocked
#define pc putchar_unlocked
#else
#define gc _getchar_nolock
#define pc _putchar_nolock
#endif
int T=1;
using namespace std;
inline bool blank(const char x){return !(x^32)||!(x^10)||!(x^13)||!(x^9);}
template<typename Tp>inline void read(Tp &x){x=0;re bool z=true;re char a=gc();for(;!isdigit(a);a=gc())if(a=='-')z=false;for(;isdigit(a);a=gc())x=(x<<1)+(x<<3)+(a^48);x=(z?x:~x+1);}
inline void read(double &x){x=0.0;re bool z=true;re double y=0.1;re char a=gc();for(;!isdigit(a);a=gc())if(a=='-')z=false;for(;isdigit(a);a=gc())x=x*10+(a^48);if(a!='.')return x=z?x:-x,void();for(a=gc();isdigit(a);a=gc(),y/=10)x+=y*(a^48);x=(z?x:-x);}
inline void read(ld &x){x=0.0;re bool z=true;re ld y=0.1;re char a=gc();for(;!isdigit(a);a=gc())if(a=='-')z=false;for(;isdigit(a);a=gc())x=x*10+(a^48);if(a!='.')return x=z?x:-x,void();for(a=gc();isdigit(a);a=gc(),y/=10)x+=y*(a^48);x=(z?x:-x);}
inline void read(char &x){for(x=gc();blank(x)&&(x^-1);x=gc());}
inline void read(char *x){re char a=gc();for(;blank(a)&&(a^-1);a=gc());for(;!blank(a)&&(a^-1);a=gc())*x++=a;*x=0;}
inline void read(string &x){x="";re char a=gc();for(;blank(a)&&(a^-1);a=gc());for(;!blank(a)&&(a^-1);a=gc())x+=a;}
template<typename T,typename ...Tp>inline void read(T &x,Tp &...y){read(x),read(y...);}
template<typename T>inline void read(T *begin,T *end){re T *i;if(begin<end)for(i=begin;i<end;++i)read(*i);else for(i=begin-1;i>=end;--i)read(*i);}
template<typename Tp>inline void write(Tp x){if(!x)return pc(48),void();if(x<0)pc('-'),x=~x+1;re int len=0;re char tmp[64];for(;x;x/=10)tmp[++len]=x%10+48;while(len)pc(tmp[len--]);}
inline void write(const double x){re int a=6;re double b=x,c=b;if(b<0)pc('-'),b=-b,c=-c;re double y=5*powl(10,-a-1);b+=y,c+=y;re int len=0;re char tmp[64];if(b<1)pc(48);else for(;b>=1;b/=10)tmp[++len]=floor(b)-floor(b/10)*10+48;while(len)pc(tmp[len--]);pc('.');for(c*=10;a;a--,c*=10)pc(floor(c)-floor(c/10)*10+48);}
inline void write(const ld x){re int a=6;re ld b=x,c=b;if(b<0)pc('-'),b=-b,c=-c;re ld y=5*powl(10,-a-1);b+=y,c+=y;re int len=0;re char tmp[64];if(b<1)pc(48);else for(;b>=1;b/=10)tmp[++len]=floor(b)-floor(b/10)*10+48;while(len)pc(tmp[len--]);pc('.');for(c*=10;a;a--,c*=10)pc(floor(c)-floor(c/10)*10+48);}
inline void write(const pair<int,double>x){re int a=x.first;if(a<7){re double b=x.second,c=b;if(b<0)pc('-'),b=-b,c=-c;re double y=5*powl(10,-a-1);b+=y,c+=y;re int len=0;re char tmp[64];if(b<1)pc(48);else for(;b>=1;b/=10)tmp[++len]=floor(b)-floor(b/10)*10+48;while(len)pc(tmp[len--]);a&&(pc('.'));for(c*=10;a;a--,c*=10)pc(floor(c)-floor(c/10)*10+48);}else printf("%.*lf",a,x.second);}
inline void write(const pair<int,ld>x){re int a=x.first;if(a<7){re ld b=x.second,c=b;if(b<0)pc('-'),b=-b,c=-c;re ld y=5*powl(10,-a-1);b+=y,c+=y;re int len=0;re char tmp[64];if(b<1)pc(48);else for(;b>=1;b/=10)tmp[++len]=floor(b)-floor(b/10)*10+48;while(len)pc(tmp[len--]);a&&(pc('.'));for(c*=10;a;a--,c*=10)pc(floor(c)-floor(c/10)*10+48);}else printf("%.*Lf",a,x.second);}
inline void write(const char x){pc(x);}
inline void write(const bool x){pc(x?49:48);}
inline void write(char *x){fputs(x,stdout);}
inline void write(const char *x){fputs(x,stdout);}
inline void write(const string &x){fputs(x.c_str(),stdout);}
template<typename T,typename ...Tp> inline void write(T x,Tp ...y){write(x),write(y...);}
template<typename T>inline void write(T *begin,T *end,const char c=' '){re T *i;for(i=begin;i<end;++i)write(*i,c);}
template<typename T>inline void init(T *begin,T *end,const T& val=T()){re T* i;for(i=begin;i<end;++i)*i=val;}
template<typename T>inline T max(T *begin,T *end){re T *ans,*i;for(i=begin;i<end;++i)if(i==begin||*ans<*i)ans=i;return *ans;}
template<typename T>inline T min(T *begin,T *end){re T *ans,*i;for(i=begin;i<end;++i)if(i==begin||*i<*ans)ans=i;return *ans;}
template<typename T>inline T calc_sum(T *begin,T *end,const T& val=T()){re T ans=val,*i;for(i=begin;i<end;++i)ans+=*i;return ans;}
template<typename T>inline bool is_equal(T *begin,T *end,const T& val=T()){re T *i;for(i=begin;i<end;++i)if(*i!=val)return false;return true;}

ll mod=0;
const int MAXN=0;
const int N=MAXN+10;
//#define DEBUG
#define more_text

struct frac{
	vl p,q;
	frac():p(0ll),q(1ll){}
	frac(vl P,vl Q):p(P),q(Q){}
	frac(vl P):p(P),q(1ll){}
	~frac(){}
};
inline void write(const frac& f){if(f.q!=1)write(f.p,'/',f.q);else write(f.p);}
frac operator+(frac a,frac b){vl p=a.p*b.q+b.p*a.q,q=a.q*b.q;vl gcd=__gcd(p,q);p/=gcd,q/=gcd;if(q<0)p=-p,q=-q;return frac(p,q);}
frac operator-(frac a,frac b){b.p=-b.p;return a+b;}
frac operator*(frac a,frac b){vl p=a.p*b.p,q=a.q*b.q;vl gcd=__gcd(p,q);p/=gcd,q/=gcd;if(q<0)p=-p,q=-q;return frac(p,q);}
frac operator/(frac a,frac b){swap(b.p,b.q);return a*b;}

struct poly{
	frac fact[22];
	poly(){}
	~poly(){}
};
inline void write(poly P){bool flag=false;for(int i=21;i>=1;--i){if(P.fact[i].p==0)continue;if(flag&&P.fact[i].p>0)write('+');if(P.fact[i].p==-1)write('-');if((P.fact[i].p>0&&(P.fact[i].p!=1||P.fact[i].q!=1))||(P.fact[i].p<0&&(P.fact[i].p!=-1||P.fact[i].q!=1)))write(P.fact[i]);write('n');if(i!=1)write("^",i);flag=true;}if(flag&&P.fact[0].p>0)write('+');if(P.fact[0].p!=0||(!flag&&P.fact[0].p==0))write(P.fact[0]);}
poly operator+(poly a,poly b){poly ans;for(int i=0;i<=21;++i)ans.fact[i]=a.fact[i]+b.fact[i];return ans;}
poly operator-(poly a,poly b){poly ans;for(int i=0;i<=21;++i)ans.fact[i]=a.fact[i]-b.fact[i];return ans;}
poly operator*(frac a,poly b){poly ans;for(int i=0;i<=21;++i)ans.fact[i]=a*b.fact[i];return ans;}
poly operator*(poly a,frac b){poly ans;for(int i=0;i<=21;++i)ans.fact[i]=b*a.fact[i];return ans;}
poly operator/(poly a,frac b){poly ans;for(int i=0;i<=21;++i)ans.fact[i]=a.fact[i]/b;return ans;}
poly operator*(poly a,poly b){poly ans;for(int i=0;i<=21;++i)for(int j=0,k=i;j<=i;++j,--k)ans.fact[i]=ans.fact[i]+a.fact[j]*b.fact[k];return ans;}

vl C(vl a,vl b){vl ans=1;for(vl i=a+1;i<=b;++i)ans*=i;for(vl i=1;i<=b-a;++i)ans/=i;return ans;}

vl k;

poly S[21],n_1,pow_poly;

void print(poly P){
	vl M=P.fact[0].q;
	for(int i=1;i<=k+1;++i)M=M*P.fact[i].q/__gcd(M,P.fact[i].q);
	
	write(M);
	
	for(int i=k+1;i>=0;--i)write(' ',M*P.fact[i]);
}

void SOLVE(int step){
	S[0].fact[1]=1;n_1.fact[1]=n_1.fact[0]=1;pow_poly=n_1;
	
	read(k);
	for(vl i=1;i<=k;++i){
		pow_poly=pow_poly*n_1;
		S[i]=pow_poly;S[i].fact[0]=S[i].fact[0]-frac(1,1);
		for(vl j=0;j<i;++j)S[i]=S[i]-C(j,i+1)*S[j];
		S[i]=S[i]/(i+1);
	}
	
	print(S[k]);write('\n');
	if(step<T-1)write('\n');
}
/*
Input:

Output:

Outline:

*/
int main(){
	#ifdef DEBUG
	freopen("test.in","r",stdin);
	freopen("test.out","w",stdout);
	#endif
	#ifdef more_text
	read(T);
	#endif 
	for(int i=0;i<T;++i)SOLVE(i);
	#ifdef DEBUG
	fclose(stdin);fclose(stdout);
	#endif
	return 0;
}