BZOJ3309 DZY Loves Math

BZOJ3309 DZY Loves Math

参考自这位大佬的推导思路以及这位大佬的实现思路

题目简述

定义\(f(n)\)为\(n\)所含质因子的最大幂指数，如\(f(18)=f(2*3^2)=2\)
输入T组数据，求\(\sum_{i=1}^{n} \sum_{j=1}^{m} f(gcd(i,j))\)的值
\(T<=10000,1<=a,b<=10^7\)

公式推导(略)

\[\sum_{i=1}^{n} \sum_{j=1}^{m} f(gcd(i,j))\]
\[=\sum_{T=1}^{min(n,m)} \lfloor \frac nT \rfloor \lfloor \frac mT \rfloor \sum_{d|T} f(d) \mu(\frac Td)\]

• 令\(g(T)=\sum_{d|T}f(d) \mu(\frac Td)\)
  则如今的目标是在线性时间内求出\(g(x)\)的值
  设\(T=p_1^{a_1}p_2^{a_2}…p_k^{a_k}\)，\(d=p_1^{b_1}p_2^{b_2}…p_k^{b_k}\)，\(c_i=a_i-b_i\)
  则\(\frac Td=p_1^{c_1}p_2^{c_2}…p_k^{c_k}\)
  要使\(\mu(\frac Td) \neq0\)，则\(c_1,c_2,…,c_k \leq1\)

分类讨论\(g(T)\)的值

1. \(a_1,a_2,…,a_k\)并不全都相等
   此时设\(a_x\)最大，对于全部包含\(p_x^{a_x}\)的\(d\)均有\(f(d)=a_x\)
   对于这些\(d\)，\(p_i(i \neq x)\)的指数的选择均有两个，为\(a_i-1,a_i\)，则一共能组合出\(2^{k-1}\)个不一样的\(d\)
   设其中\(d_1,d_2\)的质因子中除\(p_y\)一项以外的指数均相等，则\(\mu(d_1)=-\mu(d_2)\)
   而这样\(\mu\)值相反的\(d\)共有\(2^{k-2}\)组，各自相加后的和就为0即
   \[[\sum_{d|T,p_x^{a_x} |d} f(d) \mu(\frac Td) ]=0\]
   以后对于剩下取\(p_x^{a_x-1}\)的\(d\) 执行一样的分析，发现全部以\(a_q\)为最大值的的\(d\)的和均为0，因此此时\(g(T)=0\)

2. 令\(A=a_1=a_2=…=a_k\)即\(a_1,a_2,…,a_k\)均相等
   此时\(f(d)=\begin {cases} A(b_1,b_2,…,b_k中至少有一个值为A)\\ A-1(b_1=b_2=…=b_k=A-1) \end {cases}\)
   假设\(b_1=b_2=…=b_k=A-1\)的时\(d\) 值\(f(d)\)也为A
   则\(g(T)\)相似于上一种状况的推导，其的和也为0
   固然，最后\(g(T)\)要减去某特殊\(d\)多加的\(1*\mu(\frac Td)\)
   因此此状况下\(g(T)=-(-1)^{k}=(-1)^{k+1}\)

而线性求g(T),则须要记录每一个点的最小质因子及个数，除去该最小质因子后的值
这里就感性理解一下吧OvO

#include <iostream>
#include <cstdio>
#include <cstring>
#define N 10000010
using namespace std;
int n,prim[N],pcnt,at[N],la[N],cnt[N];//at[i]记录i的最小质因子，la[i]=记录i除以全部at[i]后的值，cnt[i]记录at[i]的数量 
long long sum[N],g[N],f[N];
void eular()
{
    g[1]=0;
    for(int i=2;i<=N-10;i++)
    {
        if(!at[i])
        {
            prim[++pcnt]=i;
            at[i]=pcnt;
            g[i]=la[i]=cnt[i]=1;
        }
        for(int j=1;j<=at[i];j++)
        {
            if((long long)i*prim[j]>N-10)break;
            at[i*prim[j]]=j;
            if(j==at[i])//prim[j]为i的最小质因子，即也为prim[j]*i的最小质因子 
            {
                la[i*prim[j]]=la[i];
                cnt[i*prim[j]]=cnt[i]+1;
                if(la[i]==1)g[i*prim[j]]=1;
                else g[i*prim[j]]=(cnt[i]+1==cnt[la[i]]?-g[la[i]]:0);
            }
            else//prim[j]为新出现的质因子，个数为1 
            {
                la[i*prim[j]]=i;
                cnt[i*prim[j]]=1;
                g[i*prim[j]]=(cnt[i]==1?-g[i]:0);
            }
        }
        sum[i]=sum[i-1]+g[i];
    }
}
long long mb(int n,int m)
{
    long long ans=0;
    if(n>m)swap(n,m);
    for(int l=1,r;l<=n;l=r+1)
    {
        r=min(n/(n/l),m/(m/l));
        ans+=(sum[r]-sum[l-1])*(n/l)*(m/l);
    }
    return ans;
}
int main(){
    eular();
    int a,b,T;
    cin>>T;
    while(T--)
    {
        scanf("%d%d",&a,&b);
        printf("%lld\n",mb(a,b));   
    }
    return 0;
}

\(\color{white}{神秘数字：4372wyx}\)