雅礼集训2019Day1 T3—Math(矩阵快速幂优化dp)

描述

给出 n n n m m m x x x,你需要求出下列式子的值:

∑ ( ∑ i ≤ m k i ) = n ∏ i ≤ m s i n ( k i ∗ x ) \sum_{(\sum_{i≤m}{k_i})=n}\prod_{i≤m}sin(k_i∗x) (imki)=nimsin(kix)
其中 k i k_i ki为正整数,由于答案非常大,你只需要输出答案(保证不为 0)的正负(如果是负数输出负号,否则输出正号)和从左往右第一个非 0 数位上的数字即可。


Solution

设?(?, ?)表示 k 的序列大小为 m,k 之和为 n 时的答案。
我们对 k m k_m km进行分类讨论来得到?(?, ?)的转移:

  1. k m k_m km=1,则?(? − 1, ? − 1) ∗ sin(?) −→ ?(?, ?)
  2. k m k_m km> 1,考虑对???(? ∗ ?)进行下列变换:
    ???(? ∗ ?) = ???(?) ∗ ???((? − 1) ∗ ?) + ???(?) ∗ ???((? − 1) ∗ ?)
    又???(?) ∗ ???(? ∗ ?) = ???(?) ∗ (???(?) ∗ ???((? − 1) ∗ ?) − ???(?) ∗ ???((? − 1) ∗ ?))
    = ???(?) ∗ ???(?) ∗ ???((? − 1) ∗ ?) + (???2
    (?) − 1) ∗ ???((? − 1) ∗ ?)
    = ???(?) ∗ (???(?) ∗ ???((? − 1) ∗ ?) + ???(?) ∗ ???((? − 1) ∗ ?)) − ???((? − 1) ∗ ?)
    = ???(?) ∗ ???(? ∗ ?) − ???((? − 1) ∗ ?)
    所以???(? ∗ ?) = 2 ???(?) ∗ ???((? − 1) ∗ ?) − ???((? − 2) ∗ ?)
    则这部分的贡献为:
    2?(? − 1, ?) ∗ ???(?) − ?(? − 2, ?) ⟶ \longrightarrow ?(?, ?)
    显然可以用矩阵乘法优化这个过程。
    所以复杂度为 O ( m 3 l o g n ) O(m^3logn) O(m3logn)

#include<bits/stdc++.h>
using namespace std;
const int RLEN=1<<20|1;
#define ll long long
inline char gc(){
	static char ibuf[RLEN],*ib,*ob;
	(ib==ob)&&(ob=(ib=ibuf)+fread(ibuf,1,RLEN,stdin));
	return (ib==ob)?EOF:*ib++;
}
inline int read(){
	char ch=gc();
	int res=0,f=1;
	while(!isdigit(ch))f^=ch=='-',ch=gc();
	while(isdigit(ch))res=(res+(res<<2)<<1)+(ch^48),ch=gc();
	return f?res:-res;
}
const int N=62;
int siz;
struct mat{
	double a[N][N];
	mat(){memset(a,0,sizeof(a));}
	inline const double *const operator [](const int &oset)const{return a[oset];}
	inline double *const operator [](const int &oset){return a[oset];}
	friend inline mat operator *(const mat &a,const mat &b){
		mat c;
		for(int i=1;i<=siz;i++)
			for(int j=1;j<=siz;j++)
				for(int k=1;k<=siz;k++)	
					c[i][k]+=a[i][j]*b[j][k];
		return c;
	}
}A,B;
inline mat ksm(mat a,int b,mat res){
	for(;b;b>>=1,a=a*a)if(b&1)res=res*a;
	return res;
} 
int T;
int n,m;
double x;
int main(){
	cin>>T;
	while(T--){
		cin>>m>>n>>x;siz=m*2;
		A=mat(),B=mat();
		for(int i=m+1;i<=m*2;i++){
			A[i][i-m]=1;
			A[i-m][i]=-1;
			A[i][i]=2*cos(x);
			if(i<m*2)A[i][i+1]=sin(x);
		}
		B[1][1]=sin(x);
		B[1][m+1]=sin(2*x);
		B[1][m+2]=sin(x)*sin(x);
		B=ksm(A,n-1,B);
		double res=B[1][m];
		if(res<0)putchar('-'),res=-res;
		else putchar('+');
		while(res<1)res*=10;
		while(res>=10)res/=10;
		putchar((int)(floor(res))^48);puts("");
	}
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值