I've been struggling with this for over a month now. The coding challenge on HR is similar to the official Project Euler 233, but instead of finding $f(N)=420$, your code may need to solve for $f(N)=12$. Also, instead of the result being the sum of the $n$'s, it is just the count of the $n$'s.

Following is the official challenge:

Let $f(n)$ be the number of points with integer coordinates that are on a circle passing through $(0,0)$, $(0,n)$, $(n,0)$ and $(n,n)$.

It can be shown that $f(10000)=36$.

Given two integers $N$ and $m$, what is the number of all positive integers $n \le N$ such that $f(n)=4m$?

So, I must be going about this the wrong way, and that is basically what the HR community told me, but I really do not see any other way to generate the correct solution.

To use one of the smallest examples I can think of, if $N=13$... then we know that:

$f(5) = 12$, $f(10) = 12$, and $f(13) = 12$, so the answer would be $3$.

I can see that $f(n) = 12$ when the number has one and only one prime of the form $4k+1$.

So my approach would generally involve producing all primes $4k+1 \le N$ in order to determine the answer. But the problem with that approach is that when $f(n)=12$, $N$ can be as large as $5 \times 10^{10}$. That is far too many primes to generate for the code to run in under a couple seconds.

I can only imagine that there must be some way to come up with the answer $3$ if you know that $N=13$ without generating all primes $\le 13$, but I am really clueless.

I guess my actual question would be that: "Is there a way to find the solution without knowing all primes $\le N$ or would you, in fact, need all primes $\le N$ to solve for $f(n)=12$?"

If the answer to that question is "Yes, there is a way to find the solution without knowing all primes $\le N$", then my next question would be "can you offer some clues and/or point me to online material I can research that could help me recognize this solution?"

Thanks in advance!

  • 1
    $\begingroup$ Please include the definition of $f(n)$ and other details necessary to understand the statement of your question. You're far more likely to get assistance if people don't have to do an internet search to even parse your query. $\endgroup$ – Greg Martin Jun 12 at 15:26
  • 1
    $\begingroup$ Welcome to MSE! Please use Mathjax. $\endgroup$ – Vineet Jun 12 at 15:28

I found a PDF, "COUNTING PRIMES IN RESIDUE CLASSES" , which shows what appears to be a fast way to count the number of primes congruent to $l$ modulo $k$ up to $x$.

So I guess I'm saying this is kind of an answer (but I still have questions)

The equations are quite complex and I have not yet translated them into layman programmer terms in order to produce a C (my language of choice) algorithm.

I think what I really need is a person who speaks both "math" and "code" to translate this math for me, and then I can implement the algorithm and use it to solve the overall challenge.

I think I might understand $P_2(x,y,l)$ but I am unsure about a few elements in the following sections:

From - 2. Proof of Theorem 1

We now explain the method we used to compute $π(x, k, l)$ for large values of $x$. From now on, we assume that $k$ is fixed and we write $π(x, l)$ instead of $π(x, k, l)$. Let $y$ be a real positive number and let $T(x, y, l)$ be the set of positive integers $n$ such that \begin{cases} n ≤ x,\\ n ≡ l \text{(mod k)},\\ p | n ⇒ p > y \end{cases} Assume that $y$ is such that $x^{1/3} ≤ y ≤ x^{1/2}$. Then each element $n$ of $T(x, y, l)$ has at most two (not necessarily distinct) prime factors. Thus we can split this set into three disjoint subsets $T_0(x, y, l)$, $T_1(x, y, l)$, and $T_2(x, y, l)$, according to the number of (not necessarily distinct) prime factors.

Then the function $P_2(x,y,l)$ is later reduced to:

From - 2.1 Computation of $P_2(x,y,l)$

$$P_2(x,y,l)=\sum_{y<p\le x^{1/2}}\pi(x/p,lp^{-1})-\sum_{y<p\le x^{1/2}}\pi(p-1,lp^{-1})$$

Following are some of what I'm confused about:

I don't understand the significance of $y$.

Could the summation iteration not just say $x^{1/3} < p ≤ x^{1/2}$?

I am not sure what the term $lp^{-1}$ does exactly.

I think that any number $n^{-1}$ is the same as $\frac1n$ which would make it $\frac{l}{p}$ but that is a fraction and $pi(x,l)$ and $P_2(x,l)$ expect $l$ to be a whole number (I assume), so I get pretty lost there.

My best code translation of these summation formulas (which is very incomplete atm) is:


I assume k = 4.

I'm not defining pi(x,l) yet, as it involves more terms than just P2(x,l).

pi(x,l) calls P2(x,l) and P2(x,l) calls pi(x,l) so I declare them both
and then I attempt to define P2(x,l) with my best interpretation but I get stuck.


typedef uint64_t u64;

// The maximum number I'll be working with, 10^11, 100 billion.
#define MAX_N 100000000000

u64 *primes // all primes up to sqrt(MAX_N)
  , *p_last // the last prime up to sqrt(MAX_N)
  , *p_end // pointer that is 1 greater than p_last, so p_curr should always be < p_end

void init(); // this initializes primes up to sqrt(MAX_N) via generic sieve.

u64 P2(u64 x,u64 l); // declaration needed by pi(x,l)

// place-holder for other terms needed by pi(x,l)

u64 pi(u64 x, u64 l); // declaration needed by P2(x,l)

u64 P2(u64 x, u64 l) { // definition
  u64 sum = 0
    , cbrtx = cbrt(x)
    , sqrtx = sqrt(x)
    , *p_curr = primes
  while(*p_curr <= cbrtx) ++p_curr; // get to the first prime > cbrtx

    Note: I am not sure what lp^-1 is in this case.
          I am also not positive that you are supposed to start at cbrtx but it says
          to assume y >= to cbrtx and it says that p > y, so I think it is right.
  while(*p_curr <= sqrtx) sum += pi(x/*p_curr, ???) - pi(*p_curr - 1, ???);

  return sum;


  And then pi(x,l) would call one of 2 other functions, pi_e or pi_b:
    if(x > *p_last) {
      return pi_e(x,l), use the equations (of course I have to correctly interpret
        P2 and all the other terms in the equations in order to do this)
    } else {
      pi_b(x,l), use brute iteration through primes already sieved

int main() {
  // should print the number of primes congruent to 1 mod 4 <= 10^11
  printf("%ld\n", pi(10^11,1)); 


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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