c++ - 查找 LCM 求和的最有效算法

Question

问题:查找

n 的范围:1<= n <=

主要挑战是处理可能很大的查询 (Q)。 1 <= Q <=

到目前为止我使用的方法:

蛮力

while(Q--)
{
    int N;
    cin>>N;
    for(int i=1;i<=N;i++)
        ans += lcm(i,N)/i ;
}

复杂度:

在 中预处理和处理查询

首先，我建立了一个表，其中包含每个 N 的 euler totient 函数的值。这可以在 O(N) 中完成。

void sieve()
{
    // phi table holds euler totient function value
    // lp holds the lowest prime factor for a number
    // pr is a vector which contains prime numbers 
    phi[1]=1;
    for(int i=2;i<=MAX;i++)
    {
        if(lp[i]==0)
        {
            lp[i]=i;
            phi[i]=i-1;
            pr.push_back(i);
        }
        else
        {
            if(lp[i]==lp[i/lp[i]])
                phi[i] = phi[i/lp[i]]*lp[i];
            else phi[i] = phi[i/lp[i]]*(lp[i]-1);
        }
    for(int j=0;j<(int)pr.size()&&pr[j]<=lp[i]&&i*pr[j]<=MAX;j++)
        lp[i*pr[j]] = pr[j];
}

对于每个查询，分解 N 并将 d*phi[d] 添加到结果中。

for(int i=1;i*i<=n;i++)
{
   if(n%i==0)
   {
    // i is a factor
    sum += (n/i)*phi[n/i];
    if(i*i!=n)
        {
            // n/i is a factor too
            sum += i*phi[i];
        }
   }
}

这需要 O(sqrt(N)) 。

复杂度:O(Q*sqrt(N))

在 O(1) 中处理查询

在我上面描述的筛选方法中，我添加了一个部分来计算我们在 O(NLogN) 中需要的答案

for(int i=1;i<=MAX;++i)
{
    //MAX is 10^7
    for(int j=i;j<=MAX;j+=i)
    {
        ans[j] += i*phi[i];
    }
}

不幸的是，对于给定的约束和时间限制(1 秒)超时。

我认为这涉及到一些关于 N 的质因数分解的巧妙想法。我可以使用上面构建的 lp(最低素数)表对 O(LogN) 中的数字进行质因数分解，但我无法弄清楚如何使用因式分解得出答案。

Answer 1

您可以尝试以下算法:

lcm(i,n) / i  = i * n / i * gcd(i, n) = n / gcd(i, n)

现在应该找到数字总和 n/gcd(i, n)。

让 n = p1^i1 * p2^i2 * p3^j3 其中数字 p1, p2, ... pk 是质数。

项目数 n/gdc(i, n) 其中 gcd(i , n) == 1 是 phi[n] = n*( p1-1)*(p2-1)*...*(pk-1)/(p1*p2*...*pk)，所以加到总和 n*phi[n] 。

项目数 n/gdc(i, n) 其中 gcd(i , n) == p1 是 phi[n/p1] = ( n/p1)*(p1-1)*(p2-1)*...*(pk-1)/(p1*p2*...*pk)，所以加起来 n/p1*phi[n/p1].

项目数 n/gdc(i, n) 其中 gcd(i , n) == p1*p2 是 phi[n/(p1 *p2)] = (n/(p1*p2))*(p1-1)*(p2-1)*...*(pk-1)/(p1*p2*...*pk)，因此添加到总和 n/(p1*p2)*phi[n/(p1*p2)]。

现在答案是总和

n/(p1^j1*p2^j2*...*pk^jk)  phi[n/(p1^j1*p2^j2*...*pk^jk)]

总的来说

j1=0,...,i1  
j2=0,...,i2
....  
jk=0,...,ik

此总和中的项目总数为 i1*i2*...*ik，明显小于 O(n)。

要计算此总和，您可以使用具有自由参数初始数字、当前表示和初始表示的递归函数:

initial = {p1:i1, p2:i2, ... ,pn:in} 
current = {p1:i1, p2:i2, ... ,pn:in} 
visited = {}

int calc(n, initial, current, visited):
   if(current in visited):
      return 0
   visited add current 
   int sum = 0
   for pj in keys of current:
      if current[pj] == 0:
        continue
      current[pj]--
      sum += calc(n, initial, current)
      current[pj]++  

   mult1 = n 
   for pj in keys of current:
      mult1 /= pj^current[pj]

   mult2 = mult1
   for pj in keys of current:
      if initial[pj] == current[pj]:
         continue
      mult2 = mult2*(pj -1)/pj

   sum += milt1 * mult2
   return sum