Showing posts with label modified collatz conjecture. Show all posts
Showing posts with label modified collatz conjecture. Show all posts

Sunday, August 10, 2014

Modified Collatz Conjecture

I was always interested in mathematics, and have been tinkering around with puzzles, conjectures on and off;mostly when I saw a pen and paper lying around. Recently, I was trying out problems from Project Euler.
After poking around a few problems, I came across this - the modified collatz conjecture problem.

Collatz conjecture itself has been around for a long time and hasn't been proved/disproved yet.
The gist of the problem being - it is a sequence of steps which when executed for any given number always results in '1' after which its trapped in a infinite sequence around 1. And needless to stay, in my free time I did try my hand around proving it works for everything,

Modified Collatz Conjecture(MCC) itself follows a slightly different set of steps. The project euler problem revolves around finding the smallest number that follows a given pattern, the pattern being the sequence of steps.
I grabbed the nearest envelope lying around and started chipping away at the solution.

Solution
The brute force solution was checking if every number > N followed the sequence, by performing the steps from the Collatz Conjecture. This obviously did not scale after a point.
After a couple of hours, I had a reasonable working solution. I obtained a generic equation for the MCC. I could write an expression of the form

(Ax+B)/C

where ,
A - 2 ^ d, where 42 are from the Collatz expression (4+ 1)/3 and (2x -1)/3, and and are the total 'U' and 'd' steps in the sequence

- This is a complex expression, but essentially obtained by convolution of the constants of the conjecture operation (ax+b)/with increasing powers of 3

C - k, where is the sequence length

The solution was to solve for (Ax+B)/C = Y , where is what we would obtain after performing steps in the sequence. The smallest integer solution for xsuch that x > N is what I wanted. Substituting values for Y to get potential x worked fine for smaller solutions. But as I tried out the question's test pattern, I faced two issues

  1. My Java variables were overflowing and my computations were incorrect
  2. It took a long time to compute the correct Y.
Working with BigInteger
For the first problem, I had to move to a datatype that could handle bigger numbers. BigInteger seemed like the right datatype for handling the scale computations. It was a little cumbersome to use a function for every simple arithmetic operation when i started. But once I got used to it, I realised there were a lot of neat math functions that were wrapped into the datatype.
For example something like a % b could be accomplished by a.mod(b).

Diophantine equation and Euclid's extended algorithm
I noticed that brute force to find out the smallest integer was taking too long.
Enter Euclid's extended algorithm. It computed the GCD (F, G) for an equation of the form
+ G y = c  (Diophantine equation) and got generic expressions for and y , i.e
T,
y T,
where p, q, m, n are constants that are calculated and is any integer. After finding this equation, finding the right (such that N) by substituting the lower bound N as x, was easy.

The running time for the test case reduced from over a minute to 6ms while returning the correct solution.


The Java class and functions are given below.

The important functions
calcCollatzExpressionResult
calculateResultThroughGCDComputation
You can pass the test number to the below function. Initialize the sequence in the constructor.
computeSmallestNumWithSequence 

In the following method, you might have to replace the expression coefficients if you are using only the below code.
storeSequenceConstants





 package com.eulerproject.modifiedcollatz;  
   
 import java.math.BigInteger;  
   
 /**  
  * Solution Project Euler Problem 277 Modified Collatz Conjecture  
  * The problem is to calculate a number above a certain bound that follows a known sequence in the Modified  
  * Collatz Conjecture  
  *   
  * The solution involves creating a generic expression for k levels of the sequence and calculating the coefficients  
  * based on this general formula to get a equation of the form aX + bY = C, which we solve using euclid's extended   
  * algorithm   
  *   
  * The generic expression looks like (Ax + B)/C ,   
  * where A is the numerator sum, B is the constant sum and C is the denominator result  
  *   
  * @author Suhas Shetty  
  *  
  */  
 public class ModifiedCollatzSequence {  
   
      //Stores the coefficients for the constant term of the generic expression  
      private Integer[] constantCoeffArray;  
   
      //Stores the coefficients for the numerator multiplier term of the generic expression  
      private Integer[] numeratorCoeffArray;  
   
      private int index = 0;  
      private int upwardStepCount = 0;  
      private int downwardStepCount = 0;  
   
      // These are the values that will be calculated from the generic expression  
      private Double expressionNumeratorSum = 0d;  
      private Double expressionConstantSum = 0d;  
      private Double expressionDenominator = 0d;  
   
      private CollatzExpression upwardExpression;  
      private CollatzExpression downwardExpression;  
      private CollatzExpression bigDownwardExpression;  
   
      private String sequence = null;  
   
      /**  
       * Constructor for the class  
       *   
       * @param sequence  
       */  
      public ModifiedCollatzSequence(String sequence) {  
           this.sequence = sequence;  
           constantCoeffArray = new Integer[sequence.length()];  
           numeratorCoeffArray = new Integer[sequence.length()];  
   
           //Initializing the collatz expression coefficients  
           bigDownwardExpression = new CollatzExpression(1, 0, 3);  
           upwardExpression = new CollatzExpression(4, 2, 3);  
           downwardExpression = new CollatzExpression(2, -1, 3);  
      }  
   
      /**  
       * This method stores collatz expression coeffs into an array. These values come from the expression  
       * for the modified Collatz Conjecture. i.e  
       * We treat the Collatz operation to be performed as (ax + b)/c, where  
       * a - multiplier  
       * b - constant  
       * c - denominator - 3 (Since its constant across all the operations)  
       *   
       * For each of the steps in the sequence consisting of 'U', 'D', 'd', we store the expression coeffs into   
       * respective arrays  
       *   
       * @throws Exception   
       */  
      private void storeSequenceConstants() throws Exception {  
           for(char sequenceStep: sequence.toCharArray()){  
                switch (sequenceStep){  
                case 'U' : constantCoeffArray[index] = upwardExpression.getConstantCoefficient();  
                numeratorCoeffArray[index++] = upwardExpression.getNumeratorCoefficient();  
                upwardStepCount++; break;  
                case 'd' : constantCoeffArray[index] = downwardExpression.getConstantCoefficient();   
                numeratorCoeffArray[index++] = downwardExpression.getNumeratorCoefficient();  
                downwardStepCount++;break;  
                case 'D' : constantCoeffArray[index] = bigDownwardExpression.getConstantCoefficient();  
                numeratorCoeffArray[index++] = bigDownwardExpression.getNumeratorCoefficient();  
                break;  
                default :   
                     throw new Exception("Prefix Sequence incorrect ");  
                }  
           }  
      }  
   
      /**  
       * This method performs prefix calculation for modified collatz conjecture.  
       * Given a sequence like "DdDddUUdDD....", we first perform the modified Collatz conjecture assuming an initial x.  
       * Based on the expression, we can derive that  
       * After k steps, we will have  
       * Sum = (A (x) + B ) /C  
       *  where               
       * A = a1 ^ U * a2 ^ d * a3 ^ D , where a1 = 4, a2 = 2, a3 = 1 and U,D,d represent the number of times they occur  
       *   in the sequence  
       * B = Sigma[ Pi[multiplierCoeffArray[j]] * constantCoeffArray[i]) * 3^i], i->k-1,0 ; j->1,k-1-i   
       * C = 3 ^ k,   
       *  
       * @param boundRange  
       */  
      private Long calcCollatzExpressionResult(Long boundRange) {  
             
           int sequenceLength = sequence.length();  
           expressionConstantSum = 0d;  
           for(int outer_index=sequenceLength-1; outer_index>= 0; outer_index--){  
                Long product = 1L;  
                for(int inner_index =1; inner_index<= (sequenceLength-1-outer_index); inner_index++) {  
                     product *= numeratorCoeffArray[sequenceLength - inner_index];  
                }  
                expressionConstantSum += product * Math.pow(3, outer_index) * constantCoeffArray[outer_index];  
           }  
           expressionNumeratorSum = Math.pow(4,upwardStepCount) * Math.pow(2,downwardStepCount);  
           expressionDenominator = Math.pow(3,sequenceLength);  
           return calculateResultThroughGCDComputation(expressionNumeratorSum, -expressionDenominator,  
                     -expressionConstantSum, boundRange);  
      }  
   
      /**  
       * Once we compute Ax + B / C as described above, we have gone through the provided input sequence. We want the   
       * smallest value for x greater than a starting Bound (Z) such that Ax + B / C results in an integer Y. This would   
       * ensure that the integer sequence is valid. We find integral solutions for x and Y using the Extended Euclid  
       * algorithm for computing GCD. We obtain equations of the form x = p + q * T , Y = m + n* T, where p, q, m, n are  
       * constants obtained from the GCD and the original equation's coefficients. We can substitute any integer value  
       * for T in order to obtain a valid solution for x. We can choose a good starting point for T by using Z as x,   
       * and then T is either (starting point - 1) or (starting point + 1)  
       *   
       * @param inputCoeffA  
       * @param inputCoeffB  
       * @param inputCoeffC  
       * @param inputBound  
       * @return  
       */  
      private Long calculateResultThroughGCDComputation(Double inputCoeffA, Double inputCoeffB, Double inputCoeffC, Long inputBound){            
           //Check if equation coefficients are negative  
           boolean coeffANegative = inputCoeffA < 0;  
           boolean coeffBNegative = inputCoeffB < 0;  
   
           // Storing the coefficients of the input equation as a BigInteger  
           BigInteger coeffA = BigInteger.valueOf(Math.abs(inputCoeffA.longValue()));  
           BigInteger coeffB = BigInteger.valueOf(Math.abs(inputCoeffB.longValue()));  
           BigInteger constantC = BigInteger.valueOf(inputCoeffC.longValue());  
             
           BigInteger coeffATemp = coeffA;  
           BigInteger coeffBTemp = coeffB;  
           BigInteger lowerBound = BigInteger.valueOf(inputBound);  
   
           BigInteger currentX = new BigInteger("0");  
           BigInteger currentY = new BigInteger("1");  
           BigInteger previousX = new BigInteger("1");  
           BigInteger previousY = new BigInteger("0");  
           BigInteger temp;  
   
           BigInteger remainder;  
           BigInteger quotient;  
             
           BigInteger solution;  
           BigInteger startingValue;  
           BigInteger xVar;  
             
           // Extended Euclid Algorithm to calculate GCD as well as the expression of GCD in terms of input coeffs  
           while (!coeffBTemp.equals(BigInteger.ZERO)) {  
                quotient = coeffATemp.divide(coeffBTemp);  
                remainder = coeffATemp.mod(coeffBTemp);  
                coeffATemp = coeffBTemp;  
                coeffBTemp = remainder;  
   
                temp = currentX;  
                currentX = previousX.subtract(quotient.multiply(currentX));  
                previousX = temp;  
                temp = currentY;  
                currentY = previousY.subtract(quotient.multiply(currentY));  
                previousY = temp;        
           }  
             
           // If coefficients were negative then transferring the sign to the roots of the equation  
           if(coeffANegative){  
                previousX = BigInteger.ZERO.subtract(previousX);  
           }  
           if(coeffBNegative){  
                previousY = BigInteger.ZERO.subtract(previousY);  
           }  
   
           //The coefficient of the variable part of the generic equation for all solutions of X  
           //i.e if x = p + q * T, this represents q;  
           xVar = coeffB.divide(coeffATemp);  
   
           //Scale the X and Y values to account for the constant in the equation  
           previousX = previousX.multiply(constantC);  
           previousY = previousY.multiply(constantC);  
   
           startingValue = (lowerBound.subtract(previousX)).divide(xVar);  
             
           // Calculate the generic equation for all solutions of X using the initialValue and xVar(coeffB/GCD)  
           solution = previousX.add(xVar.multiply(startingValue));  
             
           startingValue = startingValue.add(BigInteger.ONE);  
             
           while(solution.compareTo(lowerBound) == -1) {  
                solution = previousX.add(xVar.multiply(startingValue));  
                startingValue = startingValue.subtract(BigInteger.ONE);  
           }  
           return solution.longValue();  
      }  
   
      /**  
       * This method calculates the smallest number < than the specified bound that is prefixed with the given sequence   
       *  
       * @throws Exception   
       */  
      public Long computeSmallestNumWithSequence(Long boundingInput) throws Exception{  
           if(boundingInput != null  &&  boundingInput >= 0 && sequence != null) {  
                storeSequenceConstants();  
                return calcCollatzExpressionResult(boundingInput);  
           } else{  
                throw new Exception("Computation Input is Invalid");  
           }  
      }