So this means we are looking for the largest prime q that divides 4k^2 + 1, here is where I upped my level on sieve understanding.

I found the following article: http://www.devalco.de/quadr_Sieb_4x%5E2+1.php which explains how to prime sieve 4n^2 + 1. You can read the article for the mathematical proof, I will explain the algorithm.

1. Initialise a list f where f[x] = 4x^2 + 1 from x = 0 to x = limit (limit is defined by you)

2. Initialise a second list max_prime_factor = [0] * (limit + 1), in this list max_prime_factor[x] = the maximum prime factors of 4x^2 + 1

3. Go through the list f from the first element to the last element

   1. Initialise divisor = f[x]

   2. If f[x] > 1, then we continue with the sieve, otherwise we continue going through the loop

      1. From the article we need to check if, f[x + k*divisor] is divisible by divisor where (x + k*divisor) < limit

         1. While f[x + k*divisor] is divisible by divisor, we divide it and we need to assign max_prime_factor[x + k*divisor] = the maximum between it's current value and divisor, in python it looks like this: max_prime_factor[x + k*divisor] = max(max_prime_factor[x + k*divisor], fx)

      2. From the article we also need to check if, f[-x + k*divisor] is divisible by divisor where (-x + k*divisor) < limit

         1. While f[-x + k*divisor] is divisible by divisor, we divide it and we need to assign max_prime_factor[-x + k*divisor] = the maximum between it's current value and divisor, in python it looks like this: max_prime_factor[-x + k*divisor] = max(max_prime_factor[-x + k*divisor], fx)

4. At the end the list f should be reduce to a list of only 1's, and max_prime_factor is a list where max_prime_factor[x] = the maximum prime factors of 4x^2 + 1, so all we need to do is sum the list max_prime_factor and mod 10^18 to find the last 18 digits