筛法瞎写

看见APJifengc在写min25筛发现我这个东西都不会了。赶紧复习了一把去写点题。写一个更一个。

loj6682 梦中的数论

答案显然 \(\frac 12\sum_{i=1}^n\binom{d(i)}{2}=\frac 12\sum_{i=1}^nd^2(i)-d(i)\) 。推导略，以前好像写过。

对了根号的 \(\sum_{i=1}^nd(i)\) 是个傻逼问题，我拿这个诈骗joke3579然后他被骗了。

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
const int mod=998244353;
long long n,m,sq,cnt,sum,g[200010],w[200010];
int p[200010],pos1[200010],pos2[200010];
bool v[200010];
long long find(long long x){
    return x<sq?pos1[x]:pos2[n/x];
}
void get(int n){
    for(int i=2;i<=n;i++){
        if(!v[i])p[++cnt]=i;
        for(int j=1;j<=cnt&&i*p[j]<=n;j++){
            v[i*p[j]]=true;
            if(i%p[j]==0)break;
        }
    }
}
long long s(long long x,long long y){
    if(p[y]>=x)return 0;
    long long ret=(g[find(x)]-g[find(p[y])]+mod)%mod;
    for(int i=y+1;i<=cnt&&1ll*p[i]*p[i]<=x;i++){
        for(long long j=p[i],e=1;j<=x;j=j*p[i],e++){
            ret=(ret+1ll*(e+1)*(e+1)*(s(x/j,i)+(e>1)))%mod;
        }
    }
    return ret;
}
int main(){
    scanf("%lld",&n);
    sq=sqrt(n);
    get(sq);
    for(long long l=1,r;l<=n;l=r+1){
        r=n/(n/l);
        w[++m]=n/l;
        sum=(sum+1ll*w[m]*(r-l+1))%mod;
        if(w[m]<sq)pos1[w[m]]=m;
        else pos2[n/w[m]]=m;
        g[m]=4ll*(w[m]-1)%mod;
    }
    for(int j=1;j<=cnt;j++){
        for(int i=1;i<=m&&1ll*p[j]*p[j]<=w[i];i++){
            g[i]=(g[i]-1ll*(g[find(w[i]/p[j])]-g[find(p[j-1])])%mod+mod)%mod;
        }
    }
    long long ans=1ll*(s(n,0)+1-sum+mod)%mod*((mod+1)>>1)%mod;
    printf("%lld\n",ans);
    return 0;
}

DIVCNTK

我本机开氧气极限数据1.4s。

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
#define int unsigned long long
int n,m,sq,cnt,sum,k,g[2000010],w[2000010];
int p[2000010],pos1[2000010],pos2[2000010];
bool v[2000010];
int find(int x){
    return x<sq?pos1[x]:pos2[n/x];
}
void get(int n){
    for(int i=2;i<=n;i++){
        if(!v[i])p[++cnt]=i;
        for(int j=1;j<=cnt&&i*p[j]<=n;j++){
            v[i*p[j]]=true;
            if(i%p[j]==0)break;
        }
    }
}
int s(int x,int y){
    if(p[y]>=x)return 0;
    int ret=1ll*(k+1)*(g[find(x)]-g[find(p[y])]);
    for(int i=y+1;i<=cnt&&1ll*p[i]*p[i]<=x;i++){
        for(int j=p[i],e=1;j<=x;j=j*p[i],e++){
            ret=(ret+(1ll*k*e+1)*(s(x/j,i)+(e>1)));
        }
    }
    return ret;
}
signed main(){
    int tim;scanf("%llu",&tim);
    get(100000);
    while(tim--){
        scanf("%llu%llu",&n,&k);
        sq=sqrt(n);m=0;
        for(int l=1,r;l<=n;l=r+1){
            r=n/(n/l);
            w[++m]=n/l;
            sum=(sum+1ll*w[m]*(r-l+1));
            if(w[m]<sq)pos1[w[m]]=m;
            else pos2[n/w[m]]=m;
            g[m]=(w[m]-1);
        }
        for(int j=1;j<=cnt;j++){
            for(int i=1;i<=m&&1ll*p[j]*p[j]<=w[i];i++){
                g[i]=(g[i]-(g[find(w[i]/p[j])]-g[find(p[j-1])]));
            }
        }
        int ans=(s(n,0)+1);
        printf("%llu\n",ans);
    }
    
    return 0;
}

DIVCNT3

这题卡你小点，线筛就行了。

#pragma GCC optimize(3)
#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
#define int long long
int n,m,sq,cnt,sum,k,g[1000010],w[1000010];
int p[1000010],pos1[1000010],pos2[1000010];
int low[1000010],pw[1000010],f[1000010];
bool v[1000010];
int find(int x){
    return x<sq?pos1[x]:pos2[n/x];
}
void get(int n){
    f[1]=1;
    for(int i=2;i<=n;i++){
        if(!v[i])p[++cnt]=i,low[i]=i,pw[i]=1,f[i]=4;
        for(int j=1;j<=cnt&&i*p[j]<=n;j++){
            v[i*p[j]]=true;
            if(i%p[j]==0){
                low[i*p[j]]=low[i]*p[j];
                pw[i*p[j]]=pw[i]+1;
                f[i*p[j]]=f[i/low[i]]*(3ll*pw[i*p[j]]+1);
                break;
            }
            low[i*p[j]]=p[j];pw[i*p[j]]=1;
            f[i*p[j]]=f[i]*f[p[j]];
        }
    }
    for(int i=1;i<=n;i++)f[i]+=f[i-1];
}
int s(int x,int y){
    if(p[y]>=x)return 0;
    int ret=4ll*(g[find(x)]-g[find(p[y])]);
    for(int i=y+1;i<=cnt&&1ll*p[i]*p[i]<=x;i++){
        for(int j=p[i],e=1;j<=x;j=j*p[i],e++){
            ret=(ret+(3ll*e+1)*(s(x/j,i)+(e!=1)));
        }
    }
    return ret;
}
signed main(){
    int tim;scanf("%lld",&tim);
    get(500000);
    while(tim--){
        scanf("%lld",&n);
        if(n<=500000){
            printf("%lld\n",f[n]);
            continue;
        }
        sq=sqrt(n);m=0;
        for(int l=1,r;l<=n;l=r+1){
            r=n/(n/l);
            w[++m]=n/l;
            sum=(sum+1ll*w[m]*(r-l+1));
            if(w[m]<sq)pos1[w[m]]=m;
            else pos2[n/w[m]]=m;
            g[m]=(w[m]-1);
        }
        for(int j=1;j<=cnt;j++){
            for(int i=1;i<=m&&1ll*p[j]*p[j]<=w[i];i++){
                g[i]=(g[i]-(g[find(w[i]/p[j])]-g[find(p[j-1])]));
            }
        }
        int ans=(s(n,0)+1);
        printf("%lld\n",ans);
    }
    return 0;
}

DIVCNT2

今天的成就：爆切了三个水黑题。

卡常。

#pragma GCC optimize(3)
#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
#define int unsigned long long
int n,m,sq,cnt,sum,k,g[2000010],w[2000010];
int p[2000010],pos1[2000010],pos2[2000010];
bool v[2000010];
int find(int x){
    return x<sq?pos1[x]:pos2[n/x];
}
void get(int n){
    for(int i=2;i<=n;i++){
        if(!v[i])p[++cnt]=i;
        for(int j=1;j<=cnt&&i*p[j]<=n;j++){
            v[i*p[j]]=true;
            if(i%p[j]==0)break;
        }
    }
}
int s(int x,int y){
    if(p[y]>=x)return 0;
    int ret=3*(g[find(x)]-y);
    for(int i=y+1;i<=cnt&&p[i]*p[i]<=x;i++){
        for(int j=p[i],e=1;j<=x;j=j*p[i],e++){
            ret+=(2*e+1)*(s(x/j,i)+(e!=1));
        }
    }
    return ret;
}
signed main(){
    int tim;scanf("%llu",&tim);
    get(1000000);
    while(tim--){
        scanf("%llu",&n);
        sq=sqrt(n);m=0;
        for(int l=1,r;l<=n;l=r+1){
            r=n/(n/l);
            w[++m]=n/l;
            sum=(sum+1ll*w[m]*(r-l+1));
            if(w[m]<sq)pos1[w[m]]=m;
            else pos2[n/w[m]]=m;
            g[m]=(w[m]-1);
        }
        for(int j=1;j<=cnt;j++){
            for(int i=1;i<=m&&1ll*p[j]*p[j]<=w[i];i++){
                g[i]-=(g[find(w[i]/p[j])]-j+1);
            }
        }
        int ans=(s(n,0)+1);
        printf("%llu\n",ans);
    }
    return 0;
}

DIVCNT1

真的还算筛法吗。

首先你发现这是个傻逼式子：

\[\sum_{i=1}^n\lfloor\frac ni\rfloor \]

然后你发现这题数据范围 longlong 以内所以过不了。进行一步没有什么卵用的转化：

\[-\lfloor \sqrt n\rfloor^2+\sum_{d=1}^{\lfloor \sqrt n\rfloor}\lfloor \frac nd \rfloor \]

它还是根号的。因为这玩意的图象关于 \(y=x\) 对称所以成立。

我们用一堆横纵坐标都是整数的向量来拟合这个图象。为了完美拟合，需要 Stern-Brocot Tree ~~（SB树）~~。

开个单调栈来存一堆向量拟合这个东西。显然单调栈里的元素斜率单调递减。假如说我们当前拟合到 \((x,y)\) 。那么一直执行如下过程：

取出栈顶向量，不断加上它，每步累计贡献，直到下一步会走进区域（就是 \(y=\dfrac 1x\) 和坐标轴围成的区域）内。
不断弹栈，直到加上栈顶刚好走进区域内。那么我们现在就有了两个向量，一个 \(l\) 加上刚好能进去，一个 \(r\) 刚好进不去。
用 SB 树不断二分一个中间向量 \(mid\) 。如果 \(mid\) 能走进去，那么进栈，\(r=mid\) ，接着二分。反之如果 \(r\) 的斜率小于图象在 \(x+x_mid\) 的斜率，那么之后二分的所有向量都不能走进去，直接停。反之 \(l=mid\) 继续。

边界是 \(y\le \sqrt[3]{n}\) ，直接暴力。总复杂度 \(O(\sqrt[3]n\log n)\) 。

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
#define int __int128
int read(){
    int x=0;char ch=getchar();
    while(!isdigit(ch))ch=getchar();
    while(isdigit(ch))x=(x<<3)+(x<<1)+(ch^48),ch=getchar();
    return x;
}
void print(int x){
    if(x>9)print(x/10);
    putchar(x%10+'0');
}
struct node{
    int x,y;
    node operator+(const node &s)const{
        return {x+s.x,y+s.y};
    }
}stk[1000010],l,r,p;
int ans,n,top,sq,cb;
signed main(){
    int tim=read();
    while(tim--){
        n=read();ans=top=0;
        sq=sqrt(n);cb=cbrt(n);
        p={n/sq,sq+1};
        stk[++top]={1,0};stk[++top]={1,1};
        while(1){
            l=stk[top];top--;
            while(n<(p.x+l.x)*(p.y-l.y)){
                ans+=p.x*l.y+(l.x-1)*(l.y+1)/2;
                p.x+=l.x;p.y-=l.y;
            }
            if(p.y<=cb)break;
            r=stk[top];
            while(1){
                node mid=l+r;
                if(n<(p.x+mid.x)*(p.y-mid.y)){
                    r=mid;stk[++top]=r;
                }
                else if(n*r.x<=1ll*(p.x+mid.x)*(p.x+mid.x)*r.y)break;
                else l=mid;
            }
        }
        for(int i=1;i<p.y;i++)ans+=n/i;
        ans=ans*2-sq*sq;
        print(ans);putchar('\n');
    }
    return 0;
}

posted @ 2022-11-20 10:34 gtm1514 阅读(61) 评论(4) 编辑收藏举报

会员力量，点亮园子希望

刷新页面返回顶部

gtm1514