bzoj 2820 莫比乌斯反演

题意:给定N, M,求1<=x<=N, 1<=y<=M且gcd(x, y)为质数的(x, y)有多少

难道只有我一个蒟蒻觉得难点不在莫比乌斯反演而在优化的预处理么...

思路和容斥原理后的2301很像,不在详细解释只列公式

ans=sigma(1) (1<=i<=n,1<=j<=m,gcd(i,j)=p,p为不大于min(n,m)的质数)

ans=sigma(1)(1<=i<=floor(n/p),1<=j<=floor(m/p),gcd(i,j)=1,p为不大于min(n,m)的质数)

令 f(n)表示 gcd(x,y)=n 的数对数,F(n) 表示 n能整除 gcd(x,y) 的数对数

易知 F(x) = floor(n/x) * floor(m/x)

由莫比乌斯反演得

f(x)=sigma(miu(d/x)*F(d)) (x能整除d)

     =sigma(miu(d/x)* floor(n/d)* floor(m/d))(x能整除d且 d<=min(n,m))

所以

ans = f(1) =sigma(miu(d) *floor(n/pd) * floor(m/pd)) (1<=d<=min(n,m),p为不大于min(n,m)的质数)

因为我们需要枚举质数p和d,所以每个ans的时间复杂度为 sqrt(min(n,m) )*min(n,m) ,情理之中意料之外的TLE..

我们令T=pd

则 ans=sigma(miu(T/p)*floor(n/T)*floor(m/T) )(1<=T<=min(n,m) ,p为不大于min(n,m) 的质数且p能整除T) 

我们如果能预处理出sigma(miu(T/p)) 那么就很愉快了...

按照莫比乌斯反演枚举除数的套路,我们需要对sigma(miu(T/p)) 维护前缀和,然后计算的时候直接用

那么,如何预处理sigma(miu(T/p)) 呢?

我们在线性筛的时候,我们处理一下T就可以了

令 g(i)=sigma(miu(i/p) (p为不大于min(n,m) 的质数且p能整除i)

(一) 对于质数,g(i)=1 

假设我们在处理 i*prime[j] 

(二)prime[j] 不能整除 i

         (1)miu[i]=0

                  说明i至少有一个质因子的次数大于1

                  ①i有且只有一个质因子的次数恰好为2

                     我们假设p就是那个次数为2的唯一的质因子,除了枚举出p以外,其他的莫比乌斯函数值都为0,

                    所以 g[i*prime[j]] = miu[i*prime[j]/p]

                  ②除了①以外的其他所有情况

                    显然,g[i*prime[j]]=0

         (2)miu[i]<>0

                   无论枚举哪个质数,莫比乌斯函数值都是一样的。我们用r[x]表示x的质因子个数

                   则 g[i*prime[j]]=r[i*prime[j]]*(-1)^(r[i*prime[j]]-1)

                  因为 r[i*prime[j]]=r[i]+1 ,且 g[i]= r[i]*(-1)^(r[i]-1) 

                  所以 g[i*prime[j]]=(r[i]+1)*(-1)^r[i] 

                  通过瞪眼观察法,我们很容易得到性质:g[i*prime[j]] 与 r[i]异号且绝对值大1

(三)prime[j]整除i

           说明i*prime[j]至少有一个质因子的次数大于1

           只有枚举的质数为prime[j]时,莫比乌斯函数值为miu[i]其余为0

           所以此时 g[i*prime[j]]=miu[i]

对于(二)中的(1),我们如何判断i是否有且只有一个质因子的次数恰好为2呢?

我们用一个数组f[i]表示,如果i有且只有一个质因子的次数恰好为2则f[i]=i 否则 f[i]=1

那么 f[i]=1 时有如下三种情况:

  ①i有且只有一个质因子的次数大于为2

  ②i有且不只有一个质因子的次数大于1

  ③i所有质因子的次数都为1

第三种情况通过莫比乌斯函数值我们很容易辨别,前两种情况没用根本无需辨别

接下来我们解决在线性筛中f[i]的赋值:

(一) 质数的f=1

(二) prime[j] 不能整除i

            此时 f[i*prime[j]]=f[i]

(三)prime[j]能整除i

           (1)f[i]=1 且 miu[i]!=0 则 f[i*prime[j]]=prime[j]

           (2) 其余情况都为1

线性筛预处理的工作到此全部结束..最后维护一个g[i]的前缀和即可

var
    t,n,m,i,last,c      :longint;
    flag                :array[0..10000010] of boolean;
    miu,prime,f         :array[0..10000010] of longint;
    g                   :array[0..10000010] of int64;
    ans                 :int64;
function min(a,b:longint):longint;
begin
   if a<b then exit(a) else exit(b);
end;

procedure pre_do;
var
    i,j,t:longint;
begin
   miu[1]:=1;
   t:=0;
   for i:=2 to 10000000 do
   begin
      if not flag[i] then
      begin
         inc(t);
         prime[t]:=i;
         miu[i]:=-1;
         f[i]:=1;
         g[i]:=1;
      end;

      for j:=1 to t do
        if (i*prime[j]>10000000) then break else
        begin
           flag[i*prime[j]]:=true;
           if i mod prime[j]<>0 then
           begin
              miu[i*prime[j]]:=-miu[i];
              f[i*prime[j]]:=f[i];
              if miu[i]=0 then
              begin
                 if f[i]=1 then g[i*prime[j]]:=0 else g[i*prime[j]]:=miu[i*prime[j] div f[i]];
              end else
              begin
                 if g[i]>0 then g[i*prime[j]]:=-g[i]-1 else g[i*prime[j]]:=-g[i]+1;
              end;
           end else
           begin
              miu[i*prime[j]]:=0;
              if f[i]=1 then
                if miu[i]=0 then f[i*prime[j]]:=1
                  else f[i*prime[j]]:=prime[j]
                    else f[i*prime[j]]:=1;
              g[i*prime[j]]:=miu[i];
              break;
           end;
        end;
   end;
   for i:=2 to 10000000 do g[i]:=g[i-1]+g[i];
end;

begin
   pre_do;
   read(t);
   while (t>0) do
   begin
      dec(t);
      read(n,m);
      if n>m then
      begin
         c:=n;n:=m;m:=c;
      end;
      i:=2;
      ans:=0;
      while (i<=n) do
      begin
         last:=min(n div (n div i),m div (m div i));
         ans:=ans+(g[last]-g[i-1])*(n div i)*(m div i);
         i:=last+1;
      end;
      writeln(ans);
   end;
end.
——by Eirlys

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值