package umontreal.iro.lecuyer.hups;
/* *
 * Title:        

 * Description:  

 * Copyright:    Copyright (c) 

 * Company:      

 * @author
 * @version 1.0
 */


/**
 * This class implements the sequence of Halton,
 * which is essentially a modification of Hammersley nets for producing 
 * an infinite sequence of points having low discrepancy.
 * The ith point in s dimensions is 
 * 
 * 

 * 

 * i = (ψb1(i), ψb2(i),..., ψbs(i)),
 * 

 * for 
 * i = 0, 1, 2,..., where ψb is the radical inverse function
 * in base b, defined in class {@link RadicalInverse}, and where
 * 
 * 2 = b1 <  ...  < bs are the s smallest prime numbers in 
 * increasing order.
 * 
 * 

 * A fast method is implemented to generate randomized Halton sequences, starting from an arbitrary points x0.
 * 
 * 

 * The points can be ``scrambled'' by applying a permutation to the 
 * digits of i before computing each coordinate, in the same way as for the class
 * {@link HammersleyPointSet}, for all coordinates j >=  0.
 * 
 */
public class HaltonSequence extends PointSet { 
   private int[] base;           // Vector of prime bases.
   private int[][] permutation;  // Digits permutation, for each dimension.
   private boolean permuted;     // Permute digits?
   private RadicalInverse[] radinv; // Vector of RadicalInverse's.


   /**
    * Constructs a new Halton sequence 
    * in dim dimensions.
    * 
    * @param dim dimension
    * 
    */
   public HaltonSequence (int dim)  {
      if (dim < 1)
         throw new IllegalArgumentException
            ("Halton sequence must have positive dimension dim");
      this.dim  = dim;
      numPoints = Integer.MAX_VALUE;
      base = RadicalInverse.getPrimes (dim);
   }


   /**
    * Initializes the Halton sequence starting at point x0.
    *    The dimension of x0 must be at least as large as the dimension
    *    of this object.
    * 
    * @param x0 starting point of the Halton sequence
    * 
    * 
    */
   public void init (double[] x0)  {
      radinv = new RadicalInverse[dim];
      for (int i = 0; i < dim; i++)
         radinv[i] = new RadicalInverse (base[i], x0[i]);
   }


   /**
    * Permutes the digits using Faure permutations for all coordinates.
    *   After the method is called, the coordinates ui, j are generated via
    * 
    * 

    * 

    * ui, j = ∑r=0k-1πj[ar]bj-r-1,
    * 

    * for 
    * j = 0,..., s - 1,
    *  where πj is the Faure permutation of 
    * {0,..., bj -1}.
    * 
    */
   public void addFaurePermutations() {
      permutation = new int[dim][];
      for (int i = 0; i < dim; i++) {
         permutation[i] = new int[base[i]];
         RadicalInverse.getFaurePermutation (base[i], permutation[i]);
      }
      permuted = true;
   }



   /**
    * Erases the Faure permutations: from now on, the digits will not be
    *   Faure permuted.
    * 
    */
   public void ErasePermutations() {
      permuted = false;
      permutation = null;
   }



    
   public int getNumPoints () {
      return Integer.MAX_VALUE;
   }

   public double getCoordinate (int i, int j) {
      if (radinv != null) {
         if (!permuted) {
            return radinv[j].nextRadicalInverse ();
         } else {
            throw new UnsupportedOperationException (
            "Fast radical inverse is not implemented in case of permutation");
         }
      } else {
         if (permuted)
            return RadicalInverse.permutedRadicalInverse 
                      (base[j], permutation[j], i);
         else 
            return RadicalInverse.radicalInverse (base[j], i);
      }
   }
}