hdu6584

2019-08-05

字数统计: 1.1k字 | 阅读时长≈ 6分

next hexonext butterfly volantis yearn yilia shoka indigo apollo landscape cactus matery icarus fluid material

转移自老blog

hdu6584

题目本质上是找分子分母小于n的第k小的真分数,我们直接对答案二分即可，采取分数二分的方式，找到一个区间，它包含最多一个解

$$ \begin{aligned} &答案肯定是找到一个分子分母小于n的分数\frac{p}{q}他满足下面的特征\\ &(\sum_{i=1}^n\sum_{j=1}^n[gcd(i,j)=1][\frac{i}{j} \leq \frac{p}{q}]) = k\\ &我们对左式化简\\ &=\sum_{i=1}^n\sum_{j=1}^n[gcd(i,j)=1][i \leq \frac{p}{q}j]\\ &=\sum_{i=1}^{\lfloor \frac{p}{q}j\rfloor}\sum_{j=1}^n[gcd(i,j)=1]\\ &=\sum_{i=1}^{\lfloor \frac{p}{q}j\rfloor}\sum_{j=1}^ne(gcd(i,j))\\ &=\sum_{i=1}^{\lfloor \frac{p}{q}j\rfloor}\sum_{j=1}^n(u*1)(gcd(i,j))\\ &=\sum_{i=1}^{\lfloor \frac{p}{q}j\rfloor}\sum_{j=1}^n\sum_{d|gcd(i,j)}u(d)*1(d)\\ &=\sum_{i=1}^{\lfloor \frac{p}{q}j\rfloor}\sum_{j=1}^n\sum_{d|gcd(i,j)}u(d)\\ &=\sum_{i=1}^{\lfloor \frac{p}{q}j\rfloor}\sum_{j=1}^n\sum_{d|i,d|j}u(d)\\ &=\sum_{d=1}^{n}u(d)\sum_{i=1}^{\lfloor \frac{p}{q}j\rfloor}\sum_{j=1}^n[d|i,d|j]\\ &=\sum_{d=1}^{n}u(d)\sum_{xd=1}^{\lfloor \frac{p}{q}(yd)\rfloor}\sum_{yd=1}^n[d|(xd),d|(yd)]\\ &=\sum_{d=1}^{n}u(d)\sum_{x=1}^{\lfloor\frac{\lfloor \frac{p}{q}(yd)\rfloor}{d}\rfloor}\sum_{y=1}^{\lfloor\frac{n}{d}\rfloor}1\\ &=\sum_{d=1}^{n}u(d)\sum_{y=1}^{\lfloor\frac{n}{d}\rfloor}\lfloor\frac{\lfloor \frac{p}{q}(yd)\rfloor}{d}\rfloor\\ &=\sum_{d=1}^{n}u(d)\sum_{y=1}^{\lfloor\frac{n}{d}\rfloor}{\lfloor \frac{p}{q}y\rfloor}\\ \end{aligned} $$

这里是可以求出答案的,对d分块，右边的部分采用类欧几里得算法
我们一直往下二分，直到区间足够小，最后用 Stern-Brocot Tree 或法雷序列找出答案

#include<bits/stdc++.h>
using namespace std;

typedef long long ll;

/****  * 超级积性函数线性筛 *  ****/
typedef long long ll;
const ll maxn=2e6+10;
ll no_pri[maxn]={0,1,0},pri[maxn],low[maxn];
ll PHI[maxn],DDD[maxn],XDX[maxn],MUU[maxn],SIG[maxn];
void f_ini(){
    for(ll i=2;i<maxn;i++){
        if(!no_pri[i]) low[i]=pri[++pri[0]]=i;
        for(ll j=1;pri[j]*i<maxn;j++){
            no_pri[pri[j]*i]=1;
            if(i%pri[j]==0) {
                low[pri[j]*i]=low[i]*pri[j];
                break;
            }
            else low[pri[j]*i]=pri[j];
        }
    }

    DDD[1]=PHI[1]=MUU[1]=SIG[1]=1;// 改这里
    for(ll i=1;i<=pri[0];i++){
        for(ll mul=pri[i],ct=1;mul<maxn;mul*=pri[i],ct++){
            DDD[mul]=ct+1;// 改这里
            SIG[mul]=SIG[mul/pri[i]]+mul;// 改这里
            MUU[mul]=ct==1?-1:0;// 改这里
            PHI[mul]=mul/pri[i]*(pri[i]-1);// 改这里
        }
    }

    for(ll i=2;i<maxn;i++){
        for(ll j=1;pri[j]*i<maxn;j++){
            ll x=low[i*pri[j]], y=i*pri[j]/x;
            DDD[x*y]=DDD[x]*DDD[y];
            MUU[x*y]=MUU[x]*MUU[y];
            PHI[x*y]=PHI[x]*PHI[y];
            SIG[x*y]=SIG[x]*SIG[y];
            if(i%pri[j]==0) break;
        }
    }

    for(ll i=1;i<maxn;i++) {
        DDD[i]+=DDD[i-1];
        MUU[i]+=MUU[i-1];
        PHI[i]+=PHI[i-1];
        SIG[i]+=SIG[i-1];
        XDX[i]=(DDD[i]-DDD[i-1])*i+XDX[i-1];
    }
}


struct frac{
    ll x,y;
    frac(ll x_=0,ll y_=1){
        ll gcd=__gcd(x_,y_);
        x=x_/gcd;
        y=y_/gcd;
    }
    frac operator +(const frac&rhs){
        ll lcm=y/__gcd(y,rhs.y)*rhs.y;
        return frac(x*(lcm/y)+rhs.x*(lcm/rhs.y),lcm);
    }
    frac operator /(ll k){
        ll gcd=__gcd(k,x);
        return frac(x/gcd,y*(k/gcd));
    }
    bool operator <=(const frac&rhs){
        ll lcm=y/__gcd(y,rhs.y)*rhs.y;
        return x*(lcm/y)<=rhs.x*(lcm/rhs.y);
    }
};

// a>=0 b>=0 c>0 n>=0         -> O(lg(a,c))
void calfgh(ll a,ll b,ll c,ll n,ll&f,ll&g,ll&h){
    ll A=a/c,B=b/c,s0=n+1,s1=n*(n+1)/2,s2=n*(n+1)*(2*n+1)/6;
    f=s1*A+s0*B;
    g=s2*A+s1*B;
    h=s2*A*A+s0*B*B+2*s1*A*B-2*B*f-2*A*g;// 先减掉
    a%=c,b%=c;
    ll m=(a*n+b)/c;
    if(m!=0) {
        ll ff,gg,hh; calfgh(c,c-b-1,a,m-1,ff,gg,hh);
        f+=n*m-ff;
        g+=(n*m*(n+1)-hh-ff)/2;
        h+=n*m*m-2*gg-ff;
    }
    h+=2*B*f+2*A*g;//再加上
}


ll count(frac k,int n){
    ll ret=0;
    for(int i=1,ed;i<=n;i=ed+1){
        ed=n/(n/i);
        ll a[3]; calfgh(k.x,0,k.y,n/i,a[0],a[1],a[2]);
        ret+=1ll*(MUU[ed]-MUU[i-1])*a[0];
    }
    return ret;
}

int main(){
    f_ini();
    ll t,n,k;
    scanf("%lld",&t);
    while(t--){
        scanf("%lld%lld",&n,&k);
        frac l(0,1),r(1,1);// [l,r]
        for(int ijk=0;ijk<40;ijk++){
            frac mid=(l+r)/2;
            ll ct=count(mid,n);//[0,mid]
            if(ct>=k)r=mid;
            else l=mid;
        }
        //[l,r]
        frac L(0,1),R(1,0);
        while(true){
            frac mid(L.x+R.x,L.y+R.y);
            if(mid.x<=n&&mid.y<=n&&l<=mid&&mid<=r){
                printf("%lld/%lld\n",mid.x,mid.y);
                break;
            }
            if(!(l<=mid)){
                L=mid;
            }
            if(!(mid<=r)){
                R=mid;
            }
        }
    }
}

本文作者： fightinggg
本文链接： http://fightinggg.github.io/yilia/yilia/hdu6584.html
版权声明： 本博客所有文章除特别声明外，均采用 MIT 许可协议。转载请注明出处！