【XSY2523】神社闭店之日 莫比乌斯反演

题目大意

  给你\(a_1\ldots a_n,l,c\)每次给你\(x,y\),求有多少个序列满足:长度\(\leq l\),每个元素是\([1,c]\),循环右移\(a_j(x\leq j\leq y)\)次后和原序列相同。

  \(n,q\leq 100000,l,c\leq{10}^9,lcm(a_1,\ldots a_n)\leq{10}^13\)

题解

  显然只有右移\(\gcd(a_x,a_{x+1},\ldots,a_y)\)次后和原序列相同才满足条件。

  先求出\(s=\gcd(a_x,a_{x+1},\ldots,a_y)\)

  枚举长度\(i\),答案\(ans\)
\[ \sum_{i=1}^lc^{\gcd(i,s)} \]
  先进行简单的莫比乌斯反演:
\[ \begin{align} ans&=\sum_{i=1}^lc^{\gcd(i,s)}\\ &=\sum_{d|s}c^d\sum_{i=1}^l[\gcd(i,s)=d]\\ &=\sum_{d|s}c^d\sum_{i=1}^{\frac{l}{d}}[\gcd(i,\frac{s}{d})=1]\\ &=\sum_{d|s}c^d\sum_{i|\frac{s}{d}}\mu(i)\lfloor\frac{l}{id}\rfloor\\ \end{align} \]
  根据套路,我们要枚举\(j=id\)
\[ \begin{align} ans&=\sum_{j|s}\lfloor\frac{l}{j}\rfloor\sum_{i|j}\mu(i)c^\frac{j}{i}\\ \end{align} \]
  记\(f(i)=\sum_{j|i}\mu(j)c^\frac{i}{j}\),发现\(f(i)\)\(s\)无关,所以可以先把所有\(f(i)\)预处理出来,每次直接枚举\(s\)的因子计算。

  注意到这题很有一个性质:所有\(a_i\)的lcm\(\leq {10}^{13}\),那么\(s\)一定是lcm的因子。\(\leq {10}^{13}\)的数最多有\(10752\)个因子,可以先把这些因子求出来,做一波\(因子个数\text{因子个数}^2\)的DP,询问时查表。

  时间复杂度:\(因子个数O(\text{因子个数}^2+q\log n\log a)\)

代码

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cstdlib>
#include<ctime>
#include<utility>
#include<cmath>
#include<functional>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int,int> pii;
typedef pair<ll,ll> pll;
void sort(int &a,int &b)
{
    if(a>b)
        swap(a,b);
}
void open(const char *s)
{
#ifndef ONLINE_JUDGE
    char str[100];
    sprintf(str,"%s.in",s);
    freopen(str,"r",stdin);
    sprintf(str,"%s.out",s);
    freopen(str,"w",stdout);
#endif
}
int rd()
{
    int s=0,c;
    while((c=getchar())<'0'||c>'9');
    do
    {
        s=s*10+c-'0';
    }
    while((c=getchar())>='0'&&c<='9');
    return s;
}
void put(int x)
{
    if(!x)
    {
        putchar('0');
        return;
    }
    static int c[20];
    int t=0;
    while(x)
    {
        c[++t]=x%10;
        x/=10;
    }
    while(t)
        putchar(c[t--]+'0');
}
int upmin(int &a,int b)
{
    if(b<a)
    {
        a=b;
        return 1;
    }
    return 0;
}
int upmax(int &a,int b)
{
    if(b>a)
    {
        a=b;
        return 1;
    }
    return 0;
}
ll gcd(ll a,ll b)
{
    return b?gcd(b,a%b):a;
}
ll lcm(ll a,ll b)
{
    return a/gcd(a,b)*b;
}
const ll p=998244353;
ll fp(ll a,ll b)
{
    ll s=1;
    for(;b;b>>=1,a=a*a%p)
        if(b&1)
            s=s*a%p;
    return s;
}
ll d[100010];
namespace seg
{
    struct node
    {
        int l,r,ls,rs;
        ll s;
    };
    node a[200010];
    int rt;
    int cnt;
    void build(int &p,int l,int r)
    {
        p=++cnt;
        a[p].l=l;
        a[p].r=r;
        if(l==r)
        {
            a[p].s=d[l];
            return;
        }
        int mid=(l+r)>>1;
        build(a[p].ls,l,mid);
        build(a[p].rs,mid+1,r);
        a[p].s=gcd(a[a[p].ls].s,a[a[p].rs].s);
    }
    ll query(int p,int l,int r)
    {
        if(l<=a[p].l&&r>=a[p].r)
            return a[p].s;
        int mid=(a[p].l+a[p].r)>>1;
        ll s=0;
        if(l<=mid)
            s=gcd(s,query(a[p].ls,l,r));
        if(r>mid)
            s=gcd(s,query(a[p].rs,l,r));
        return s;
    }
}
struct hashset
{
    ll v[20010];
    int w[20010];
    int t[20010];
    int h[100010];
    int n;
    void insert(ll x,int y)
    {
        int hs=x%100007;
        n++;
        v[n]=x;
        w[n]=y;
        t[n]=h[hs];
        h[hs]=n;
    }
    int query(ll x)
    {
        int hs=x%100007;
        int i;
        for(i=h[hs];i;i=t[i])
            if(v[i]==x)
                return w[i];
        return 0;
    }
};
hashset h;
struct list
{
    int v[20000010];
    int t[20000010];
    int h[20010];
    int n;
    void add(int x,int y)
    {
        n++;
        v[n]=y;
        t[n]=h[x];
        h[x]=n;
    }
};
list li;
ll a[20010];
int m=0;
ll miu[20010];
ll f[20010];
ll l,c;
ll e=1;
ll g[20010];
ll gao(int x)
{
    ll ans=0;
    int i;
    for(i=1;i<=x;i++)
        if(a[x]%a[i]==0)
            ans=(ans+l/a[i]*f[i])%p;
    ans=(ans+p)%p;
    return ans;
}
ll f1[20010];
ll f2[20010];
ll ans[20010];
int main()
{
    open("b");
    int i;
    int n,q;
    scanf("%d%lld%lld%d",&n,&l,&c,&q);
    for(i=1;i<=n;i++)
    {
        scanf("%lld",&d[i]);
        e=lcm(e,d[i]);
    }
    for(i=1;ll(i)*i<=e;i++)
        if(e%i==0)
        {
            a[++m]=i;
            if(ll(i)*i!=e)
                a[++m]=e/i;
        }
    int j;
    sort(a+1,a+m+1);
    for(i=1;i<=m;i++)
    {
        g[i]=fp(c,a[i]);
        h.insert(a[i],i);
    }
    for(i=1;i<=m;i++)
        for(j=1;j<=i;j++)
            if(a[i]%a[j]==0)
                li.add(i,j);
    for(i=1;i<=m;i++)
        if(a[i]!=1)
        {
            if(a[i]<=1000000000ll&&e%(a[i]*a[i])==0)
                f1[h.query(a[i]*a[i])]=1;
            for(j=li.h[i];j;j=li.t[j])
                if(li.v[j]!=1&&li.v[j]!=i)
                {
                    f1[i]|=f1[li.v[j]];
                    f2[i]=1;
                }
        }
    for(i=1;i<=m;i++)
        if(f1[i])
            miu[i]=0;
        else
        {
            int s=0;
            for(j=li.h[i];j;j=li.t[j])
                if(!f2[li.v[j]])
                    s++;
            miu[i]=(s&1?1:-1);
        }
    for(i=1;i<=m;i++)
    {
        f[i]=0;
        for(j=li.h[i];j;j=li.t[j])
            f[i]=(f[i]+g[li.v[j]]*miu[h.query(a[i]/a[li.v[j]])])%p;
    }
    for(i=1;i<=m;i++)
        ans[i]=gao(i);
    seg::build(seg::rt,1,n);
    for(i=1;i<=q;i++)
    {
        int x,y;
        scanf("%d%d",&x,&y);
        ll s=seg::query(seg::rt,x,y);
        ll ss=ans[h.query(s)];
        printf("%lld\n",ss);
    }
    return 0;
}

转载于:https://www.cnblogs.com/ywwyww/p/8513362.html

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值