bzoj 2820 mobius反演

学了一晚上mobius，终于A了一道了。。。。

假设枚举到i，质数枚举到p(程序里的prime[j])，要更新A=i*p的信息。

1. p|i
    这时A的素数分解式中，p这一项的次数>=2。

  考虑g(A)的求和式：

  如果枚举的质数p'不等于p，A/p'就也会有p这一项，次数>=2，这时候miu(A/p')=0

    如果枚举的质数p'=p，A/p=i，这一项就是miu(i)

    因此g(A)=miu(i)
2. p!|i (即i%p!=0）

   这时候A比i多一个质因子p，对miu(i)分情况讨论。

  2.1 miu(i)==0 (即i有大于1次的项)

           这时A除去任何一个p'都会留下i的那个大于1次的项，除非是下面这一种非常特殊的情况：

    2.1.1 i的素数分解式中，大于1次的项只有一个，且这一项为2次。记这一项为p0。

         这时除去任何一个p'!=p0都会留下这一项，但是除去p0则会得到A/p0——这个数所有的项都是1次的。因此g(A)=miu(A/p0)

          2.1.2 i的素数分解式大于1次的项不止一个 或者 大于1次的项唯一，但次数高于2次。易见g(A)=0

     2.2 miu(i)!=0 (即i全是1次) 这个时候A的项也全是1次。设r(x)为x的质因子个数。

     则可以得到g(A)=r(A)*(-1)^(r(A)-1)。因为除以任何一个p'，miu(A/p')都是一样的。

     同理g(i)=r(i)*(-1)^(r(i)-1)，且有r(A)=r(i)+1。 利用r(A)=r(i)+1可以方便地得到：g(A)和g(i)异号，且绝对值比g(i)多1。

     亦即g(A)=(g(i)>0)?-1:1 -g(i) 

    

看情况2.1.1，我们有这么个遗留问题：

如果x的大于1次的项唯一，且这一项为2次，则令f(x)为这个项，否则f(x)=1。

事实上f(x)=1包含3种情况：

1. 大于1的项不唯一

2. 大于1次的项唯一但大于2次。

3. 全为1次

1和2利用现有的结果无法区分，但事实上不需要区分。3则可以用miu(x)判出来。

好，我们来对付f(x)，仍然是线性筛，变量意义同g(x)的讨论。

1. p|i

 A由i把最小因子p的次数加1得到，显然这一项的次数>=2。

1.1 f(i)!=1 

1.1.1 如果f(i)=p，那么A中p的次数就是3次了，f(A)=1。

1.1.2 如果f(i)!=p，那么A中大于1次的项就不唯一了，仍有f(A)=1

因此f(i)!=1必然有f(A)=1

1.2 i全为1次 即f(i)=1且miu(i)!=0 这时显然f(A)=p

1.3 i不全为1次 即f(i)=1且miu(i)=0 这时显然f(A)=1

2. p!|i

A比i多一个1次的质因数p，那么应有f(A)=f(i)

//By BLADEVIL
var
    mu, prime, mindiv, g, f            :array[0..10000010] of longint;
    gs                                :array[0..10000010] of int64;
    n, m, tt                        :longint;
    ans                                :int64;
    
procedure init;
var
    i, j, a                            :longint;
begin
    mu[1]:=1;
    for i:=2 to 10000000 do 
    begin
        if mindiv[i]=0 then 
        begin
            inc(prime[0]);
            prime[prime[0]]:=i;
            mindiv[i]:=i;
            mu[i]:=-1;
            f[i]:=1;
            g[i]:=1;
        end;
        for j:=1 to prime[0] do 
        begin
            if i*prime[j]>10000000 then break;
            a:=i*prime[j];
            mindiv[a]:=prime[j];
            if i mod prime[j]<>0 then 
            begin
                mu[a]:=-mu[i];
                f[a]:=f[i];
                if mu[i]=0 then 
                begin
                    if f[i]<>1 then g[a]:=mu[a div f[i]] else g[a]:=0;
                end else 
                begin
                    if g[i]>0 then g[a]:=-g[i]-1 else g[a]:=-g[i]+1;
                end;
            end else 
            begin
                mu[a]:=0;
                if f[i]=1 then 
                    if mu[i]=0 then f[a]:=1 else f[a]:=prime[j] else 
                    f[a]:=1;
                g[a]:=mu[i];
                break;
            end;
        end;
    end;
    for i:=2 to 10000000 do gs[i]:=gs[i-1]+g[i];
end;

procedure main;
var
    k, i                            :longint;
    t, t1, t2                          :longint;
    
begin
    read(tt);
    for k:=1 to tt do 
    begin
        read(n,m);
        if n<m then 
        begin
            t:=n; n:=m; m:=t;
        end;
        ans:=0;
        i:=2;
        while i<=m do 
        begin
            t1:=n div (n div i);
            t2:=m div (m div i);
            if t1<t2 then t:=t1 else t:=t2;
            ans:=ans+(gs[t]-gs[i-1])*(n div i)*(m div i);
            i:=t+1;
        end;
        writeln(ans);
    end;

end;

begin
    init;
    main;
end.