SearcherCBC.java

/*
 * Class:        SearcherCBC
 * Description:  searches the best rank-1 lattices with respect to a given
                 discrepancy, using component-by-component (CBC) searches.
 * Environment:  Java
 * Software:     SSJ 
 * Copyright (C) 2001  Pierre L'Ecuyer and Universite de Montreal
 * Organization: DIRO, Universite de Montreal
 * @author       Richard Simard
 * @since        January 2009
 
 * SSJ is free software: you can redistribute it and/or modify it under
 * the terms of the GNU General Public License (GPL) as published by the
 * Free Software Foundation, either version 3 of the License, or
 * any later version.
 
 * SSJ is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 
 * A copy of the GNU General Public License is available at
   <a href="http://www.gnu.org/licenses">GPL licence site</a>.
 */
package umontreal.ssj.discrepancy;
 
import umontreal.ssj.hups.Rank1Lattice;
import umontreal.ssj.util.Num;

public class SearcherCBC extends Searcher {
   protected double[] bestVals; // best value of discr in dim 0, 1, ..., (s-1)
 
   private double exhaust(int s, boolean relPrime) {
      int n = disc.getNumPoints();
      gamma = disc.getGamma();
      int pos = -1;
      double err;
      double best = 0;
      int i, j;
      bestAs[0] = 1;
      bestVals[0] = -1;
 
      for (j = 1; j < s; j++) {
         best = Double.MAX_VALUE;
         pos = -1;
         for (i = 1; i < n; i++) {
            if (relPrime) {
               // Consider only values i relatively prime to n
               if (power2F) {
                  if ((i & 1) == 0)
                     continue;
               } else {
                  if (Num.gcd(n, i) != 1)
                     continue;
               }
            }
            bestAs[j] = i;
            // print (bestAs, j + 1);
            lat = new Rank1Lattice(n, bestAs, j + 1);
            err = disc.compute(lat, gamma);
            if (err < best) {
               best = err;
               pos = i;
            }
         }
         bestAs[j] = pos;
         bestVals[j] = best;
      }
 
      bestVal = best;
      return best;
   }
 
   private double random(int s, int k, boolean relPrime) {
      int n = disc.getNumPoints();
      if (k >= n)
         return exhaust(s, relPrime);
 
      gamma = disc.getGamma();
      final int nm1 = n - 1;
      int pos = -1;
      double err;
      double best = 0;
      int i, j;
      bestAs[0] = 1;
      bestVals[0] = -1;
 
      for (j = 1; j < s; j++) {
         best = Double.MAX_VALUE;
         i = 0;
         while (i < k) {
            if (power2F) {
               bestAs[j] = gen.nextInt(1, nm1);
               if (relPrime)
                  bestAs[j] |= 1;
            } else {
               do {
                  bestAs[j] = gen.nextInt(1, nm1);
               } while (relPrime && (Num.gcd(n, bestAs[j]) != 1));
            }
            i++;
            // print (bestAs, j + 1);
            lat = new Rank1Lattice(n, bestAs, j + 1);
            err = disc.compute(lat, gamma);
            if (err < best) {
               best = err;
               pos = bestAs[j];
            }
         }
         bestAs[j] = pos;
         bestVals[j] = best;
      }
 
      bestVal = best;
      return best;
   }

   public SearcherCBC(Discrepancy disc, boolean primeN) {
      super(disc, primeN);
      int s = disc.getDimension();
      bestVals = new double[s];
   }

   public double exhaust(int s) {
      return exhaust(s, false);
   }

   public double exhaustPrime(int s) {
      if (primeN)
         return exhaust(s, false); // all a_j are rel. prime to n
      else
         return exhaust(s, true);
   }

   public double random(int s, int k) {
      return random(s, k, false);
   }

   public double randomPrime(int s, int k) {
      if (primeN)
         return random(s, k, false); // all a_j are rel. prime to n
      else
         return random(s, k, true);
   }

   public double[] getBestVals() {
      return bestVals;
   }
 
}