多校冲刺NOIP模拟19 - 乘法——倍增、DP

偶耶XJX

于 2021-10-31 12:25:24 发布

阅读量210

点赞数

分类专栏：信息竞赛解题文章标签：算法数学倍增动态规划

本文链接：https://blog.csdn.net/weixin_43960287/article/details/121061712

版权

信息竞赛解题专栏收录该内容

195 篇文章 5 订阅

订阅专栏

此题不提供链接

题目描述

求 $n!$ （ $n<2^{64}$ ）转成十六进制后，去掉末尾的零的最后 16 位。

题解

设 $n!=a\cdot2^b,a\,\,is\,\,odd$ ，显然答案应该是 $a\cdot2^{b\bmod 4}$ 。

由于保留十六进制最后 16 位，相当于二进制下保留 64 位，取模 $2^{64}$ ，所以全程用 unsigned long long 自然溢出即可。剩下的问题是如何求出 n 以内所有数去除质因子 2 之后相乘。

首先，可以把所有奇数相乘，然后剩下的数除以 2 过后就变成了 $\{1,2,3,...,\lfloor\frac{n}{2}\rfloor\}$ ，再次取所有奇数相乘，把剩下的数除以 2…

所以，我们设答案为 $f (n)$ ，设 $g (n)$ 为 n 以内所有奇数相乘，那么有
$f(n)=g(n)+f(\lfloor\frac{n}{2}\rfloor)$
由此我们可以想到倍增做法，只需要每次能较快求出 $g (n)$ 即可。

由于 $g(n)=\prod_{i=0}^{\frac{n-1}{2}}(2i+1)$ ，我们把括号拆开，就变为所有由若干个 $2 i$ 相乘的项相加。由于取模 $2^{64}$ ，所以我们最多只需要关注不超过 63 个 $2 i$ 相乘的项，再大则必定为0。所以对于每个 $g (n)$ ，我们只需要求出所有 $dp(n,k)(0\le k<64)$ 表示选 k 个 $2 i$ 相乘的和，然后加起来。

再一次转换， $d p (n, k)$ 怎么求？还是考虑用倍增。假设我们已经求得所有 $dp(\frac{n}{2},k)$ ，相当于知道了前面一半的 $\sum\prod_{...}2i$ 。设后面一半为 $dp'(\frac{n}{2},k)=\sum\prod_{...}(2i+x)$ ，那么有
$dp(n,k)=\sum_{i=0}^kdp(\frac{n}{2},i)\cdot dp'(\frac{n}{2},k-i)$
很像卷积？对啊，我也不知道为什么，可它就是一个普通的状态转移方程。

观察 $dp'(\frac{n}{2},k)$ 的式子，发现里面这个括号又可以展开，变为若干个 $2 i$ 与 $x$ 的若干次方相乘的和。也就是说，我们可以利用 $dp(\frac{n}{2},k)$ 来求出 $dp'(\frac{n}{2},k)$ ：
$dp'(\frac{n}{2},k)=\sum_{i=0}^{k}dp(\frac{n}{2},i)\cdot x^{k-i}\,\,(?)$
搞出来一个貌似正确的式子，可是它连 $n = 5$ 都会求错。说明式子肯定不这么简单，后面还要乘上一个系数。经过不断尝试演算，我试算出来发现是 ${m-i\choose k-i}$ ，其中 m 是 $\frac{n}{2}$ 以内奇数的个数。
$dp'(\frac{n}{2},k)=\sum_{i=0}^{k}dp(\frac{n}{2},i)\cdot x^{k-i}\cdot{m-i\choose k-i}$
其实也很好理解，因为 ${m-i\choose k-i}=\frac{{m\choose k}{k\choose i}}{{m\choose i}}$ 。

怎么求组合数呢？由于 $k - i$ 很小，所以可以考虑用同一行的组合数递推。但是需要模意义下除法，模数不为质数不好求逆元怎么办？

由于组合数是整数，所以递推的过程中保证分母之积整除分子之积，分子的质因子2的次数肯定大于等于分母的。所以我们可以先把所有分子分母的2的都取出来，最后再乘回去。把2取出来后，剩下的分母必定与 $2^{64}$ 互质，直接预处理求逆元即可。

最后有一些小细节，比如奇数个数翻倍后多1或少1，需要用背包把dp再推一步或倒推一步。

两个倍增合在一起处理，复杂度可达到 $O(\log^3n)$ 。我没有预处理 $x$ 的次幂和2的次幂，直接使用快速幂，所以复杂度是 $O(\log^3n\log\log n)$ 。

代码

自己推了好久，发现居然就是 $\rm O\color{red}neInDark$ 考场上的做法！我还是被 $\rm O\color{red}neInDark$ 远远地甩在后面！

#include<cstdio>//JZM yyds!!
#include<cstring>
#include<iostream>
#include<algorithm>
#include<cmath>
#include<ctime>
#include<vector>
#include<queue>
#include<stack>
#include<map>
#include<set>
#define ll long long
#define uns unsigned
#define MOD 
#define MAXN 
#define INF 1e17
#define IF (it->first)
#define IS (it->second)
using namespace std;
inline ll read(){
	ll x=0;bool f=1;char s=getchar();
	while((s<'0'||s>'9')&&s>0)f^=(s=='-'),s=getchar();
	while(s>='0'&&s<='9')x=(x<<1)+(x<<3)+(s^48),s=getchar();
	return f?x:-x;
}
int pt[30],lp;
inline void print(ll x,char c='\n'){
	if(x<0)putchar('-'),x=-x;
	pt[lp=1]=x%10;
	while(x>9)x/=10,pt[++lp]=x%10;
	while(lp)putchar(pt[lp--]^48);
	putchar(c);
}
inline uns ll lowbit(uns ll x){return x&-x;}
uns ll n;
char p16[16]={'0','1','2','3','4','5','6','7','8','9','A','B','C','D','E','F'};
inline void print16(uns ll x,char c='\n'){
	pt[lp=1]=x&15ull;
	while(x>15ull)x>>=4,pt[++lp]=x&15ull;
	while(lp)putchar(p16[pt[lp--]]);
	putchar(c);
}
inline uns ll ksm(uns ll a,uns ll b){
	uns ll res=1;
	for(;b;b>>=1,a*=a)if(b&1)res*=a;
	return res;
}
uns ll inv[105];
const int w=64;
struct itn{
	uns ll f[w];
	int n;
	itn(){memset(f,0,sizeof(f)),n=0;}
};
uns ll ans;
inline int Log2(double x){
	if(x==1)return 0;
	return ((*(unsigned long long*)&x>>52)&1023)+1;
}
inline int coun(uns ll n){return Log2(lowbit(n));}
inline uns ll del(uns ll n){return n/lowbit(n);}
inline itn C(uns ll n){
	itn res;
	res.n=min(n+1ull,w+0ull),res.f[0]=1;
	uns ll k=res.n-1;
	itn cnt;cnt.n=res.n;
	for(uns ll i=1;i<=k;i++){
		res.f[i]=res.f[i-1]*del(n-i+1)*inv[del(i)];
		cnt.f[i]=cnt.f[i-1]+coun(n-i+1)-coun(i);
	}
	for(uns ll i=1;i<=k;i++)res.f[i]*=ksm(2,cnt.f[i]);
	return res;
}
inline itn solve(uns ll n){
	if(n==1){
		itn res;res.f[0]=1,res.n=1;
		return res;
	}
	uns ll m=(n>>1),ad=m-(~m&1)+1;
	if(n&1)m++;
	uns ll c=((n>>1)+1)>>1,d;
	itn a=solve(n>>1),b,res;
	b.n=min(a.n+1,w),res.n=min(a.n+b.n-1,w);
	itn ch[b.n];
	for(int i=0;i<b.n;i++)ch[i]=C(c-i);
	for(int i=b.n-1;~i;i--)for(int j=i;~j;j--)
		b.f[i]+=a.f[j]*ksm(ad,i-j)*ch[j].f[i-j];
	for(int i=res.n-1;~i;i--)for(int j=i;~j;j--)
		res.f[i]+=a.f[j]*b.f[i-j];
	c<<=1;
	while(m>c){
		c++,d=(c<<1)-1ull;
		res.n=min(res.n+1,w);
		for(int i=res.n-1;i>0;i--)res.f[i]+=res.f[i-1]*(d-1ull);
	}while(m<c){
		d=(c<<1)-1ull,c--;
		for(int i=1;i<res.n;i++)res.f[i]-=res.f[i-1]*(d-1ull);
		if(res.n>0&&!res.f[res.n-1])res.n--;
	}
	uns ll cg=0;
	for(int i=0;i<res.n;i++)cg+=res.f[i];
	ans*=cg;
	return res;
}
signed main()
{
	freopen("multiplication.in","r",stdin);
	freopen("multiplication.out","w",stdout);
	for(int i=1;i<=65;i+=2)inv[i]=ksm(i,(1ull<<63)-1);
	for(int T=read();T--;){
		scanf("%llu",&n);
		uns ll m=n;
		int e=0;
		while(m>1)e+=(m>>1)&3ull,e&=3,m>>=1;
		ans=1;
		solve(n);
		while(e--)ans<<=1;
		print16(ans);
	}
	return 0;
}