luogu1587 [NOI2016]循环之美

题目链接:luogu1587

首先是题目中的“纯循环小数”让人感觉十分清奇

对于一个分数\(\frac{x}{y}\),要满足条件首先得使得\(gcd(x,y)=1\),其次因为纯循环必然存在一个循环节的长度\(l\)使得将小数点右移\(l\)位之后两个数的小数部分相同,写出来就是这个样子
\[ \frac{x}{y}-\lfloor\frac{x}{y}\rfloor=\frac{xk^l}{y}-\lfloor\frac{xk^l}{y}\rfloor \]

两边去分母可得
\[ x-y\lfloor\frac{x}{y}\rfloor=xk^l-y\lfloor\frac{xk^l}{y}\rfloor \]
移项
\[ xk^l-x=y(\lfloor\frac{xk^l}{y}\rfloor-\lfloor\frac{x}{y}\rfloor) \]
所以最后就是\(xk^l\equiv x(mod\ y)\),由于\(gcd(x,y)=1\),所以两边同时除以\(x\)得到\(k^l\equiv 1(mod\ y)\),很显然有\(gcd(k,y)=1\)
所以最后题目求的就是
\[ \sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=1][gcd(j,k)=1] \]
如果两个都同时进行反演的话显然很难处理,于是考虑先做一个
\[ \begin{aligned} &\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=1][gcd(j,k)=1]\\=&\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=1]\sum_{d|gcd(j,k)}\mu(d)\\ =&\sum_{d|k}\mu(d)\sum_{i=1}^n\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}[gcd(i,jd)=1]\\ =&\sum_{d|k}\mu(d)\sum_{i=1}^n\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}[gcd(i,j)=1][gcd(i,d)=1] \end{aligned} \]
这有什么意义呢?注意到前后两个式子形式大体一致,于是可以记\(f(n,m,k)=\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=1][gcd(j,k)=1]\),则\(f(n,m,k)=\sum_{d|k}\mu(d)f(\lfloor\frac{m}{d}\rfloor,n,d)\)
边界条件的话首先\(f(0,m,k)=f(n,0,k)=0\)
接下来\(f(n,m,1)=\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=1]\),这个随便反演一下就有\(\sum_{d=1}^n\mu(i)\lfloor\frac{n}{d}\rfloor\lfloor\frac{m}{d}\rfloor\),直接数论分块,但是要注意\(n\leq 10^9\),于是考虑杜教筛求\(\mu\)的和
一个小技巧就是我们只会去计算\(\mu(d)\neq 0\)时的值,其它的直接剪枝即可

#include<iostream>
#include<string>
#include<string.h>
#include<stdio.h>
#include<algorithm>
#include<math.h>
#include<vector>
#include<queue>
#include<stack>
#include<map>
#include<set>
using namespace std;
typedef long long ll;
#define lowbit(x) (x)&(-x)
#define sqr(x) (x)*(x)
#define rep(i,a,b) for (register int i=a;i<=b;i++)
#define per(i,a,b) for (register int i=a;i>=b;i--)
#define fir first
#define sec second
#define maxd 1000000007
#define eps 1e-8
#define pb(a) push_back(a)
#define mp(a,b) make_pair(a,b)
const int N=10000000;
struct node{
    int n,m,k;
};
bool operator <(node p,node q)
{
    if (p.n!=q.n) return p.n<q.n;
    else if (p.m!=q.m) return p.m<q.m;
    else return p.k<q.k;
}
map<node,ll> ans;
map<int,ll> sumu;
int n,m,k,tot=0,mu[N+10],pri[1001000],len;
bool nopri[N+10];
ll sum[N+10];
vector<int> fac;

int read()
{
    int x=0,f=1;char ch=getchar();
    while ((ch<'0') && (ch>'9')) {if (ch=='-') f=-1;ch=getchar();}
    while ((ch>='0') && (ch<='9')) {x=x*10+(ch-'0');ch=getchar();}
    return x*f;
}

void sieve()
{
    mu[1]=1;
    rep(i,2,N)
    {
        if (!nopri[i]) {pri[++tot]=i;mu[i]=-1;}
        int j;
        for (j=1;j<=tot && i*pri[j]<=N;j++)
        {
            nopri[i*pri[j]]=1;
            if (i%pri[j]==0) break;
            else mu[i*pri[j]]-=mu[i];
        }
    }
    rep(i,1,N) sum[i]=sum[i-1]+mu[i];
}

ll query(int x)
{
    if (x<=N) return sum[x];
    if (sumu[x]) return sumu[x];
    ll ans=1;int l,r;
    for (l=2;l<=x;l=r+1)
    {
        r=x/(x/l);
        ans-=1ll*(r-l+1)*query(x/l);
    }
    sumu[x]=ans;
    return ans;
}

ll calc(int n,int m)
{
    if (n>m) swap(n,m);
    //cout << n << " " << m << " ";
    int l,r;ll ans=0;
    for (l=1;l<=n;l=r+1)
    {
        r=min(n/(n/l),m/(m/l));
        ans+=1ll*(n/l)*(m/l)*(query(r)-query(l-1));
    }
    //cout << ans << endl;
    return ans;
}

ll solve(int n,int m,int k)
{
    if ((!n) || (!m)) return 0;
    node now=(node){n,m,k};
    if (ans[now]) return ans[now];
    if (k==1) ans[now]=calc(n,m);
    else
    {
        int i;
        for (i=0;i<len && fac[i]<=k;i++)
        {
            if ((k%fac[i]==0) && (mu[fac[i]]))
                ans[now]+=solve(m/fac[i],n,fac[i])*mu[fac[i]];
        }
    }
    return ans[now];
}

int main()
{
    n=read();m=read();k=read();
    sieve();
    rep(i,1,k)
        if (k%i==0) fac.pb(i);len=fac.size();
    printf("%lld",solve(n,m,k));
    return 0;
}

转载于:https://www.cnblogs.com/encodetalker/p/11129977.html

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值