6 Ocak 2018

Euler Projesi 233. Soru

Bir çemberde kafes noktaları

F(N), (0,0), (N, 0), (0, N) ve (N, N) noktalarından geçen bir daire üzerinde olan tam sayı koordinatlarına sahip noktaların sayısı olsun.

F (10000) = 36 olduğu gösterilebilir.

F(N) = 420 olan tüm pozitif N ≤ 1011 tam sayılarının toplamı nedir?

Cevap: 271204031455541309

C++:
N=400000;
A=vector(N); 
 
highestpower(i)=if(i%2==1,0,1+highestpower(i/2))

for(i=1,N,\
ifact=factor(i);\
L=matsize(ifact);\
L=L[1];\
isok=1;\
for(j=1,L,if(ifact[j,1]%4==1,isok=0;));\
if(isok==1,j=i;while(j<=N,A[j]+=i;j+=2^highestpower(j);)))

B=vector(N);

for(i=1,N,\
j=i;\
while(j>0,B[i]+=A[j];j-=2^highestpower(j);))

count=0;

n=100000000000;
for(a=2,70,\
for(b=2,70,\
if(a*b==105,\
x=(a-1)/2;y=(b-1)/2;\
p1=5;\
while(p1^(x+y)<=n,\
if(isprime(p1),\
p2=p1+4;\
while(p1^x*p2^(y)<=n,\
if(isprime(p2),\
number=p1^x*p2^y;\
count+=number*B[n\number];);\
p2+=4;));\
p1+=4;))))

for(a=2,10,\
for(b=2,10,\
for(c=2,10,\
if(a*b*c==105,\
x=(a-1)/2;y=(b-1)/2;z=(c-1)/2;\
p1=5;\
while(p1^(x+y+z)<=n,\
if(isprime(p1),\
p2=p1+4;\
while(p1^x*p2^(y+z)<=n,\
if(isprime(p2),\
p3=p2+4;\
while(p1^x*p2^y*p3^z<=n,\
if(isprime(p3),\
number=p1^x*p2^y*p3^z;\
count+=number*B[n\number];);\
p3+=4;\
));\
p2+=4;));\
p1+=4;)))))
count

Python:
import math

nmax = 5000000

primes = nmax*[True]

primes[0:2] = [False,False]

for x in range(2, int(math.sqrt(len(primes)) + 1)):
    if primes[x]:
        for y in range(2 * x, len(primes), x):
            primes[y] = False

p = [x for x in range(len(primes)) if primes[x] and x % 4 == 1]

q = [x for x in range(len(primes)) if primes[x] and x % 4 != 1]

n = 10 ** 11

s = 0

def cnt(x, i = 0):
    global s
    s += x
    for j in range(i, len(q)):
        if x * q[j] > n:
            break
        cnt(x * q[j], j)

for a in p:
    if a ** 3 > n:
        break
    for b in p:
        if b != a:
            if a ** 3 * b ** 2 > n:
                break
            for c in p:
                if c != a and c != b:
                    if a ** 3 * b ** 2 * c > n:
                        break
                    cnt(a ** 3 * b ** 2 * c)

for a in p:
    if a ** 7 > n:
        break
    for b in p:
        if b != a:
            if a ** 7 * b ** 3 > n:
                break
            cnt(a ** 7 * b ** 3)

for a in p:
    if a ** 10 > n:
        break
    for b in p:
        if b != a:
            if a ** 10 * b ** 2 > n:
                break
            cnt(a ** 10 * b ** 2)

print(s)

Haskell:
module Main where

import Utils

import Data.List

main :: IO ()
main = print . foldl' (+) 0 $ solve (10^11)

solve :: Integer -> [Integer]
solve maxn = concatMap (multiples maxn) (baseNums maxn)

coolPrimes :: [Integer]
coolPrimes = filter ((== 1) . (`mod` 4)) primes

uncoolPrimes :: [Integer]
uncoolPrimes = filter ((/= 1) . (`mod` 4)) primes

uncoolNums :: Integer -> [Integer]
uncoolNums maxd = go 1 uncoolPrimes
   where
   go d pps@(p:ps)
      | d  > maxd = []
      | dp > maxd = [d]
      | otherwise = go dp pps ++ go d ps
      where
      dp = d * p

expPatterns :: [[Integer]]
expPatterns = concatMap permutations [[1,2,3],[3,7],[2,10]]

baseNums :: Integer -> [Integer]
baseNums maxn = go 1 coolPrimes =<< expPatterns
   where
   go n _ _
      | n > maxn    = []
   go n _ []        = [n]
   go n pps@(p:ps) ees@(e:es)
      | minn > maxn = []
      | otherwise   = go (n * p^e) ps es ++ go n ps ees
      where
      minn = n * product (zipWith (^) pps ees)

multiples :: Integer -> Integer -> [Integer]
multiples maxn n = map (* n) (uncoolNums (maxn `div` n))