【BZOJ4652】循环之美(NOI2016)-杜教筛

本文详细解析了一道关于纯循环小数的算法题,利用杜教筛算法进行求解。通过对题目条件的数学抽象,将问题转化为求特定条件下两数的最大公约数等于1的计数问题,并给出了具体的实现代码。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

测试地址:循环之美
做法:本题需要用到杜教筛。
首先,思考一下纯循环小数这个条件有什么意义。注意到,如果一个数 w w 的循环节长度为x,那么我们用 kxw k x w 减去 w w 可以得到一个整数u,那么 w w 就可以表示成ukx1,注意到分数能表达成这种形式是满足题目条件的充要条件,那么条件就变成一个分数能不能表示成这种形式,也就是存不存在分母 w|kx1 w | k x − 1 x1 x ≥ 1 的情况。
注意到 w|kx1 w | k x − 1 ,也就意味着 kxmodw=1 k x mod w = 1 ,注意到这个东西可能成立当且仅当 gcd(w,k)=1 gcd ( w , k ) = 1 ,即 w,k w , k 互质,那么我们就把题目条件转化成了一个可以列式的条件了。
综上所述,这题要求的式子就是:
ans=ni=1mj=1[gcd(i,j)=1][gcd(j,k)=1] a n s = ∑ i = 1 n ∑ j = 1 m [ gcd ( i , j ) = 1 ] [ gcd ( j , k ) = 1 ]
接下来有几种推式子的方法,由于我太菜了,只会比较简单的不能通过全部测试点的方法(在BZOJ的总时限下勉强能过),推法如下:
常规把 [gcd(j,k)=1] [ gcd ( j , k ) = 1 ] 处理掉,得到:
ans=d|kμ(d)ni=1d|j[gcd(i,j)=1] a n s = ∑ d | k μ ( d ) ∑ i = 1 n ∑ d | j [ gcd ( i , j ) = 1 ]
把用 jd j d 替换 j j ,得:
ans=d|kμ(d)i=1nj=1md[gcd(i,j)=1][gcd(i,d)=1]
我们发现后面式子的形式和我们一开始列出的式子非常相似,那么我们令 ans=g(n,m,k) a n s = g ( n , m , k ) ,那么有:
g(n,m,k)=d|kμ(d)g(md,n,d) g ( n , m , k ) = ∑ d | k μ ( d ) g ( ⌊ m d ⌋ , n , d )
边界为 g(0,m,k)=g(n,0,k)=0 g ( 0 , m , k ) = g ( n , 0 , k ) = 0 ,而 g(n,m,1)=ni=1mj=1[gcd(i,j)=1] g ( n , m , 1 ) = ∑ i = 1 n ∑ j = 1 m [ gcd ( i , j ) = 1 ] 这个式子我们在其他题目里推了好多遍了,直接杜教筛处理即可,这样我们就可以递归计算 ans a n s 了。这个方法可以拿到很高的分数:92分,但不幸的是,由于这个方法不是很优或者本人太菜,这里给的代码并不能过最大的两个点,有待学习。
以下是本人代码:

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll hashsiz=7000003;
ll n,m,k,limit;
ll prime[1100010],mu[1100010],summu[1100010];
ll hashlist[7000010]={0},hashval[7000010],inhashlist=0;
bool vis[1100010]={0};

void calc_mu()
{
    mu[1]=1;
    prime[0]=0;
    for(ll i=2;i<=limit;i++)
    {
        if (!vis[i])
        {
            prime[++prime[0]]=i;
            mu[i]=-1;
        }
        for(ll j=1;j<=prime[0]&&i*prime[j]<=limit;j++)
        {
            vis[i*prime[j]]=1;
            if (i%prime[j]==0)
            {
                mu[i*prime[j]]=0;
                break;
            }
            mu[i*prime[j]]=-mu[i];
        }
    }
    summu[0]=0;
    for(int i=1;i<=limit;i++)
        summu[i]=summu[i-1]+mu[i];
}

ll hashfind(ll x)
{
    ll pos=x%hashsiz;
    while(hashlist[pos]&&hashlist[pos]!=x)
    {
        pos++;
        if (pos>=hashsiz) pos-=hashsiz;
    }
    if (hashlist[pos]==x) return pos;
    else return -1;
}

void hashinsert(ll x,ll val)
{
    if (inhashlist>=hashsiz) return;
    inhashlist++;
    ll pos=x%hashsiz;
    while(hashlist[pos]&&hashlist[pos]!=x)
    {
        pos++;
        if (pos>=hashsiz) pos-=hashsiz;
    }
    hashlist[pos]=x;
    hashval[pos]=val;
}

ll sum_mu(ll n)
{
    if (n<=limit) return summu[n];
    int pos=hashfind(n);
    if (pos!=-1) return hashval[pos];
    ll ans=1;
    for(ll i=n;i>=2;i=n/(n/i+1))
    {
        ll l=max(n/(n/i+1)+1,2ll),r=i;
        ans-=(r-l+1)*sum_mu(n/i);
    }
    hashinsert(n,ans);
    return ans;
}

ll sum(ll n,ll m)
{
    ll ans=0;
    for(ll i=min(n,m);i>=1;i=max(n/(n/i+1),m/(m/i+1)))
    {
        ll l=max(n/(n/i+1),m/(m/i+1))+1,r=i;
        ans+=(sum_mu(r)-sum_mu(l-1))*(n/i)*(m/i);
    }
    return ans;
}

ll calc(ll n,ll m,ll k)
{
    ll ans=0;
    if (!n||!m) return 0;
    if (k==1) return sum(n,m);
    for(ll i=1;i*i<=k;i++)
        if (k%i==0)
        {
            ans+=mu[i]*calc(m/i,n,i);
            if (i*i!=k) ans+=mu[k/i]*calc(m/(k/i),n,k/i);
        }
    return ans;
}

int main()
{
    scanf("%lld%lld%lld",&n,&m,&k);
    limit=1;
    while(limit*limit*limit<=max(n,m)) limit++;
    limit*=limit;
    limit=max(limit,k);
    calc_mu();
    printf("%lld",calc(n,m,k));

    return 0;
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值