ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/jsr166/jsr166/src/test/loops/TSPExchangerTest.java
(Generate patch)

Comparing jsr166/src/test/loops/TSPExchangerTest.java (file contents):
Revision 1.1 by dl, Sun Aug 7 19:25:55 2005 UTC vs.
Revision 1.19 by jsr166, Sun Sep 13 16:28:14 2015 UTC

# Line 1 | Line 1
1   /*
2 < * Written by Bill Scherer and Doug Lea with assistance from members
3 < * of JCP JSR-166 Expert Group and released to the public domain. Use,
4 < * modify, and redistribute this code in any way without
5 < * acknowledgement.
2 > * Written by Doug Lea and Bill Scherer with assistance from members
3 > * of JCP JSR-166 Expert Group and released to the public domain, as
4 > * explained at http://creativecommons.org/publicdomain/zero/1.0/
5   */
6  
7 + import java.util.*;
8   import java.util.concurrent.*;
9   import java.util.concurrent.atomic.*;
10   import java.util.concurrent.locks.*;
11  
12 + /**
13 + * A parallel Traveling Salesperson Problem (TSP) program based on a
14 + * genetic algorithm using an Exchanger.  A population of chromosomes is
15 + * distributed among "subpops".  Each chromosomes represents a tour,
16 + * and its fitness is the total tour length.
17 + *
18 + * A set of worker threads perform updates on subpops. The basic
19 + * update step is:
20 + * <ol>
21 + *   <li>Select a breeder b from the subpop
22 + *   <li>Create a strand of its tour with a random starting point and length
23 + *   <li>Offer the strand to the exchanger, receiving a strand from
24 + *       another subpop
25 + *   <li>Combine b and the received strand using crossing function to
26 + *       create new chromosome c.
27 + *   <li>Replace a chromosome in the subpop with c.
28 + * </ol>
29 + *
30 + * This continues for a given number of generations per subpop.
31 + * Because there are normally more subpops than threads, each worker
32 + * thread performs small (randomly sized) run of updates for one
33 + * subpop and then selects another. A run continues until there is at
34 + * most one remaining thread performing updates.
35 + *
36 + * See below for more details.
37 + */
38   public class TSPExchangerTest {
39 <    // Set SLS true to use as default the settings in Scherer, Lea, and
40 <    // Scott paper. Otherwise smaller values are used to speed up testing
41 <    static final boolean SLS = false;
42 <
43 <    static final int DEFAULT_THREADS      = SLS?    32:     8;
44 <    static final int DEFAULT_CITIES       = SLS?   100:    50;
45 <    static final int DEFAULT_POPULATION   = SLS?  1000:   500;
46 <    static final int DEFAULT_BREEDERS     = SLS?   200:   100;
47 <    static final int DEFAULT_GENERATIONS  = SLS? 20000: 10000;
39 >    static final int NCPUS = Runtime.getRuntime().availableProcessors();
40 >
41 >    /** Runs start with two threads, increasing by two through max */
42 >    static final int DEFAULT_MAX_THREADS = Math.max(4, NCPUS + NCPUS/2);
43 >
44 >    /** The number of replication runs per thread value */
45 >    static final int DEFAULT_REPLICATIONS = 3;
46 >
47 >    /** If true, print statistics in SNAPSHOT_RATE intervals */
48 >    static boolean verbose = true;
49 >    static final long SNAPSHOT_RATE = 10000; // in milliseconds
50 >
51 >    /**
52 >     * The problem size. Each city is a random point. The goal is to
53 >     * find a tour among them with smallest total Euclidean distance.
54 >     */
55 >    static final int DEFAULT_CITIES = 144;
56 >
57 >    // Tuning parameters.
58 >
59 >    /**
60 >     * The number of chromosomes per subpop. Must be a power of two.
61 >     *
62 >     * Smaller values lead to faster iterations but poorer quality
63 >     * results
64 >     */
65 >    static final int DEFAULT_SUBPOP_SIZE = 32;
66 >
67 >    /**
68 >     * The number of iterations per subpop. Convergence appears
69 >     * to be roughly proportional to #cities-squared
70 >     */
71 >    static final int DEFAULT_GENERATIONS = DEFAULT_CITIES * DEFAULT_CITIES;
72 >
73 >    /**
74 >     * The number of subpops. The total population is #subpops * subpopSize,
75 >     * which should be roughly on the order of #cities-squared
76 >     *
77 >     * Smaller values lead to faster total runs but poorer quality
78 >     * results
79 >     */
80 >    static final int DEFAULT_NSUBPOPS = DEFAULT_GENERATIONS / DEFAULT_SUBPOP_SIZE;
81 >
82 >    /**
83 >     * The minimum length for a random chromosome strand.
84 >     * Must be at least 1.
85 >     */
86 >    static final int MIN_STRAND_LENGTH = 3;
87 >
88 >    /**
89 >     * The probability mask value for creating random strands,
90 >     * that have lengths at least MIN_STRAND_LENGTH, and grow
91 >     * with exponential decay 2^(-(1/(RANDOM_STRAND_MASK + 1)
92 >     * Must be 1 less than a power of two.
93 >     */
94 >    static final int RANDOM_STRAND_MASK = 7;
95 >
96 >    /**
97 >     * Probability control for selecting breeders.
98 >     * Breeders are selected starting at the best-fitness chromosome,
99 >     * with exponentially decaying probability
100 >     * 1 / (subpopSize >>> BREEDER_DECAY).
101 >     *
102 >     * Larger values usually cause faster convergence but poorer
103 >     * quality results
104 >     */
105 >    static final int BREEDER_DECAY = 1;
106 >
107 >    /**
108 >     * Probability control for selecting dyers.
109 >     * Dyers are selected starting at the worst-fitness chromosome,
110 >     * with exponentially decaying probability
111 >     * 1 / (subpopSize >>> DYER_DECAY)
112 >     *
113 >     * Larger values usually cause faster convergence but poorer
114 >     * quality results
115 >     */
116 >    static final int DYER_DECAY = 1;
117  
118 +    /**
119 +     * The set of cities. Created once per program run, to
120 +     * make it easier to compare solutions across different runs.
121 +     */
122 +    static CitySet cities;
123  
124      public static void main(String[] args) throws Exception {
125 <        int maxThreads = DEFAULT_THREADS;
125 >        int maxThreads = DEFAULT_MAX_THREADS;
126          int nCities = DEFAULT_CITIES;
127 <        int pSize = DEFAULT_POPULATION;
128 <        int nBreeders = DEFAULT_BREEDERS;
129 <        int numGenerations = DEFAULT_GENERATIONS;
127 >        int subpopSize = DEFAULT_SUBPOP_SIZE;
128 >        int nGen = nCities * nCities;
129 >        int nSubpops = nCities * nCities / subpopSize;
130 >        int nReps = DEFAULT_REPLICATIONS;
131  
31        // Parse and check args
32        int argc = 0;
132          try {
133 +            int argc = 0;
134              while (argc < args.length) {
135                  String option = args[argc++];
136 <                if (option.equals("-b"))
37 <                    nBreeders = Integer.parseInt(args[argc]);
38 <                else if (option.equals("-c"))
136 >                if (option.equals("-c")) {
137                      nCities = Integer.parseInt(args[argc]);
138 +                    nGen = nCities * nCities;
139 +                    nSubpops = nCities * nCities / subpopSize;
140 +                }
141 +                else if (option.equals("-p"))
142 +                    subpopSize = Integer.parseInt(args[argc]);
143                  else if (option.equals("-g"))
144 <                    numGenerations = Integer.parseInt(args[argc]);
145 <                else if (option.equals("-p"))
146 <                    pSize = Integer.parseInt(args[argc]);
147 <                else
144 >                    nGen = Integer.parseInt(args[argc]);
145 >                else if (option.equals("-n"))
146 >                    nSubpops = Integer.parseInt(args[argc]);
147 >                else if (option.equals("-q")) {
148 >                    verbose = false;
149 >                    argc--;
150 >                }
151 >                else if (option.equals("-r"))
152 >                    nReps = Integer.parseInt(args[argc]);
153 >                else
154                      maxThreads = Integer.parseInt(option);
155                  argc++;
156              }
157          }
49        catch (NumberFormatException e) {
50            reportUsageErrorAndDie();
51            System.exit(0);
52        }
158          catch (Exception e) {
159              reportUsageErrorAndDie();
160          }
161  
162 <        // Display runtime parameters
163 <        System.out.print("TSPExchangerTest -b " + nBreeders);
164 <        System.out.print(" -c " + nCities);
165 <        System.out.print(" -g " + numGenerations);
166 <        System.out.print(" -p " + pSize);
167 <        System.out.print(" max threads " + maxThreads);
168 <        System.out.println();
169 <
170 <        // warmup
171 <        System.out.print("Threads: " + 2 + "\t");
172 <        oneRun(2,
173 <               nCities,
174 <               pSize,
175 <               nBreeders,
176 <               numGenerations);
177 <        Thread.sleep(100);
178 <
179 <        int k = 4;
180 <        for (int i = 2; i <= maxThreads;) {
181 <            System.out.print("Threads: " + i + "\t");
182 <            oneRun(i,
183 <                   nCities,
184 <                   pSize,
185 <                   nBreeders,
186 <                   numGenerations);
82 <            Thread.sleep(100);
83 <            if (i == k) {
84 <                k = i << 1;
85 <                i = i + (i >>> 1);
86 <            }
87 <            else
88 <                i = k;
162 >        System.out.print("TSPExchangerTest");
163 >        System.out.print(" -c " + nCities);
164 >        System.out.print(" -g " + nGen);
165 >        System.out.print(" -p " + subpopSize);
166 >        System.out.print(" -n " + nSubpops);
167 >        System.out.print(" -r " + nReps);
168 >        System.out.print(" max threads " + maxThreads);
169 >        System.out.println();
170 >
171 >        cities = new CitySet(nCities);
172 >
173 >        if (false && NCPUS > 4) {
174 >            int h = NCPUS/2;
175 >            System.out.printf("Threads: %4d Warmup\n", h);
176 >            oneRun(h, nSubpops, subpopSize, nGen);
177 >            Thread.sleep(500);
178 >        }
179 >
180 >        int maxt = (maxThreads < nSubpops) ? maxThreads : nSubpops;
181 >        for (int j = 0; j < nReps; ++j) {
182 >            for (int i = 2; i <= maxt; i += 2) {
183 >                System.out.printf("Threads: %4d Replication: %2d\n", i, j);
184 >                oneRun(i, nSubpops, subpopSize, nGen);
185 >                Thread.sleep(500);
186 >            }
187          }
188      }
189  
190 <    private static void reportUsageErrorAndDie() {
191 <        System.out.print("usage: TSPExchangerTest [-b #breeders] [-c #cities]");
192 <        System.out.println(" [-g #generations]");
193 <        System.out.println(" [-p population size] [ #threads]");
190 >    static void reportUsageErrorAndDie() {
191 >        System.out.print("usage: TSPExchangerTest");
192 >        System.out.print(" [-c #cities]");
193 >        System.out.print(" [-p #subpopSize]");
194 >        System.out.print(" [-g #generations]");
195 >        System.out.print(" [-n #subpops]");
196 >        System.out.print(" [-r #replications]");
197 >        System.out.print(" [-q <quiet>]");
198 >        System.out.print(" #threads]");
199 >        System.out.println();
200          System.exit(0);
201      }
202  
203 <    static void oneRun(int nThreads,
204 <                       int nCities,
205 <                       int pSize,
206 <                       int nBreeders,
207 <                       int numGenerations)
208 <        throws Exception {
209 <        CyclicBarrier runBarrier = new CyclicBarrier(nThreads + 1);
210 <        Population p = new Population(nCities, pSize, nBreeders, nThreads,
211 <                                      numGenerations, runBarrier);
212 <
213 <        // Run the test
214 <        long startTime = System.currentTimeMillis();
215 <        runBarrier.await(); // start 'em off
216 <        runBarrier.await(); // wait 'til they're done
217 <        long stopTime = System.currentTimeMillis();
218 <        long elapsed = stopTime - startTime;
219 <        long rate = (numGenerations * 1000) / elapsed;
220 <        double secs = (double)elapsed / 1000.0;
203 >    /**
204 >     * Performs one run with the given parameters.  Each run completes
205 >     * when there are fewer than 2 active threads.  When there is
206 >     * only one remaining thread, it will have no one to exchange
207 >     * with, so it is terminated (via interrupt).
208 >     */
209 >    static void oneRun(int nThreads, int nSubpops, int subpopSize, int nGen)
210 >        throws InterruptedException {
211 >        Population p = new Population(nThreads, nSubpops, subpopSize, nGen);
212 >        ProgressMonitor mon = null;
213 >        if (verbose) {
214 >            p.printSnapshot(0);
215 >            mon = new ProgressMonitor(p);
216 >            mon.start();
217 >        }
218 >        long startTime = System.nanoTime();
219 >        p.start();
220 >        p.awaitDone();
221 >        long stopTime = System.nanoTime();
222 >        if (mon != null)
223 >            mon.interrupt();
224 >        p.shutdown();
225 >        //        Thread.sleep(100);
226  
227 <        // Display results
228 <        System.out.print(LoopHelpers.rightJustify((int)p.bestFitness()) +
229 <                         " fitness");
121 <        System.out.print(LoopHelpers.rightJustify(rate) + " gen/s \t");
122 <        System.out.print(secs + "s elapsed");
123 <        System.out.println();
227 >        long elapsed = stopTime - startTime;
228 >        double secs = (double) elapsed / 1000000000.0;
229 >        p.printSnapshot(secs);
230      }
231  
232 +    /**
233 +     * A Population creates the subpops, subpops, and threads for a run
234 +     * and has control methods to start, stop, and report progress.
235 +     */
236      static final class Population {
237 <        final Chromosome[] individuals;
238 <        final Exchanger<Chromosome> x;
239 <        final CitySet cities;
240 <        final int[] dyers;
241 <        final int[] breeders;
242 <        final CyclicBarrier generationBarrier;
133 <        final Thread[] threads;
134 <        final boolean[] doneMating;
135 <        final ReentrantLock matingBarrierLock;
136 <        final Condition matingBarrier;
137 <        final LoopHelpers.SimpleRandom[] rngs;
237 >        final Worker[] threads;
238 >        final Subpop[] subpops;
239 >        final Exchanger<Strand> exchanger;
240 >        final CountDownLatch done;
241 >        final int nGen;
242 >        final int subpopSize;
243          final int nThreads;
139        volatile int matingBarrierCount;
140
141        // action to run between each generation
142        class BarrierAction implements Runnable {
143            public void run() {
144                prepareToBreed();
145                resetMatingBarrier();
146            }
147        }
244  
245 <        Population(int nCities,
150 <                   int pSize,
151 <                   int nBreeders,
152 <                   int nThreads,
153 <                   int nGen,
154 <                   CyclicBarrier runBarrier) {
245 >        Population(int nThreads, int nSubpops, int subpopSize, int nGen) {
246              this.nThreads = nThreads;
247 <            // rngs[nThreads] is for global actions; others are per-thread
248 <            this.rngs = new LoopHelpers.SimpleRandom[nThreads+1];
249 <            for (int i = 0; i < rngs.length; ++i)
250 <                rngs[i] = new LoopHelpers.SimpleRandom();
251 <            this.cities = new CitySet(nCities, rngs[nThreads]);
252 <            this.individuals = new Chromosome[pSize];
253 <            for (int i = 0; i < individuals.length; i++)
254 <                individuals[i] = new Chromosome(cities, nCities,
255 <                                                rngs[nThreads]);
256 <            this.doneMating = new boolean[nThreads];
257 <            this.dyers = new int[nBreeders];
258 <            this.breeders = new int[nBreeders];
259 <
260 <            this.x = new Exchanger();
261 <
262 <            this.matingBarrierLock = new ReentrantLock();
263 <            this.matingBarrier = matingBarrierLock.newCondition();
264 <
265 <            BarrierAction ba = new BarrierAction();
266 <            this.generationBarrier = new CyclicBarrier(nThreads, ba);
267 <            ba.run(); // prepare for first generation
268 <
269 <            this.threads = new Thread[nThreads];
270 <            for (int i = 0; i < nThreads; i++) {
271 <                Runner r = new Runner(i, this, runBarrier, nGen);
272 <                threads[i] = new Thread(r);
273 <                threads[i].start();
274 <            }
275 <        }
276 <        
277 <        double averageFitness() {
278 <            double total = 0;
279 <            for (int i = 0; i < individuals.length; i++)
280 <                total += individuals[i].fitness;
281 <            return total/(double)individuals.length;
282 <        }
283 <
284 <        double bestFitness() {
285 <            double best = individuals[0].fitness;
286 <            for (int i = 0; i < individuals.length; i++)
287 <                if (individuals[i].fitness < best)
288 <                    best = individuals[i].fitness;
289 <            return best;
290 <        }
291 <
292 <        void resetMatingBarrier() {
293 <            matingBarrierCount = nThreads - 1;
294 <        }
295 <
296 <        void awaitMatingBarrier(int tid) {
297 <            doneMating[tid] = true; // heuristically set before lock
298 <            matingBarrierLock.lock();
299 <            try {
300 <                int m = matingBarrierCount--;
301 <                if (m < 1) {
302 <                    for (int i = 0; i < doneMating.length; ++i)
303 <                        doneMating[i] = false;
304 <                    Thread.interrupted(); // clear
305 <                    matingBarrier.signalAll();
306 <                } else {
307 <                    doneMating[tid] = true;
308 <                    if (m == 1 && nThreads > 2) {
309 <                        for (int j = 0; j < doneMating.length; ++j) {
310 <                            if (!doneMating[j]) {
311 <                                threads[j].interrupt();
312 <                                break;
313 <                            }
314 <                        }
315 <                    }
316 <                    try {
317 <                        do {
318 <                            matingBarrier.await();
319 <                        } while (matingBarrierCount >= 0);
320 <                    } catch(InterruptedException ie) {}
230 <                }
231 <            } finally {
232 <                matingBarrierLock.unlock();
233 <            }
247 >            this.nGen = nGen;
248 >            this.subpopSize = subpopSize;
249 >            this.exchanger = new Exchanger<Strand>();
250 >            this.done = new CountDownLatch(nThreads - 1);
251 >
252 >            this.subpops = new Subpop[nSubpops];
253 >            for (int i = 0; i < nSubpops; i++)
254 >                subpops[i] = new Subpop(this);
255 >
256 >            this.threads = new Worker[nThreads];
257 >            int maxExchanges = nGen * nSubpops / nThreads;
258 >            for (int i = 0; i < nThreads; ++i) {
259 >                threads[i] = new Worker(this, maxExchanges);
260 >            }
261 >
262 >        }
263 >
264 >        void start() {
265 >            for (int i = 0; i < nThreads; ++i) {
266 >                threads[i].start();
267 >            }
268 >        }
269 >
270 >        /** Stop the tasks */
271 >        void shutdown() {
272 >            for (int i = 0; i < threads.length; ++ i)
273 >                threads[i].interrupt();
274 >        }
275 >
276 >        void threadDone() {
277 >            done.countDown();
278 >        }
279 >
280 >        /** Wait for tasks to complete */
281 >        void awaitDone() throws InterruptedException {
282 >            done.await();
283 >        }
284 >
285 >        int totalExchanges() {
286 >            int xs = 0;
287 >            for (int i = 0; i < threads.length; ++i)
288 >                xs += threads[i].exchanges;
289 >            return xs;
290 >        }
291 >
292 >        /**
293 >         * Prints statistics, including best and worst tour lengths
294 >         * for points scaled in [0,1), scaled by the square root of
295 >         * number of points. This simplifies checking results.  The
296 >         * expected optimal TSP for random points is believed to be
297 >         * around 0.76 * sqrt(N). For papers discussing this, see
298 >         * http://www.densis.fee.unicamp.br/~moscato/TSPBIB_home.html
299 >         */
300 >        void printSnapshot(double secs) {
301 >            int xs = totalExchanges();
302 >            long rate = (xs == 0) ? 0L : (long) ((secs * 1000000000.0) / xs);
303 >            Chromosome bestc = subpops[0].chromosomes[0];
304 >            Chromosome worstc = bestc;
305 >            for (int k = 0; k < subpops.length; ++k) {
306 >                Chromosome[] cs = subpops[k].chromosomes;
307 >                if (cs[0].fitness < bestc.fitness)
308 >                    bestc = cs[0];
309 >                int w = cs[cs.length-1].fitness;
310 >                if (cs[cs.length-1].fitness > worstc.fitness)
311 >                    worstc = cs[cs.length-1];
312 >            }
313 >            double sqrtn = Math.sqrt(cities.length);
314 >            double best = bestc.unitTourLength() / sqrtn;
315 >            double worst = worstc.unitTourLength() / sqrtn;
316 >            System.out.printf("N:%4d T:%8.3f B:%6.3f W:%6.3f X:%9d R:%7d\n",
317 >                              nThreads, secs, best, worst, xs, rate);
318 >            //            exchanger.printStats();
319 >            //            System.out.print(" s: " + exchanger.aveSpins());
320 >            //            System.out.print(" p: " + exchanger.aveParks());
321          }
322 +    }
323  
324 <        void prepareToBreed() {
324 >    /**
325 >     * Worker threads perform updates on subpops.
326 >     */
327 >    static final class Worker extends Thread {
328 >        final Population pop;
329 >        final int maxExchanges;
330 >        int exchanges;
331 >        final RNG rng = new RNG();
332 >
333 >        Worker(Population pop, int maxExchanges) {
334 >            this.pop = pop;
335 >            this.maxExchanges = maxExchanges;
336 >        }
337  
338 <            // Calculate statistics
339 <            double totalFitness = 0;
340 <            double worstFitness = 0;
341 <            double bestFitness = individuals[0].fitness;
342 <
343 <            for (int i = 0; i < individuals.length; i++) {
344 <                totalFitness += individuals[i].fitness;
345 <                if (individuals[i].fitness > worstFitness)
346 <                    worstFitness = individuals[i].fitness;
347 <                if (individuals[i].fitness < bestFitness)
348 <                    bestFitness = individuals[i].fitness;
349 <            }
350 <            
351 <            double[] lifeNorm = new double[individuals.length];
352 <            double lifeNormTotal = 0;
253 <            double[] deathNorm = new double[individuals.length];
254 <            double deathNormTotal = 0;
255 <            for (int i = 0; i < individuals.length; i++) {
256 <                deathNorm[i] = (individuals[i].fitness - bestFitness)
257 <                    / (worstFitness - bestFitness + 1) + .05;
258 <                deathNorm[i] = (deathNorm[i] * deathNorm[i]);
259 <                lifeNorm[i] = 1.0 - deathNorm[i];
260 <                lifeNormTotal += lifeNorm[i];
261 <                deathNormTotal += deathNorm[i];
262 <            }
263 <            
264 <            double deathScale = deathNormTotal / (double)0x7FFFFFFF;
265 <            double lifeScale = lifeNormTotal / (double)0x7FFFFFFF;
266 <            
267 <            int nSub = breeders.length / nThreads;
268 <            LoopHelpers.SimpleRandom random = rngs[nThreads];
269 <
270 <            // Select breeders (need to be distinct)
271 <            for (int i = 0; i < nSub; i++) {
272 <                boolean newBreeder;
273 <                int lucky;
274 <                do {
275 <                    newBreeder = true;
276 <                    double choice = lifeScale * (double)random.next();
277 <                    for (lucky = 0; lucky < individuals.length; lucky++) {
278 <                        choice -= lifeNorm[lucky];
279 <                        if (choice <= 0)
280 <                            break;
281 <                    }
282 <                    for (int j = 0; j < i; j++)
283 <                        if (breeders[j] == lucky)
284 <                            newBreeder = false;
285 <                } while (!newBreeder);
286 <                breeders[i] = lucky;
287 <            }
288 <            
289 <            // Select dead guys (need to be distinct)
290 <            for (int i = 0; i < nSub; i++) {
291 <                boolean newDead;
292 <                int victim;
293 <                do {
294 <                    newDead = true;
295 <                    double choice = deathScale * (double)random.next();
296 <                    for (victim = 0; victim < individuals.length; victim++) {
297 <                        choice -= deathNorm[victim];
298 <                        if (choice <= 0)
299 <                            break;
338 >        /**
339 >         * Repeatedly, find a subpop that is not being updated by
340 >         * another thread, and run a random number of updates on it.
341 >         */
342 >        public void run() {
343 >            try {
344 >                int len = pop.subpops.length;
345 >                int pos = (rng.next() & 0x7FFFFFFF) % len;
346 >                while (exchanges < maxExchanges) {
347 >                    Subpop s = pop.subpops[pos];
348 >                    AtomicBoolean busy = s.busy;
349 >                    if (!busy.get() && busy.compareAndSet(false, true)) {
350 >                        exchanges += s.runUpdates();
351 >                        busy.set(false);
352 >                        pos = (rng.next() & 0x7FFFFFFF) % len;
353                      }
354 <                    for (int j = 0; j < i; j++)
355 <                        if (dyers[j] == victim)
356 <                            newDead = false;
357 <                } while (!newDead);
358 <                dyers[i] = victim;
354 >                    else if (++pos >= len)
355 >                        pos = 0;
356 >                }
357 >                pop.threadDone();
358 >            } catch (InterruptedException fallthrough) {
359              }
307
360          }
361 +    }
362  
363 <        
364 <        void nextGeneration(int tid, int matings)
365 <            throws InterruptedException, BrokenBarrierException {
366 <
367 <            int firstChild = ((individuals.length * tid) / nThreads);
368 <            int lastChild = ((individuals.length * (tid + 1)) / nThreads);
369 <            int nChildren =  lastChild - firstChild;
370 <            
371 <            int firstSub = ((breeders.length * tid) / nThreads);
372 <            int lastSub = ((breeders.length * (tid + 1)) / nThreads);
373 <            int nSub = lastSub - firstSub;
374 <            
375 <            Chromosome[] children = new Chromosome[nChildren];
376 <
377 <            LoopHelpers.SimpleRandom random = rngs[tid];
378 <
379 <            for (int i = 0; i < nSub; i++) {
380 <                Chromosome parent = individuals[breeders[firstSub + i]];
381 <                Chromosome offspring = new Chromosome(parent);
382 <                int k = 0;
383 <                while (k < matings && matingBarrierCount > 0) {
384 <                    try {
385 <                        Chromosome other = x.exchange(offspring);
386 <                        offspring = offspring.reproduceWith(other, random);
387 <                        ++k;
388 <                    } catch (InterruptedException to) {
389 <                        break;
390 <                    }
363 >    /**
364 >     * A Subpop maintains a set of chromosomes.
365 >     */
366 >    static final class Subpop {
367 >        /** The chromosomes, kept in sorted order */
368 >        final Chromosome[] chromosomes;
369 >        /** The parent population */
370 >        final Population pop;
371 >        /** Reservation bit for worker threads */
372 >        final AtomicBoolean busy;
373 >        /** The common exchanger, same for all subpops */
374 >        final Exchanger<Strand> exchanger;
375 >        /** The current strand being exchanged */
376 >        Strand strand;
377 >        /** Bitset used in cross */
378 >        final int[] inTour;
379 >        final RNG rng;
380 >        final int subpopSize;
381 >
382 >        Subpop(Population pop) {
383 >            this.pop = pop;
384 >            this.subpopSize = pop.subpopSize;
385 >            this.exchanger = pop.exchanger;
386 >            this.busy = new AtomicBoolean(false);
387 >            this.rng = new RNG();
388 >            int length = cities.length;
389 >            this.strand = new Strand(length);
390 >            this.inTour = new int[(length >>> 5) + 1];
391 >            this.chromosomes = new Chromosome[subpopSize];
392 >            for (int j = 0; j < subpopSize; ++j)
393 >                chromosomes[j] = new Chromosome(length, rng);
394 >            Arrays.sort(chromosomes);
395 >        }
396 >
397 >        /**
398 >         * Run a random number of updates.  The number of updates is
399 >         * at least 1 and no more than subpopSize.  This
400 >         * controls the granularity of multiplexing subpop updates on
401 >         * to threads. It is small enough to balance out updates
402 >         * across tasks, but large enough to avoid having runs
403 >         * dominated by subpop selection. It is randomized to avoid
404 >         * long runs where pairs of subpops exchange only with each
405 >         * other.  It is hardwired because small variations of it
406 >         * don't matter much.
407 >         *
408 >         * @param g the first generation to run
409 >         */
410 >        int runUpdates() throws InterruptedException {
411 >            int n = 1 + (rng.next() & ((subpopSize << 1) - 1));
412 >            for (int i = 0; i < n; ++i)
413 >                update();
414 >            return n;
415 >        }
416 >
417 >        /**
418 >         * Chooses a breeder, exchanges strand with another subpop, and
419 >         * crosses them to create new chromosome to replace a chosen
420 >         * dyer.
421 >         */
422 >        void update() throws InterruptedException {
423 >            int b = chooseBreeder();
424 >            int d = chooseDyer(b);
425 >            Chromosome breeder = chromosomes[b];
426 >            Chromosome child = chromosomes[d];
427 >            chooseStrand(breeder);
428 >            strand = exchanger.exchange(strand);
429 >            cross(breeder, child);
430 >            fixOrder(child, d);
431 >        }
432 >
433 >        /**
434 >         * Chooses a breeder, with exponentially decreasing probability
435 >         * starting at best.
436 >         * @return index of selected breeder
437 >         */
438 >        int chooseBreeder() {
439 >            int mask = (subpopSize >>> BREEDER_DECAY) - 1;
440 >            int b = 0;
441 >            while ((rng.next() & mask) != mask) {
442 >                if (++b >= subpopSize)
443 >                    b = 0;
444 >            }
445 >            return b;
446 >        }
447 >
448 >        /**
449 >         * Chooses a chromosome that will be replaced, with
450 >         * exponentially decreasing probability starting at
451 >         * worst, ignoring the excluded index.
452 >         * @param exclude index to ignore; use -1 to not exclude any
453 >         * @return index of selected dyer
454 >         */
455 >        int chooseDyer(int exclude) {
456 >            int mask = (subpopSize >>> DYER_DECAY)  - 1;
457 >            int d = subpopSize - 1;
458 >            while (d == exclude || (rng.next() & mask) != mask) {
459 >                if (--d < 0)
460 >                    d = subpopSize - 1;
461 >            }
462 >            return d;
463 >        }
464 >
465 >        /**
466 >         * Select a random strand of b's.
467 >         * @param breeder the breeder
468 >         */
469 >        void chooseStrand(Chromosome breeder) {
470 >            int[] bs = breeder.alleles;
471 >            int length = bs.length;
472 >            int strandLength = MIN_STRAND_LENGTH;
473 >            while (strandLength < length &&
474 >                   (rng.next() & RANDOM_STRAND_MASK) != RANDOM_STRAND_MASK)
475 >                strandLength++;
476 >            strand.strandLength = strandLength;
477 >            int[] ss = strand.alleles;
478 >            int k = (rng.next() & 0x7FFFFFFF) % length;
479 >            for (int i = 0; i < strandLength; ++i) {
480 >                ss[i] = bs[k];
481 >                if (++k >= length) k = 0;
482 >            }
483 >        }
484 >
485 >        /**
486 >         * Copies current strand to start of c's, and then appends all
487 >         * remaining b's that aren't in the strand.
488 >         * @param breeder the breeder
489 >         * @param child the child
490 >         */
491 >        void cross(Chromosome breeder, Chromosome child) {
492 >            for (int k = 0; k < inTour.length; ++k) // clear bitset
493 >                inTour[k] = 0;
494 >
495 >            // Copy current strand to c
496 >            int[] cs = child.alleles;
497 >            int ssize = strand.strandLength;
498 >            int[] ss = strand.alleles;
499 >            int i;
500 >            for (i = 0; i < ssize; ++i) {
501 >                int x = ss[i];
502 >                cs[i] = x;
503 >                inTour[x >>> 5] |= 1 << (x & 31); // record in bit set
504 >            }
505 >
506 >            // Find index of matching origin in b
507 >            int first = cs[0];
508 >            int j = 0;
509 >            int[] bs = breeder.alleles;
510 >            while (bs[j] != first)
511 >                ++j;
512 >
513 >            // Append remaining b's that aren't already in tour
514 >            while (i < cs.length) {
515 >                if (++j >= bs.length) j = 0;
516 >                int x = bs[j];
517 >                if ((inTour[x >>> 5] & (1 << (x & 31))) == 0)
518 >                    cs[i++] = x;
519 >            }
520 >
521 >        }
522 >
523 >        /**
524 >         * Fixes the sort order of a changed Chromosome c at position k.
525 >         * @param c the chromosome
526 >         * @param k the index
527 >         */
528 >        void fixOrder(Chromosome c, int k) {
529 >            Chromosome[] cs = chromosomes;
530 >            int oldFitness = c.fitness;
531 >            c.recalcFitness();
532 >            int newFitness = c.fitness;
533 >            if (newFitness < oldFitness) {
534 >                int j = k;
535 >                int p = j - 1;
536 >                while (p >= 0 && cs[p].fitness > newFitness) {
537 >                    cs[j] = cs[p];
538 >                    j = p--;
539                  }
540 <                if (k != 0)
541 <                    children[i] = offspring;
542 <                else {
543 <                    // No peers, so we mate with ourselves
544 <                    for ( ; i < nSub - 1; i += 2) {
545 <                        int cur = firstSub + i;
546 <                        Chromosome bro = individuals[breeders[cur]];
346 <                        Chromosome sis = individuals[breeders[cur + 1]];
347 <                        
348 <                        children[i] = bro.breedWith(sis, matings, random);
349 <                        children[i+1] = sis.breedWith(bro, matings, random);
350 <                    }
351 <                    
352 <                    // Not even a sibling, so we go to asexual reproduction
353 <                    if (i < nSub)
354 <                        children[i] = individuals[breeders[firstSub + 1]];
355 <                    break;
540 >                cs[j] = c;
541 >            } else if (newFitness > oldFitness) {
542 >                int j = k;
543 >                int n = j + 1;
544 >                while (n < cs.length && cs[n].fitness < newFitness) {
545 >                    cs[j] = cs[n];
546 >                    j = n++;
547                  }
548 <
548 >                cs[j] = c;
549              }
359
360            awaitMatingBarrier(tid);
361
362            // Kill off dead guys
363            for (int i = 0; i < nSub; i++) {
364                individuals[dyers[firstSub + 1]] = children[i];
365            }
366
367            generationBarrier.await();
550          }
551      }
552  
553 <    static final class Chromosome {
554 <        private final CitySet cities;
555 <        private final int[] alleles;  
556 <        private final int length;    
557 <        public  double fitness; // immutable after publication
558 <        
559 <        // Basic constructor - gets randomized genetic code
560 <        Chromosome(CitySet cities, int length,
561 <                   LoopHelpers.SimpleRandom random) {
562 <            this.length = length;
563 <            this.cities = cities;
564 <            // Initialize alleles to a random shuffle
553 >    /**
554 >     * A Chromosome is a candidate TSP tour.
555 >     */
556 >    static final class Chromosome implements Comparable {
557 >        /** Index of cities in tour order */
558 >        final int[] alleles;
559 >        /** Total tour length */
560 >        int fitness;
561 >
562 >        /**
563 >         * Initializes to random tour.
564 >         */
565 >        Chromosome(int length, RNG random) {
566              alleles = new int[length];
567              for (int i = 0; i < length; i++)
568                  alleles[i] = i;
569              for (int i = length - 1; i > 0; i--) {
570 +                int idx = (random.next() & 0x7FFFFFFF) % alleles.length;
571                  int tmp = alleles[i];
388                int idx = random.next() % length;
572                  alleles[i] = alleles[idx];
573                  alleles[idx] = tmp;
574              }
575              recalcFitness();
576          }
577 <        
578 <        // Copy constructor - clones parent's genetic code
579 <        Chromosome(Chromosome clone) {
580 <            length = clone.length;
581 <            cities = clone.cities;
399 <            fitness = clone.fitness;
400 <            alleles = new int[length];
401 <            System.arraycopy(clone.alleles, 0, alleles, 0, length);
402 <        }
403 <        
404 <        int getAllele(int offset) {
405 <            return alleles[offset % length];
406 <        }
407 <        void setAllele(int offset, int v) {
408 <            alleles[offset % length] = v;
577 >
578 >        public int compareTo(Object x) { // to enable sorting
579 >            int xf = ((Chromosome) x).fitness;
580 >            int f = fitness;
581 >            return ((f == xf) ? 0 :((f < xf) ? -1 : 1));
582          }
583 <        
583 >
584          void recalcFitness() {
585 <            fitness = cities.distanceBetween(alleles[0], alleles[length-1]);
586 <            for (int i = 1; i < length; i++) {
587 <                fitness += cities.distanceBetween(alleles[i-1], alleles[i]);
588 <            }
589 <        }
590 <        
591 <        Chromosome breedWith(Chromosome partner, int n,
592 <                             LoopHelpers.SimpleRandom random) {
593 <            Chromosome offspring = new Chromosome(this);
594 <            for (int i = 0; i < n; i++)
595 <                offspring = offspring.reproduceWith(partner, random);
596 <            return offspring;
597 <        }
598 <        
599 <        Chromosome reproduceWith(Chromosome other,
600 <                                 LoopHelpers.SimpleRandom random) {
601 <            Chromosome child = new Chromosome(this);
602 <            int coStart = random.next() % length;
603 <            int coLen = 3;
604 <            while (1 == (random.next() & 1) && (coLen < length))
605 <                coLen++;
606 <            int cPos, pPos;
607 <            
608 <            int join = other.getAllele(coStart);
609 <            child.alleles[0] = join;
610 <            
611 <            for (pPos = 0; alleles[pPos] != join; pPos++)
612 <                ;
613 <            
614 <            for (cPos = 1; cPos < coLen; cPos++)
615 <                child.setAllele(cPos, other.getAllele(coStart + cPos));
616 <            
617 <            for (int i = 0; i < length; i++) {
618 <                boolean found = false;
619 <                int allele = getAllele(pPos++);
620 <                for (int j = 0; j < coLen; j++) {
621 <                    if (found = (child.getAllele(j) == allele))
622 <                        break;
623 <                }
451 <                if (!found)
452 <                    child.setAllele(cPos++, allele);
453 <            }
454 <            
455 <            child.recalcFitness();
456 <            return child;
585 >            int[] a = alleles;
586 >            int len = a.length;
587 >            int p = a[0];
588 >            long f = cities.distanceBetween(a[len-1], p);
589 >            for (int i = 1; i < len; i++) {
590 >                int n = a[i];
591 >                f += cities.distanceBetween(p, n);
592 >                p = n;
593 >            }
594 >            fitness = (int) (f / len);
595 >        }
596 >
597 >        /**
598 >         * Returns tour length for points scaled in [0, 1).
599 >         */
600 >        double unitTourLength() {
601 >            int[] a = alleles;
602 >            int len = a.length;
603 >            int p = a[0];
604 >            double f = cities.unitDistanceBetween(a[len-1], p);
605 >            for (int i = 1; i < len; i++) {
606 >                int n = a[i];
607 >                f += cities.unitDistanceBetween(p, n);
608 >                p = n;
609 >            }
610 >            return f;
611 >        }
612 >
613 >        /**
614 >         * Checks that this tour visits each city.
615 >         */
616 >        void validate() {
617 >            int len = alleles.length;
618 >            boolean[] used = new boolean[len];
619 >            for (int i = 0; i < len; ++i)
620 >                used[alleles[i]] = true;
621 >            for (int i = 0; i < len; ++i)
622 >                if (!used[i])
623 >                    throw new Error("Bad tour");
624          }
625 <        
625 >
626      }
627 +
628      /**
629 <     * A collection of (x,y) points that represent cities
629 >     * A Strand is a random sub-sequence of a Chromosome.  Each subpop
630 >     * creates only one strand, and then trades it with others,
631 >     * refilling it on each iteration.
632 >     */
633 >    static final class Strand {
634 >        final int[] alleles;
635 >        int strandLength;
636 >        Strand(int length) { alleles = new int[length]; }
637 >    }
638 >
639 >    /**
640 >     * A collection of (x,y) points that represent cities.
641       */
642      static final class CitySet {
643 <        final int XMAX = 1000;
465 <        final int YMAX = 1000;
643 >
644          final int length;
645 <        final int xPts[];
646 <        final int yPts[];
647 <        final double distances[];
648 <        
649 <        CitySet(int n, LoopHelpers.SimpleRandom random) {
645 >        final int[] xPts;
646 >        final int[] yPts;
647 >        final int[][] distances;
648 >
649 >        CitySet(int n) {
650              this.length = n;
651 <            xPts = new int[n];
652 <            yPts = new int [n];
651 >            this.xPts = new int[n];
652 >            this.yPts = new int[n];
653 >            this.distances = new int[n][n];
654 >
655 >            RNG random = new RNG();
656              for (int i = 0; i < n; i++) {
657 <                xPts[i] = random.next() % XMAX;
658 <                yPts[i] = random.next() % YMAX;
657 >                xPts[i] = (random.next() & 0x7FFFFFFF);
658 >                yPts[i] = (random.next() & 0x7FFFFFFF);
659              }
660 <            distances = new double[n * n];
660 >
661              for (int i = 0; i < n; i++) {
662                  for (int j = 0; j < n; j++) {
663 <                    double dX = (double)(xPts[i] - xPts[j]);
664 <                    double dY = (double)(yPts[i] - yPts[j]);
665 <                    distances[i + j * n] = Math.hypot(dX, dY);
663 >                    double dx = (double) xPts[i] - (double) xPts[j];
664 >                    double dy = (double) yPts[i] - (double) yPts[j];
665 >                    double dd = Math.hypot(dx, dy) / 2.0;
666 >                    long ld = Math.round(dd);
667 >                    distances[i][j] = (ld >= Integer.MAX_VALUE) ?
668 >                        Integer.MAX_VALUE : (int) ld;
669                  }
670              }
671          }
672 <        
673 <        // Retrieve the cached distance between a pair of cities
674 <        double distanceBetween(int idx1, int idx2) {
675 <            return distances[idx1 + idx2 * length];
672 >
673 >        /**
674 >         * Returns the cached distance between a pair of cities.
675 >         */
676 >        int distanceBetween(int i, int j) {
677 >            return distances[i][j];
678 >        }
679 >
680 >        // Scale ints to doubles in [0,1)
681 >        static final double PSCALE = (double) 0x80000000L;
682 >
683 >        /**
684 >         * Returns distance for points scaled in [0,1). This simplifies
685 >         * checking results.  The expected optimal TSP for random
686 >         * points is believed to be around 0.76 * sqrt(N). For papers
687 >         * discussing this, see
688 >         * http://www.densis.fee.unicamp.br/~moscato/TSPBIB_home.html
689 >         */
690 >        double unitDistanceBetween(int i, int j) {
691 >            double dx = ((double) xPts[i] - (double) xPts[j]) / PSCALE;
692 >            double dy = ((double) yPts[i] - (double) yPts[j]) / PSCALE;
693 >            return Math.hypot(dx, dy);
694          }
695 +
696      }
697  
698 <    static final class Runner implements Runnable {
699 <        final Population pop;
700 <        final CyclicBarrier b;
701 <        final int nGen;
702 <        final int tid;
703 <        static final boolean verbose = false;
704 <        
705 <        Runner(int tid, Population p, CyclicBarrier b, int n) {
706 <            this.tid = tid;
707 <            this.pop = p;
708 <            this.b = b;
709 <            this.nGen = n;
698 >    /**
699 >     * Cheap XorShift random number generator
700 >     */
701 >    static final class RNG {
702 >        /** Seed generator for XorShift RNGs */
703 >        static final Random seedGenerator = new Random();
704 >
705 >        int seed;
706 >        RNG(int seed) { this.seed = seed; }
707 >        RNG()         { this.seed = seedGenerator.nextInt() | 1; }
708 >
709 >        int next() {
710 >            int x = seed;
711 >            x ^= x << 6;
712 >            x ^= x >>> 21;
713 >            x ^= x << 7;
714 >            seed = x;
715 >            return x;
716          }
717 <        
717 >    }
718 >
719 >    static final class ProgressMonitor extends Thread {
720 >        final Population pop;
721 >        ProgressMonitor(Population p) { pop = p; }
722          public void run() {
723 +            double time = 0;
724              try {
725 <                b.await();
726 <                for (int i = 0; i < nGen; i++) {
727 <                    if (verbose && 0 == tid && 0 == i % 1000) {
728 <                        System.out.print("Gen " + i + " fitness:");
515 <                        System.out.print(" best=" + (int)pop.bestFitness());
516 <                        System.out.println("; avg=" + (int)pop.averageFitness());
517 <                    }
518 <                    int matings = (((nGen - i) * 2) / (nGen)) + 1;
519 <                    pop.nextGeneration(tid, matings);
725 >                while (!Thread.interrupted()) {
726 >                    sleep(SNAPSHOT_RATE);
727 >                    time += SNAPSHOT_RATE;
728 >                    pop.printSnapshot(time / 1000.0);
729                  }
730 <                b.await();
522 <            }
523 <            catch (InterruptedException e) {
524 <                e.printStackTrace(System.out);
525 <                System.exit(0);
526 <            }
527 <            catch (BrokenBarrierException e) {
528 <                e.printStackTrace(System.out);
529 <                System.exit(0);
530 <            }
730 >            } catch (InterruptedException ie) {}
731          }
732      }
733   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines