題目大意

有函數\(f(x)\)，\(f(0)=0,f(1)=1,f(x)=f(x-1)+f(x-2)\)
\(t\)(\(t\leq1000\))組詢問，每次給定\(n,m\)(\(n,m\leq10^6\))，求：
\[\prod_{i=1}^{n}\prod_{j=1}^{m}f(gcd(i,j))\]

題解

這個人（點這里）講得很清楚\(\color{white}{\text{shing太強了}}\)
假設\(n\)是\(n,m\)中較小的那個，講一些在得出答案\(=\prod_{d=1}^{n}f(d)^{\sum_{k=1}^{\lfloor \frac{n}ze8trgl8bvbq\rfloor}{\mu(k)\lfloor \frac{n}{dk}\rfloor\lfloor \frac{m}{dk}\rfloor}}\)之后的事
設\(p=dk\)則有答案\(=\prod_{d=1}^{n}f(d)^{\sum_{d\mid p}^{n}{\mu(\frac{p}ze8trgl8bvbq)\lfloor \frac{n}{p}\rfloor\lfloor \frac{m}{p}\rfloor}}\)
把冪次上的\(\sum\)拿下來變成\(\prod\)，答案\(=\prod_{d=1}^{n}\prod_{d\mid p}^{n}f(d)^{{\mu(\frac{p}ze8trgl8bvbq)\lfloor \frac{n}{p}\rfloor\lfloor \frac{m}{p}\rfloor}}=\prod_{d=1}^{n}(\prod_{d\mid p}^{n}f(d)^{\mu(\frac{p}ze8trgl8bvbq)} )^{\lfloor \frac{n}{p}\rfloor\lfloor \frac{m}{p}\rfloor}\)
括號里的部分\(g(d)=\prod_{d\mid p}^{n}f(d)^{\mu(\frac{p}ze8trgl8bvbq)}\)與\(n,m\)無關可以預處理，將每個\(f(a)^{\mu(b)}\)乘到\(g(a*b)\)，預處理的時間復雜度為\(\frac{n}{1}+\frac{n}{2}+\frac{n}{3}+...+\frac{n}{n}\)，大約是\(n*log n\)
接下來數論分塊計算\(\prod_{d=1}^{n}g(d)^{\lfloor \frac{n}{p}\rfloor\lfloor \frac{m}{p}\rfloor}\)就行了，這部分的時間復雜度是\(\Theta(T(\sqrt n +\sqrt m))\)

代碼

#include<algorithm> #include<cmath> #include<cstdio> #include<cstdlib> #include<cstring> #include<ctime> #include<iomanip> #include<iostream> #include<map> #include<queue> #include<set> #include<stack> #include<vector> #define rep(i,x,y) for(register int i=(x);i<=(y);++i) #define dwn(i,x,y) for(register int i=(x);i>=(y);--i) #define maxn 1000010 #define lim 1000000 #define LL long long using namespace std; int read() {int x=0,f=1;char ch=getchar();while(!isdigit(ch)&&ch!='-')ch=getchar();if(ch=='-')f=-1,ch=getchar();while(isdigit(ch))x=(x<<1)+(x<<3)+ch-'0',ch=getchar();return x*f; } void write(int x) {if(x==0){putchar('0'),putchar('\n');return;}int f=0;char ch[20];if(x<0)putchar('-'),x=-x;while(x)ch[++f]=x%10+'0',x/=10;while(f)putchar(ch[f--]);putchar('\n');return; } const LL mod=1e9+7; int n,m,t,f[maxn],pf[maxn],g[maxn],pg[maxn],p[maxn],no[maxn],cnt,mu[maxn]; int mul(int x,LL y){int res=1;while(y){if(y&1ll)res=(LL)res*(LL)x%mod;x=(LL)x*(LL)x%mod,y>>=1;}return res;} int main() {no[1]=mu[1]=1;rep(i,2,lim){if(!no[i])p[++cnt]=i,mu[i]=-1;for(int j=1;j<=cnt&&i*p[j]<=lim;j++){no[i*p[j]]=1;if(i%p[j]==0){mu[i*p[j]]=0;break;}else mu[i*p[j]]=-mu[i];}}f[1]=pf[1]=g[1]=1;rep(i,2,lim)f[i]=(f[i-1]+f[i-2])%mod,pf[i]=mul(f[i],mod-2),g[i]=1;rep(i,1,lim){for(int j=i;j<=lim;j+=i){if(mu[j/i]==1)g[j]=(LL)g[j]*(LL)f[i]%mod;else if(mu[j/i]==-1)g[j]=(LL)g[j]*(LL)pf[i]%mod;}}g[0]=pg[0]=1;rep(i,1,lim)g[i]=(LL)g[i-1]*(LL)g[i]%mod,pg[i]=mul(g[i],mod-2);t=read();while(t--){n=read(),m=read();if(n>m)swap(n,m);int ans=1;for(int l=1,r=0;l<=n;l=r+1){r=min(n/(n/l),m/(m/l));ans=(LL)ans*mul((LL)g[r]*(LL)pg[l-1]%mod,(LL)(n/l)*(LL)(m/l))%mod;}write(ans);}return 0; }

轉載于:https://www.cnblogs.com/xzyf/p/10443257.html

總結

以上是生活随笔為你收集整理的并不对劲的bzoj4816:loj2000:p3704[SDOI2017]数字表格的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。