1
$\begingroup$

I want to optimize an algorithm for calculating $g(n)=\sum^n_{k=2} p(\Omega(k))\Omega(k)$ where $$ p(n) = \begin{cases} 1 &\text{if $n$ is odd} \\ -1 &\text{if $n$ is even} \end{cases}$$ and $\Omega(k)$ is the prime factor counting function or more exactly, I want to optimize calculating all its values up to $n$

My current best try is

def Omega(k):
    Factor = 0
    b=1
    while(k>1):
        b += 1
        if (k % b == 0):
            while k % b == 0:
                k /= b
                Factor += 1
    return Factor
sum = 0
for i in range(n)
   w = Omega(i)
   sum += (w % 2 * 2 - 1) * w
   #using it

First, I am trying to calculate the current algorithm's computational complexity. we calculate $\Omega(k) ~~n$ times, I haven't managed to calculate $\Omega(k)$'s computational complexity but it seems to be $O(\sqrt k \log k)$ thus, I think the entire program is $O(n^{1.5} \log n)$.
Am I correct?
Can anyone help optimize my algorithm?

Ideas:

  1. Calculating $\Omega$ less. $\Omega$ is an Completely Additive function, i.e. for any $a,b$, $\Omega(ab)=\Omega(a)+\Omega(b)$. I've tried to use this to optimize the calculation, but have not managed to.
  2. Optimizing the prime factorizing algorithm. I'm unsure about this.
    The reason I'm interested in this is to graph the function. My current code is the first thing that went into my head, i.e., very nonoptimal.
    This function is related to Liouville function.
$\endgroup$
7
  • $\begingroup$ If you calculate Phi(k) for every k, that shows you haven't thought about this one bit. There is no need to find for example all integers with exactly two prime factors if you only need to know how many there are. $\endgroup$
    – gnasher729
    Jan 6, 2021 at 14:23
  • $\begingroup$ What programming competition is this coming from? $\endgroup$
    – gnasher729
    Jan 6, 2021 at 14:24
  • $\begingroup$ This isn't coming from a programming competition, I am just interested in this function. also, $\Omega(4)= 2 = 1+1 = \Omega(2)+\Omega(2)$, am I missing something? $\endgroup$
    – razivo
    Jan 6, 2021 at 14:28
  • $\begingroup$ I've edited in the reason I'm interested in this. $\endgroup$
    – razivo
    Jan 6, 2021 at 14:49
  • $\begingroup$ I'm calculating every $g(k)$ for every $k$ because this is what I need to graph it $\endgroup$
    – razivo
    Jan 6, 2021 at 15:38

1 Answer 1

2
$\begingroup$

Algorithm in the question

Although it is correct, it is rather slow. In function Omega(k), b goes from 1 to the largest prime divisor of $k$, which is of order $\dfrac {\pi^2n} {12\,\log n}$ on average for $1 \le k \le n$ as mentioned by the On-Line Encyclopedia of Integer Sequences. So, the running time of computing $g(n)$ is $O(\dfrac {n^2}{\log n})$. It takes more than an hour to compute $g(10^6)$ on my computer.

A Faster Algorithm

Instead of computing directly the number of prime factors of each number $k$, we can add 1 to $\Omega(k)$ whenever $k$ is a multiple of a power of $p$ for all primes $p$. For example, since 72 is a multiple of 2, 4, 8, 3 and 9 but no other powers of primes, we will get $$\Omega(72)=(1+1+1)+(1+1)=5.$$

Here is the improved algorithm.

from itertools import accumulate

def g_up_to(n):
    # return [g(0), g(1), g(2), ..., g(n)], where g(0)=g(1)=0.

    # Omega[k] for k >= 1 is the total number of prime factors
    # of k, counting multiplicity. Omega[0] = 0.
    Omega = [0] * (n+1)
    for k in range(2, n+1):
        if Omega[k] == 0:  # k is a prime.
            power = k
            while power <= n:
                m = power
                while m <= n:
                    Omega[m] += 1
                    m += power
                power *= k
    return list(accumulate((x % 2 * 2 - 1) * x for x in Omega))

print(g_up_to(9))  # print $g[0], ..., g[9]
# [0, 0, 1, 2, 0, 1, -1, 0, 3, 1]
print(g_up_to(10**6)[-5:])  # print g[999996], ..., g[1000000]
# [-14178, -14180, -14184, -14177, -14189]

On average, each $O(1)$ number of operations will yield $1$ prime factor. Thus, the running time of the algorithm above is $$O(1)\cdot O\left(\sum_{k=1}^n\Omega(k)\right)=O(n\log\log n),$$ where $\log\log n$ is the average of $\Omega(i)$ for $1\le i\le n$, as mentioned in Wikipedia. Since $\log\log n$ grows very slowly, this algorithm can be considered as an almost linear algorithm.

It takes less than 1 second to run g_up_to(10**6) on my computer.

$\endgroup$
1
  • $\begingroup$ Just determining all the primes no less than $n$ will need $O(n\log\log n)$ time. Furthermore, the hidden constant multipliers for both $O(n\log\log n)$ are within a factor of 2 to each other in all reasonable implementations. So my algorithm above should be optimal within a factor of 2, if we have to (implicitly) determine all the primes no less than $n$ during the algorithm. $\endgroup$
    – John L.
    Jan 27, 2021 at 16:45

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.