Euler Projesi 271. Soru

Modüler Küpler, kısım 1

Bir pozitif n sayısı için $1<x<n$ ve $x^3\equiv 1$ mod n olmak üzere x tam sayılarının toplamı S(n) ile tanımlansın.

n=91 iken x için 8 olası değer mevcut:  9, 16, 22, 29, 53, 74, 79, 81.
Yani, S(91)=9+16+22+29+53+74+79+81=363.

S(13082761331670030) kaçtır?
Cevap: 4617456485273129588 Python: #!/usr/bin/python import itertools def getx(a,pr): k = len(a) x = [0]*k for i in range(k): x[i] = a[i] for j in range(i): x[i] = r[i][j] * (x[i] - x[j]) % pr[i]; res = 0 mul = 1 for i in range(k): res += x[i] * mul mul *= pr[i] return res def find_obr(a,mod): return a ** (mod - 2) % mod pr = [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43] r = [[ find_obr(a,mod) for a in pr] for mod in pr] T = [[i for i in range(1,p) if (i**3) % p == 1] for p in pr] print sum([getx(a,pr) for a in apply(itertools.product,T)]) - 1 C#: int[] primes = new int[] { 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43 }; List<long> oknums = new List<long>(new long[]{121, 151}); long n = 2 * 3 * 5 * 7; for (int i = 4; i < primes.Length; i++) { oknums.Add(1); List<long> newoknums = new List<long>(); long newprime = primes[i]; long oldn = n; n *= newprime; for (int k = 0; k < oknums.Count; k++) { long num = oknums[k]; while (num < n) { if (num != 1) { //http://www.codeproject.com/KB/cs/biginteger.aspx BigInteger bi = new BigInteger(num); bi = bi.modPow(3, n); if (bi == 1 && !newoknums.Contains(num)) newoknums.Add(num); } num += oldn; } } oknums = newoknums; } for (int i = 0; i < oknums.Count; i++) { sum += oknums[i]; } Console.WriteLine(sum); C/C++: #include<stdio.h> #define int64 long long int // solves n*x==a (mod m), iterations int solve(int n, int a, int m) { int x=0,y=0; while(x<m){x++;y+=n;if(y>m)y-=m;if(y==a)return x;} fprintf(stderr,"solve: no solutions in %d x = %d (mod %d)\n",n,a,m); } int64 A=1ll*2*3*5*11*17*23*29*41; int p[]={7,13,19,31,37,43}; int m[]={2, 3, 7, 5,10, 6}; // solves m*m+m+1==0 (mod p) int64 x[6][43][43]; int main() { int i,j,ii,jj; int64 a,b,c,sum=0; for (j=0;j<6;j++)for(a=1;a<p[j];a++)for(b=1;b<p[j];b++) x[j][a][b]=solve(a,b,p[j]); for(i=1;i<3*3*3*3*3*3;i++){ a=A;b=1; for(ii=i,j=0;j<6;j++,ii/=3)if(!(ii%3))a*=p[j]; for(ii=i,j=0;j<6;j++,ii/=3)if(jj=ii%3){ c=m[j];if(jj>1)c=p[j]-1-c; c-=b; c%=p[j];if(c<0)c+=p[j]; b+=a*x[j][a%p[j]][c]; a*=p[j]; } printf("%lld\n",b); sum+=b; } printf("%lld\n",sum); } Haskell: -- 13082761331670030 -- = product.take 14 $ primes -- = product [2,3,5,7,11,13,17,19,23,29,31,37,41,43] import Number (primes) import Modulo (congruence) cbrt p = filter (\q -> mod (q^3) p == 1)[1..p-1] p270 ps = pred.sum $ map (congruence ps) (sequence.map cbrt $ ps) main = print.p270.take 14 $ primes -- congruence :: (Integral a) => [a] -> [a] -> a -- congruence ms as = (`mod` m).sum.zipWith3 add as vs $ us -- where m = product ms -- us = [div m mi | mi <- ms] -- vs = zipWith invMod us ms -- add a v u = a * v * u -- -- invMod :: Integral a => a -> a -> a -- invMod n = fst.gcd' n -- -- gcd' :: Integral a => a -> a -> (a, a) -- gcd' m n | m < n = let (a,b) = gcd' n m in (b,a) -- | r == 0 = (1,1-q) -- | otherwise = let (a,b) = gcd' n r -- in (b,a-b*q) -- where (q,r) = divMod m n Maple: pb271:=proc(nth,P,N) local s,i; if nth=15 then return(N) fi: s:=[msolve(x**3=1,ithprime(nth))]; return(add(pb271(nth+1,P*ithprime(nth),chrem([N,rhs(s[i][1])],[P,ithprime(nth)])),i=1..nops(s))) end: pb271(2,2,1)-1; Mathematica: Total@ReplaceAll[ List @@ Reduce[x^3 == 1, x, Modulus -> 13082761331670030], Equal[x, y_] -> y] - 1 Ruby: def extended_euclid(a, b) return [0,1] if a % b == 0 x,y = extended_euclid(b, a % b) return y, x - y * (a / b) end factors = [2,3,5,7,11,13,17,19,23,29,31,37,41,43] roots = [1] prod = 1 factors.each_with_index do |f,i| new_roots = [] roots.each do |root| (1...f).each do |x| if (x**3) % f == 1 s1,s2 = extended_euclid(prod, f) new_roots << (root*s2*f + x*s1*prod) % (prod*f) end end end prod *= f roots = new_roots end p roots.inject(:+)-1