线性时间欧拉全函数计算

Linear time Euler's Totient Function calculation

本文关键字:函数 计算 时间 线性      更新时间:2023-10-16

经过一段时间的浏览,我找到了这个使用Eratostenes筛计算线性时间内欧拉φ值的代码。但我没有理解这个代码中使用的主要逻辑,特别是内部for循环中做了什么,以及这个循环中计算φ值的想法。如果有人能帮助我理解这段代码,那会很有帮助。

#define MAXN 3000000
int phi[MAXN + 1], prime[MAXN/10], sz;
bitset <MAXN + 1> mark;
for (int i = 2; i <= MAXN; i++ ){
   if(!mark[i]){
      phi[i] = i-1;
      prime[sz++]= i;
   }
    for (int j=0; j<sz && prime[j]*i <= MAXN; j++ ){
            mark[prime[j]*i]=1;
            if(i%prime[j]==0){
                  phi[i*prime[j]] = phi[i]*prime[j];
                  break;
            }
            else phi[i*prime[j]] = phi[i]*(prime[j]-1 );
   }
}

编码不在这里也不在那里,但算法本身就非常出色。不过,这与埃拉托斯梯尼之筛关系不大。这种方法有点让人想起Sundaram筛,因为它系统地产生素数的倍数,以标记复合物;它比Sundaram的要好,因为每个复合物都只划掉一次(没有透支)。同样地,每个phi值被精确地计算和分配一次。

如果先更改一下代码,就会更容易理解:

enum {  N = 3000000  };
vector<unsigned> primes;
vector<unsigned> phi(N + 1);        // indexed directly with numbers, hence the +1
vector<bool> is_composite(N + 1);   // indexed directly with numbers, hence the +1
phi[1] = 1;   // phi[0] is 0 already; phi[1] needs to be set explicitly
for (unsigned i = 2; i <= N; ++i)
{
   if (not is_composite[i])  // it's a prime
   {
      phi[i] = i - 1;        // a prime is coprime to all numbers before it
      primes.push_back(i);
   }
   for (auto current_prime: primes)
   {
      unsigned ith_multiple = i * current_prime;
      if (ith_multiple > N)
         break;
      is_composite[ith_multiple] = true;
      if (i % current_prime)  // i and current_prime are coprime -> phi is multiplicative in this case
      {
         phi[ith_multiple] = phi[i] * (current_prime - 1);  // (current_prime - 1) == phi[current_prime]
      }
      else
      {
         phi[ith_multiple] = phi[i] * current_prime;  // based on phi(p^(k+1)) / phi(p^k) == p
         break;
      }
   }
}

phi(k)值的计算必须考虑三种不同的情况:

(1) k是素数:phi(k)=k-1,这是平凡的

(2) k=m*p,其中m和p互素且p素:φ(k)=φ(m*p)=φ

(3) k=m*p,其中m=m'*p^n和m',p互素和p素:phi(k)=phi(m)*p

这两个相关的数学事实在维基百科关于欧拉总函数的文章中的欧拉乘积公式下列出。情况1在外循环中处理,情况2是内循环中条件的then分支,情况3是同样终止内循环的else分支。情况2和3根据较小数字的预先存在的phi值构建复合物的phi,所有这些最终都源于素数的phi(设置在外循环中)。

该算法的精彩之处在于它安排要做的工作的方式,即每个值只计算一次,并且在需要时已经计算过了。它为复合物实现的递归是基于通过分解出其最小素因子来有效地拆分每个复合物:m=m'*p。这有助于根据上述情况2和3计算复合物的phi,并为生成没有重复的复合物提供了一个简单的规则。除其他外,外环向内环呈现所有可能的m'(尽管对于i>N/2,内环从不取任何,外环旋转仅用于收集剩余素数)。然后内环产生所有的复合物,这些复合物是m'和一个不超过其自身素数中最小素数的素数的乘积。

该代码已根据从OEIS页面检索到的用于phi函数的前100000个phi值的列表进行了验证。它是为了展示算法的工作原理而写的;在程序中使用之前,必须对其进行审查、调整和加固(特别是为了防止溢出)。

附录-如果有人无意中跳过DanaJ的评论:is_composite[]数组是多余的,因为给定i的非复合性(素性)可以通过在外循环中测试phi[i]的零来确定。原因是复合材料的phi值向上传播——当i是它们的因子之一时,它们在早期迭代中得到计算。另一种推理方式是,只有在计算和存储相应的phi值时,is_composite[m]才会被设置为true,并且这些值永远不可能为零。因此,外环中的测试变为if (phi[i] == 0)。实施者(而不是"纯粹的"鉴赏家)可能希望更仔细地研究DanaJ的评论。。。