Programming Praxis

6 Responses to “ABC Conjecture”

1. Paul said

   September 18, 2012 at 5:38 PM

   A python version. I left out code for gcd and prime factors

   def rad(n):
       return reduce(operator.mul, set(ma.prime.factors3(n)))
   
   def q(a, b, c):
       return math.log(c) / math.log(rad(a * b * c))
   
   def find_abc(limit=1000):
       pfac = [0] + [set(ma.prime.factors3(n)) for n in range(1, limit)]
       num = 0
       for a in range(1, limit - 2):
           for b in range(a + 1, limit - a):
               c = a + b
               if c < limit and not pfac[a] & pfac[b]:
                   qabc = q(a, b, c)
                   if qabc > 1:
                       num += 1
                       print "{0:3d}  {1:5.3f}  {2:3d}  {3:3d}   {4:3d}".format(num, qabc, a, b, c)
   find_abc(1000)
   
2. DGel said

   September 19, 2012 at 11:49 AM

   Haskell version of the suggested solution. Uses primes and factorization algorithm taken from the haskellwiki

   import Data.Function (fix)
   import Data.List (nub, sortBy, sort)
    
   primesTMEf = 2 : g (fix g) 
     where 
       g xs = 3 : gaps 5 (joinT [[x*x, x*x+2*x..] | x <- xs]) 
       joinT ((x:xs):t) = x : union xs (joinT (pairs t))  
       pairs ((x:xs):ys:t) = (x : union xs ys) : pairs t 
       gaps k s@(x:xs) | k<x  = k:gaps (k+2) s    -- equivalent to
                       | True =   gaps (k+2) xs   --  [k,k+2..]`minus`s, k<=head s
   
   union (x:xs) (y:ys) = case (compare x y) of 
              LT -> x : union  xs  (y:ys)
              EQ -> x : union  xs     ys 
              GT -> y : union (x:xs)  ys
   union  xs     []    = xs
   union  []     ys    = ys
   
   primeFactors n = factor primesTMEf n
     where 
       factor ps@(p:pt) n | p*p > n      = [n]               
                          | rem n p == 0 = p : factor ps (quot n p) 
                          | otherwise    =     factor pt n
   
   rad = product . nub . primeFactors
   q a b c = log (fromIntegral c) / log (fromIntegral $ rad (a * b * c))
   
   triples :: Int -> [(Int, Int, Int)]
   triples limit = triples' [(2,1),(3,1)]
       where
           triples' [] = []
           triples' ((a,b):xs)
               | c <= limit = (a,b,c) : triples' ([(2 * a - b, a), (2 * a + b, a), (a + 2 * b, b)] ++ xs)
               | otherwise = triples' xs
               where 
                   c = a + b
   
   main = putStr $ unlines . map show . sort $ [(qval, a, b, c) | (a, b, c) <- triples 1000, let qval = q a b c, qval > 1]
   
3. A piece of the abc conjecture | jverkamp.com said

   September 19, 2012 at 1:06 PM

   […] make a long story shorter, there was a challenge on Programming Praxis that intrigued me, which was to write code that given a upper bound on c would generate a list of […]
4. Catalin Cristu (@catalin_c) said

   September 19, 2012 at 2:14 PM

   A Java solution:

   import java.util.*;
   
   public class ABCConjecture {
   
       public static void main(String[] args) {
           List<ABCTriple> triples = abc(1000);
           Collections.sort(triples, new Comparator<ABCTriple>() {
               public int compare(ABCTriple t1, ABCTriple t2) {
                   if (t1.q < t2.q) return 1;
                   return -1;
               }
           });
           for (ABCTriple t : triples) {
               System.out.println(t);
           }
       }
   
       public static List<ABCTriple> abc(final int n) {
           List<ABCTriple> result = new LinkedList<ABCTriple>();
           for (int a = 1; a <= n - 1; a++) {
               for (int b = a + 1; b <= n; b++) {
                   int c = a + b;
                   if (c < n && coPrime(a, b, c)) {
                       double q = q(a, b, c);
                       if (q > 1.0) {
                           result.add(new ABCTriple(a, b, c, q));
                       }
                   }
               }
           }
           return result;
       }
   
       public static double q(int a, int b, int c) {
           return Math.log(c) / Math.log(rad(a * b * c));
       }
   
       public static int rad(int n) {
           return product(unique(primeFactors(n)));
       }
   
       private static boolean coPrime(Integer... numbers) {
           for (int i = 0; i < numbers.length - 1; i++) {
               for (int j = i + 1; j < numbers.length; j++) {
                   if (gcd(numbers[i], numbers[j]) != 1) return false;
               }
           }
           return true;
       }
   
       private static int gcd(int a, int b) {
           int max = Math.max(a, b);
           int min = Math.min(a, b);
           while (min != 0) {
               int r = max % min;
               max = min;
               min = r;
           }
           return max;
       }
   
       private static int product(Collection<Integer> numbers) {
           int result = 1;
           for (int i : numbers) {
               result *= i;
           }
           return result;
       }
   
       private static Set<Integer> unique(Collection<Integer> numbers) {
           return new HashSet<Integer>(numbers);
       }
   
       private static List<Integer> primeFactors(int n) {
           List<Integer> primeFactors = new LinkedList<Integer>();
           for (int p : PRIMES) {
               if (p * p > n) break;
               while (n % p == 0) {
                   primeFactors.add(p);
                   n /= p;
               }
           }
           if (n > 1) primeFactors.add(n);
           return primeFactors;
       }
   
       private static List<Integer> PRIMES = sieve(1000);
   
       private static List<Integer> sieve(int n) {
           List<Integer> primes = new LinkedList<Integer>();
           BitSet sieve = new BitSet(n + 1);
           for (int i = 2; i * i <= n; i = sieve.nextClearBit(i + 1)) {
               primes.add(i);
               for (int j = i + i; j <= n; j += i) {
                   sieve.set(j);
               }
           }
           return primes;
       }
   
   }
   
   class ABCTriple {
       final int a;
       final int b;
       final int c;
       final double q;
   
       ABCTriple(int a, int b, int c, double q) {
           this.a = a;
           this.b = b;
           this.c = c;
           this.q = q;
       }
   
       public String toString() {
           return q + " [" + a + ", " + b + ", " + c + "]";
       }
   }
   
5. JP said

   September 21, 2012 at 5:02 AM

   I went with quick and easy, using Racket’s excellent for*/list macro to generate the final step. Interestingly, although it just brute forces the problem rather than the “smarter” solution posted here, it actually runs about 20% faster on the n=1000 case than the code posted here (on the same machine and the same environment). Not entirely sure why, although I am curious if anyone has any thoughts.

   In any case, you can see my version on my blog.
6. Tadewos Abiye said

   September 29, 2012 at 6:57 AM

   This algorithm is still to slow to solve http://projecteuler.net/problem=127