[GRASS-SVN] r51622 - grass-addons/grass7/vector/v.net.salesman.opt
svn_grass at osgeo.org
svn_grass at osgeo.org
Mon May 14 10:37:59 EDT 2012
Author: mmetz
Date: 2012-05-14 07:37:56 -0700 (Mon, 14 May 2012)
New Revision: 51622
Added:
grass-addons/grass7/vector/v.net.salesman.opt/evolution.c
grass-addons/grass7/vector/v.net.salesman.opt/local_proto.h
grass-addons/grass7/vector/v.net.salesman.opt/optimize.c
grass-addons/grass7/vector/v.net.salesman.opt/tour.c
grass-addons/grass7/vector/v.net.salesman.opt/v.net.salesman.opt.html
Removed:
grass-addons/grass7/vector/v.net.salesman.opt/v.net.salesman.html
Modified:
grass-addons/grass7/vector/v.net.salesman.opt/Makefile
grass-addons/grass7/vector/v.net.salesman.opt/main.c
Log:
v.net.salesman with tour optimization
Modified: grass-addons/grass7/vector/v.net.salesman.opt/Makefile
===================================================================
--- grass-addons/grass7/vector/v.net.salesman.opt/Makefile 2012-05-14 14:31:54 UTC (rev 51621)
+++ grass-addons/grass7/vector/v.net.salesman.opt/Makefile 2012-05-14 14:37:56 UTC (rev 51622)
@@ -1,7 +1,7 @@
MODULE_TOPDIR = ../..
-PGM = v.net.salesman
+PGM = v.net.salesman.opt
LIBES = $(VECTORLIB) $(GISLIB)
DEPENDENCIES = $(VECTORDEP) $(GISDEP)
Added: grass-addons/grass7/vector/v.net.salesman.opt/evolution.c
===================================================================
--- grass-addons/grass7/vector/v.net.salesman.opt/evolution.c (rev 0)
+++ grass-addons/grass7/vector/v.net.salesman.opt/evolution.c 2012-05-14 14:37:56 UTC (rev 51622)
@@ -0,0 +1,1010 @@
+#include <stdlib.h>
+#include <string.h>
+#include <grass/gis.h>
+#include <grass/vector.h>
+#include <grass/glocale.h>
+#include "local_proto.h"
+
+struct tsp_tour
+{
+ int *cycle; /* order of cities */
+ double cost; /* total cost */
+ int used;
+ int opt_done;
+};
+
+struct tour_cost
+{
+ int t;
+ double cost;
+};
+
+/* common variables */
+int *tcused, *tmpcycle;
+
+/* helper functions */
+
+int init_tour(struct tsp_tour *tour);
+int init_ga_1(struct tsp_tour *tour, int ntours);
+int init_ga_2(struct tsp_tour *tour, int ntours);
+int ga_rem_dupl_tours(struct tsp_tour *tour, int ntours);
+int ga_recombine(struct tsp_tour *tour, int p1, int p2, int c);
+int ga_recombine2(struct tsp_tour *tour, int p1, int p2, int c);
+int cmp_tour_cost(const void *pa, const void *pb);
+
+/* genetic algorithm:
+ * initialization
+ * loop over
+ * 1. natural selection
+ * 2. recombination
+ * 3. mutation
+ * until there is no better solution */
+
+int ga_opt(int ntours, int nelim, int nopt, int ngen, int *best_cycle)
+{
+ int i, j, k, t;
+ double best_cost, worst_cost;
+ int no_better, gen, best_tour, new_child, ncidx;
+ int ntours_left, nparents;
+ int nunused, *tunused;
+ struct tsp_tour *tour, *tp;
+ struct tour_cost *tc, *td;
+ int init_method = 2;
+
+ int max_chain_length, chain_length;
+ int optiter, success;
+
+ tour = G_malloc(ntours * sizeof(struct tsp_tour));
+ tc = G_malloc(ntours * sizeof(struct tour_cost));
+ td = G_malloc(ntours * sizeof(struct tour_cost));
+ tunused = G_malloc(ntours * sizeof(int));
+ tcused = G_malloc(ncities * sizeof(int));
+ tmpcycle = G_malloc(ncities * sizeof(int));
+
+ /* initialize */
+ G_verbose_message(_("Generating %d tours"), ntours);
+ for (i = 0; i < ntours; i++) {
+ tp = &(tour[i]);
+ init_tour(tp);
+ tunused[i] = -1;
+ }
+ if (init_method == 1) {
+ if (init_ga_1(tour, ntours) < ntours)
+ G_fatal_error(_("Method 1 failed to create %d tours"), ntours);
+ }
+ else if (init_method == 2) {
+ if (init_ga_2(tour, ntours) < ntours)
+ G_fatal_error(_("Method 2 failed to create %d tours"), ntours);
+ }
+ else {
+ G_fatal_error(_("Method %d is not implemented"), init_method);
+ }
+
+ /* tour costs */
+ best_cost = -1;
+ best_tour = -1;
+ for (i = 0; i < ntours; i++) {
+ tp = &(tour[i]);
+ tp->cycle[ncities] = tp->cycle[0];
+ tp->cost = 0;
+ for (j = 0; j < ncities; j++)
+ tp->cost += cost_cache[tp->cycle[j]][tp->cycle[j + 1]];
+
+ tc[i].cost = tp->cost;
+ tc[i].t = i;
+
+ if (best_cost < 0) {
+ best_cost = tp->cost;
+ best_tour = i;
+ }
+ else if (best_cost > tp->cost) {
+ best_cost = tp->cost;
+ best_tour = i;
+ }
+ }
+
+ no_better = gen = 0;
+
+ while (!no_better && gen < ngen) {
+
+ G_message(_("%d. Generation"), gen + 1);
+
+ /* 1. natural selection */
+ G_verbose_message(_("Natural Selection"));
+
+ /* remove duplicate tours */
+ ga_rem_dupl_tours(tour, ntours);
+
+ ntours_left = nunused = 0;
+ for (i = 0; i < ntours; i++) {
+ if (!tour[i].used) {
+ tour[i].opt_done = 0;
+ tunused[nunused++] = i;
+ continue;
+ }
+
+ tc[ntours_left].cost = tour[i].cost;
+ tc[ntours_left].t = i;
+ ntours_left++;
+ }
+ if (nunused)
+ G_debug(1, "%d tours: %d duplicates removed", ntours, nunused);
+
+ /* sort tours ascending by cost */
+ qsort((void *)tc, ntours_left, sizeof(struct tour_cost), cmp_tour_cost);
+
+ /* sort by cost difference to previous tour */
+ for (t = 1; t < ntours_left; t++) {
+ td[t - 1].cost = tour[tc[t].t].cost - tour[tc[t - 1].t].cost;
+ td[t - 1].t = tc[t].t;
+ }
+ qsort((void *)td, ntours_left - 1, sizeof(struct tour_cost), cmp_tour_cost);
+ G_debug(3, "Smallest difference: %.3f", td[0].cost);
+
+ /* remove nelim tours */
+ k = 0;
+ while (ntours_left > ntours - nelim || td[k].cost < 0.0001) {
+ tour[td[k].t].used = 0;
+ tour[td[k].t].opt_done = 0;
+ tunused[nunused++] = td[k].t;
+ ntours_left--;
+ k++;
+ }
+ G_debug(1, "%d tours: %d unused, %d left", ntours, nunused, ntours_left);
+ if (ntours_left < 2)
+ G_fatal_error(_("Diversity loss"));
+
+ nunused = 0;
+ for (i = 0, j = 0; i < ntours; i++) {
+ if (!tour[i].used) {
+ tunused[nunused++] = i;
+ continue;
+ }
+ tc[j].cost = tour[i].cost;
+ tc[j++].t = i;
+ }
+ if (j != ntours_left)
+ G_fatal_error(_("Wrong number of remaining tours"));
+
+ /* sort tours ascending by cost */
+ qsort((void *)tc, ntours_left, sizeof(struct tour_cost), cmp_tour_cost);
+
+ /* 2. recombination */
+ G_verbose_message(_("Recombination"));
+ nparents = ntours_left;
+
+ new_child = -1;
+ ncidx = nunused - 1;
+ for (i = 0; i < nparents; i++) {
+ for (j = i + 1; j < nparents; j++) {
+
+ if (ntours_left < ntours) {
+ if (ncidx < 0)
+ G_fatal_error(_("1 ncidx too small"));
+
+ if (tunused[ncidx] < 0)
+ G_fatal_error(_("1 tunused too small"));
+
+ if (tunused[ncidx] >= ntours)
+ G_fatal_error(_("1 tunused too large"));
+
+ new_child = tunused[ncidx];
+ }
+ else {
+ worst_cost = -1;
+ for (k = 0; k < nunused; k++) {
+ if (worst_cost < tour[tunused[k]].cost) {
+ worst_cost = tour[tunused[k]].cost;
+ new_child = tunused[k];
+ }
+ }
+ }
+
+ if (ga_recombine(tour, tc[i].t, tc[j].t, new_child)) {
+ if (ntours_left < ntours) {
+ ncidx--;
+ ntours_left++;
+ }
+ }
+
+ if (ntours_left < ntours) {
+ if (ncidx < 0)
+ G_fatal_error(_("2 ncidx too small"));
+
+ if (tunused[ncidx] < 0)
+ G_fatal_error(_("2 tunused too small"));
+
+ if (tunused[ncidx] >= ntours)
+ G_fatal_error(_("2 tunused too large"));
+
+ new_child = tunused[ncidx];
+ }
+ else {
+ worst_cost = -1;
+ for (k = 0; k < nunused; k++) {
+ if (worst_cost < tour[tunused[k]].cost) {
+ worst_cost = tour[tunused[k]].cost;
+ new_child = tunused[k];
+ }
+ }
+ }
+
+ if (ga_recombine(tour, tc[j].t, tc[i].t, new_child)) {
+ if (ntours_left < ntours) {
+ ncidx--;
+ ntours_left++;
+ }
+ }
+ }
+ }
+ if (ntours_left < ntours) {
+ G_debug(1, "Failed to generate enough new tours (old: %d, new: %d)", nparents, ntours_left);
+ }
+
+ /* 3. optimization */
+ G_verbose_message(_("Mutation"));
+ for (i = 0; i < ntours; i++) {
+ success = 1;
+
+ if (!tour[i].used)
+ continue;
+
+ /*
+ if (tour[i].opt_done)
+ continue;
+ */
+
+ /* this can break out of the current local minimum */
+ optimize_tour_chains(4, 2, 2, tour[i].cycle, tcused, ncities, 0, 0, 1);
+
+ for (j = 0; j < ncities; j++) {
+ optimize_nbrs(j, ncities, tour[i].cycle);
+ tcused[j] = 1;
+ }
+
+ /* this would approximate a local minimum,
+ * but is slow and will be done later on with the best tour anyway
+ * disabled */
+#if 0
+ while (success)
+ success = optimize_tour_chains(4, 4, 1, tour[i].cycle, tcused, ncities, 0, 1, 1);
+#endif
+
+ tp = &(tour[i]);
+ tp->cycle[ncities] = tp->cycle[0];
+ tp->cost = 0;
+ for (j = 0; j < ncities; j++)
+ tp->cost += cost_cache[tp->cycle[j]][tp->cycle[j + 1]];
+
+ tour[i].opt_done = 1;
+ }
+ /* debug */
+ for (t = 0; t < ntours; t++) {
+ if (!tour[i].used)
+ continue;
+
+ for (i = 0; i < ncities; i++)
+ tcused[i] = 0;
+ for (i = 0; i < ncities; i++) {
+ if (tcused[tour[t].cycle[i]])
+ G_fatal_error(_("Duplicate city"));
+ tcused[tour[t].cycle[i]] = 1;
+ }
+ }
+
+ /* tour costs */
+ ntours_left = 0;
+ for (i = 0; i < ntours; i++) {
+ if (!tour[i].used)
+ continue;
+
+ tp = &(tour[i]);
+
+ tc[ntours_left].cost = tp->cost;
+ tc[ntours_left].t = i;
+ ntours_left++;
+ }
+
+ /* sort tours ascending by cost */
+ qsort((void *)tc, ntours_left, sizeof(struct tour_cost), cmp_tour_cost);
+
+ for (i = 0; i < ntours_left; i++) {
+ if (!tour[tc[i].t].used)
+ continue;
+ G_debug(3, "%d. tour, cost %.3f", i + 1, tc[i].cost);
+ }
+
+ if (best_cost > tour[tc[0].t].cost + 0.0001) {
+ no_better = 0;
+ }
+ else
+ no_better++;
+
+ best_cost = tour[tc[0].t].cost;
+ best_tour = tc[0].t;
+
+ G_verbose_message(_("Best tour: %d, best cost: %.3f"), best_tour, best_cost);
+
+ if (ntours_left == nparents)
+ break;
+
+ gen++;
+ }
+
+ for (i = 0; i < ncities; i++)
+ best_cycle[i] = tour[best_tour].cycle[i];
+ best_cycle[ncities] = best_cycle[0];
+
+ /* brute force optimization of the best tour */
+ G_message(_("Optimizing the best tour"));
+
+ for (j = 0; j < ncities; j++) {
+ optimize_nbrs(j, ncities, best_cycle);
+ }
+
+ optiter = 1;
+ max_chain_length = MAX_CHAIN_LENGTH < ncities / 2 ? MAX_CHAIN_LENGTH : ncities / 2;
+ for (chain_length = max_chain_length; chain_length > 0; chain_length--) {
+ success = 0;
+
+ G_message(_("%d. iteration..."), optiter++);
+
+ if (max_chain_length >= 5)
+ optimize_tour_chains(5, 5, 1,
+ best_cycle, tcused, ncities, 1, 1, 0);
+ else
+ optimize_tour_chains(max_chain_length, max_chain_length, 1,
+ best_cycle, tcused, ncities, 1, 1, 0);
+
+ success = optimize_tour_chains(max_chain_length, 1, 1, best_cycle, tcused, ncities, 1, 1, 1);
+
+ if (!success)
+ break;
+ }
+
+ return gen;
+}
+
+
+int init_tour(struct tsp_tour *tour)
+{
+ tour->cycle = (int *)G_malloc((ncities + 1) * sizeof(int));
+
+ if (tour->cycle == NULL)
+ G_fatal_error(_("Out of memory"));
+
+ tour->cost = -1;
+ tour->used = 0;
+ tour->opt_done = 0;
+
+ return 1;
+}
+
+
+int init_ga_1(struct tsp_tour *tour, int ntours)
+{
+ int t, i;
+ struct tsp_tour *tp;
+ double cost, tmpcost;
+ int start1, start2;
+ int nxt1, nxt2;
+ int tnc;
+
+ start1 = start2 = -1;
+ nxt1 = nxt2 = 0;
+ t = 0;
+
+ while (t < ntours) {
+ tp = &(tour[t]);
+
+ tnc = 0;
+ for (i = 0; i < ncities; i++)
+ tcused[i] = 0;
+
+ if (t == 0) {
+ /* the first tour is the standard tour */
+
+ cost = start1 = -1;
+ for (i = 0; i < ncities; i++) {
+ tmpcost = costs[i][ncities - 2].cost;
+ if (tmpcost > cost) {
+ cost = tmpcost;
+ start1 = i;
+ }
+ }
+ start2 = costs[start1][ncities - 2].city;
+
+ /* add these 2 cities to array */
+ add_city(start1, -1, &tnc, tp->cycle, tcused);
+ add_city(start2, 0, &tnc, tp->cycle, tcused);
+ }
+ else {
+ if (nxt2 <= nxt1)
+ nxt2 = nxt1 + 1;
+ if (nxt2 >= ncities) {
+ nxt1++;
+ nxt2 = nxt1 + 1;
+ }
+ if (nxt1 == ncities - 1)
+ return t;
+ if (nxt2 >= ncities)
+ return t;
+ if ((nxt1 == start1 && nxt2 == start2) || (nxt1 == start2 && nxt2 == start1)) {
+ nxt2++;
+ if (nxt2 >= ncities) {
+ nxt1++;
+ nxt2 = nxt1 + 1;
+ }
+ if (nxt1 == ncities - 1)
+ return t;
+ if (nxt2 >= ncities)
+ return t;
+ }
+ /* add these 2 cities to array */
+ tnc = 0;
+ add_city(nxt1, -1, &tnc, tp->cycle, tcused);
+ add_city(nxt2, 0, &tnc, tp->cycle, tcused);
+ }
+
+ build_tour(tp->cycle, tcused, &tnc, 2, 1);
+
+ /* debug */
+ for (i = 0; i < ncities; i++)
+ tcused[i] = 0;
+ for (i = 0; i < ncities; i++) {
+ if (tcused[tp->cycle[i]])
+ G_fatal_error(_("Duplicate city"));
+ tcused[tp->cycle[i]] = 1;
+ }
+
+ tp->used = 1;
+ t++;
+ }
+
+ return t;
+}
+
+int init_ga_2(struct tsp_tour *tour, int ntours)
+{
+ int t, i;
+ int city, city1;
+ int tnc;
+ struct tsp_tour *tp;
+
+ /* ntours must be <= ncities */
+ if (ntours > ncities)
+ G_fatal_error(_("The number of tours (%d) must not be larger than the number of cities (%d)"),
+ ntours, ncities);
+
+ /* use each city as start
+ * go to nearest unused city
+ * from this city
+ * go to nearest unused city
+ * etc until all cities are in the cycle */
+ for (t = 0; t < ntours; t++) {
+ tp = &(tour[t]);
+
+ tnc = 0;
+ for (i = 0; i < ncities; i++)
+ tcused[i] = 0;
+
+ city = t;
+ tp->cycle[0] = city;
+ tcused[city] = 1;
+ tnc++;
+
+ while (tnc < ncities) {
+ /* get nearest unused city */
+ city1 = -1;
+ for (i = 0; i < ncities - 1; i++) {
+ if (tcused[costs[city][i].city])
+ continue;
+ city1 = costs[city][i].city;
+ break;
+ }
+ if (city1 == -1)
+ G_fatal_error(_("No unused city left"));
+ tp->cycle[tnc] = city1;
+ city = city1;
+ tcused[city] = 1;
+ tnc++;
+ }
+ tp->used = 1;
+
+ /* debug */
+ for (i = 0; i < ncities; i++)
+ tcused[i] = 0;
+ for (i = 0; i < ncities; i++) {
+ if (tcused[tp->cycle[i]])
+ G_fatal_error(_("Duplicate city"));
+ tcused[tp->cycle[i]] = 1;
+ }
+ }
+
+ return ntours;
+}
+
+
+int ga_rem_dupl_tours(struct tsp_tour *tour, int ntours)
+{
+ int i, j, k;
+ int n_left, start1, start2, nxt1, nxt2;
+ int tequal;
+
+ /* remove duplicate tours */
+ n_left = ntours;
+
+ for (i = 0; i < ntours; i++) {
+
+ if (!tour[i].used)
+ continue;
+
+ start1 = 0;
+ for (j = i + 1; j < ntours; j++) {
+ tequal = 1;
+
+ if (!tour[j].used)
+ continue;
+
+ start2 = -1;
+ for (k = 0; k < ncities; k++) {
+ if (tour[j].cycle[k] == tour[i].cycle[start1]) {
+ start2 = k;
+ break;
+ }
+ }
+ if (start2 == -1)
+ G_fatal_error(_("start2 not found"));
+
+ if (tour[i].cycle[start1] != tour[j].cycle[start2])
+ G_fatal_error(_("No common start city"));
+
+ /* move this up */
+ if (fabs(tour[i].cost - tour[j].cost) > 0.01)
+ continue;
+
+ for (k = 0; k < ncities; k++) {
+ nxt1 = wrap_into(start1 + k, ncities);
+ nxt2 = wrap_into(start2 + k, ncities);
+ if (tour[i].cycle[nxt1] != tour[j].cycle[nxt2])
+ tequal = 0;
+ break;
+ }
+ if (tequal) {
+ /* remove the tour */
+ G_debug(1, "Tours %d and %d are identical", i, j);
+ tour[j].used = 0;
+ n_left--;
+ }
+ }
+ }
+
+ return n_left;
+}
+
+/* recombination similar to
+ * http://www.gcd.org/sengoku/docs/arob98.pdf
+ * here: not the longest (most nodes) subtour, but the shortest (costs) total tour */
+int ga_subtour(struct tsp_tour *tour, int p1, int p2, int c,
+ int next_i, int *start_i)
+{
+ int i, j;
+ int prev1, prev2, nxt1, nxt2, startsub, endsub;
+ int cncyc;
+ int have_p1, have_p2;
+ int st_max, st_length, st_i, st_j;
+ struct tsp_tour *parent1, *parent2, *child;
+
+ G_debug(1, "ga_subtour()");
+
+ parent1 = &(tour[p1]);
+ parent2 = &(tour[p2]);
+ child = &(tour[c]);
+ cncyc = 0;
+
+ child->used = 0;
+
+ /* find subtour */
+ st_max = st_length = st_i = st_j = -1;
+
+ for (i = next_i; i < ncities; i++) {
+ for (j = 0; j < ncities; j++)
+ tcused[j] = 0;
+
+ prev1 = wrap_into(i - 1, ncities);
+ nxt1 = wrap_into(i + 1, ncities);
+
+ st_length = -1;
+ j = 0;
+
+ while (parent2->cycle[j] != parent1->cycle[i]) {
+ j++;
+ if (j >= ncities)
+ G_fatal_error(_("Missing city"));
+ }
+
+ if (parent2->cycle[j] != parent1->cycle[i])
+ G_fatal_error(_("City mismatch"));
+
+ prev2 = wrap_into(j - 1, ncities);
+ nxt2 = wrap_into(j + 1, ncities);
+
+ /* makes only sense if prev[1,2] are not equal or nxt[1,2] are not equal */
+ if (parent2->cycle[prev2] != parent1->cycle[nxt1] &&
+ (parent2->cycle[prev2] != parent1->cycle[prev1] &&
+ parent2->cycle[nxt2] != parent1->cycle[nxt1])) {
+
+ /* fake recombination */
+ have_p1 = have_p2 = 1;
+ tmpcycle[0] = parent1->cycle[i];
+ tcused[parent1->cycle[i]] = 1;
+ startsub = endsub = 0;
+ st_length = 1;
+
+ while (have_p1 || have_p2) {
+
+ if (have_p1 && tcused[parent1->cycle[nxt1]])
+ have_p1 = 0;
+
+ if (have_p1) {
+ endsub = wrap_into(endsub + 1, ncities);
+ tmpcycle[endsub] = parent1->cycle[nxt1];
+ tcused[parent1->cycle[nxt1]] = 1;
+ st_length++;
+
+ nxt1 = wrap_into(nxt1 + 1, ncities);
+ if (tcused[parent1->cycle[nxt1]])
+ have_p1 = 0;
+ }
+
+ if (have_p2 && tcused[parent2->cycle[prev2]])
+ have_p2 = 0;
+
+ if (have_p2) {
+ startsub = wrap_into(startsub - 1, ncities);
+ tmpcycle[startsub] = parent2->cycle[prev2];
+ tcused[parent2->cycle[prev2]] = 1;
+ st_length++;
+
+ prev2 = wrap_into(prev2 - 1, ncities);
+ if (tcused[parent2->cycle[prev2]])
+ have_p2 = 0;
+ }
+ }
+ }
+ if (st_max <= st_length) {
+ st_max = st_length;
+ st_i = i;
+ st_j = j;
+ }
+ if (st_max > 4)
+ break;
+ }
+ *start_i = i;
+
+ if (st_max <= 0)
+ return 0;
+
+ G_debug(3, "longest subtour %d at %d", st_max, st_i);
+
+ for (i = 0; i < ncities; i++)
+ tcused[i] = 0;
+
+ /* recombination for selected subtour */
+ i = st_i;
+ j = st_j;
+ nxt1 = wrap_into(i + 1, ncities);
+ prev2 = wrap_into(j - 1, ncities);
+ have_p1 = have_p2 = 1;
+ tmpcycle[0] = parent1->cycle[i];
+ tcused[parent1->cycle[i]] = 1;
+ startsub = endsub = 0;
+
+ while (have_p1 || have_p2) {
+
+ if (have_p1 && tcused[parent1->cycle[nxt1]])
+ have_p1 = 0;
+
+ if (have_p1) {
+ endsub = wrap_into(endsub + 1, ncities);
+ tmpcycle[endsub] = parent1->cycle[nxt1];
+ tcused[parent1->cycle[nxt1]] = 1;
+
+ nxt1 = wrap_into(nxt1 + 1, ncities);
+ if (tcused[parent1->cycle[nxt1]])
+ have_p1 = 0;
+ }
+
+ if (have_p2 && tcused[parent2->cycle[prev2]])
+ have_p2 = 0;
+
+ if (have_p2) {
+ startsub = wrap_into(startsub - 1, ncities);
+ tmpcycle[startsub] = parent2->cycle[prev2];
+ tcused[parent2->cycle[prev2]] = 1;
+
+ prev2 = wrap_into(prev2 - 1, ncities);
+ if (tcused[parent2->cycle[prev2]])
+ have_p2 = 0;
+ }
+ }
+ while (startsub != endsub) {
+ if (cncyc >= ncities)
+ G_fatal_error(_("Too many cities"));
+
+ child->cycle[cncyc++] = tmpcycle[startsub];
+ startsub = wrap_into(startsub + 1, ncities);
+ }
+ if (cncyc >= ncities)
+ G_fatal_error(_("Too many cities"));
+
+ if (startsub != endsub)
+ G_fatal_error(_("startsub != endsub"));
+
+ child->cycle[cncyc++] = tmpcycle[startsub];
+
+#if 0
+ /* this method is better, but the standard method */
+
+ build_tour(child->cycle, tcused, &cncyc, 2, 0);
+
+#else
+ /* this method is faster and true recombination */
+
+ /* add unused cities */
+ i = nxt1;
+ while (cncyc < ncities) {
+ int found = 0;
+
+ i = wrap_into(i + 1, ncities);
+
+ if (tcused[parent1->cycle[i]]) {
+ continue;
+ }
+ G_debug(3, "cncyc: %d", cncyc);
+
+ if (tcused[parent1->cycle[nxt1]]) {
+ child->cycle[cncyc++] = parent1->cycle[i];
+ tcused[parent1->cycle[i]] = 1;
+ continue;
+ }
+
+ prev1 = wrap_into(i - 1, ncities);
+ nxt1 = wrap_into(i + 1, ncities);
+
+ found = 0;
+ j = 0;
+
+ while (parent2->cycle[j] != parent1->cycle[i]) {
+ j++;
+ if (j >= ncities)
+ G_fatal_error(_("Missing city"));
+ }
+
+ if (parent2->cycle[j] != parent1->cycle[i])
+ G_fatal_error(_("City mismatch"));
+
+ prev2 = wrap_into(j - 1, ncities);
+ nxt2 = wrap_into(j + 1, ncities);
+
+ if (tcused[parent2->cycle[prev2]]) {
+ child->cycle[cncyc++] = parent1->cycle[i];
+ tcused[parent1->cycle[i]] = 1;
+ continue;
+ }
+
+ /* makes only sense if prev[1,2] are not equal or nxt[1,2] are not equal */
+ if (parent2->cycle[prev2] != parent1->cycle[nxt1] &&
+ (parent2->cycle[prev2] != parent1->cycle[prev1] &&
+ parent2->cycle[nxt2] != parent1->cycle[nxt1])) {
+
+ found = 1;
+
+ /* actual recombination */
+ have_p1 = have_p2 = 1;
+ tmpcycle[0] = parent1->cycle[i];
+ tcused[parent1->cycle[i]] = 1;
+ startsub = endsub = 0;
+
+ while (have_p1 || have_p2) {
+
+ if (have_p1 && tcused[parent1->cycle[nxt1]])
+ have_p1 = 0;
+
+ if (have_p1) {
+ endsub = wrap_into(endsub + 1, ncities);
+ tmpcycle[endsub] = parent1->cycle[nxt1];
+ tcused[parent1->cycle[nxt1]] = 1;
+
+ nxt1 = wrap_into(nxt1 + 1, ncities);
+ if (tcused[parent1->cycle[nxt1]])
+ have_p1 = 0;
+ }
+
+ if (have_p2 && tcused[parent2->cycle[prev2]])
+ have_p2 = 0;
+
+ if (have_p2) {
+ startsub = wrap_into(startsub - 1, ncities);
+ tmpcycle[startsub] = parent2->cycle[prev2];
+ tcused[parent2->cycle[prev2]] = 1;
+
+ prev2 = wrap_into(prev2 - 1, ncities);
+ if (tcused[parent2->cycle[prev2]])
+ have_p2 = 0;
+ }
+ }
+ }
+
+ if (found) {
+ /* add tmpcycle to child cycle */
+ G_debug(3, "add subtour, startsub %d, endsub %d", startsub, endsub);
+
+ while (startsub != endsub) {
+ if (cncyc >= ncities)
+ G_fatal_error(_("Too many cities"));
+
+ child->cycle[cncyc++] = tmpcycle[startsub];
+ startsub = wrap_into(startsub + 1, ncities);
+ }
+ if (cncyc >= ncities)
+ G_fatal_error(_("Too many cities"));
+
+ if (startsub != endsub)
+ G_fatal_error(_("startsub != endsub"));
+
+ child->cycle[cncyc++] = tmpcycle[startsub];
+ G_debug(3, "cncyc: %d", cncyc);
+ }
+ else {
+ if (cncyc >= ncities)
+ G_fatal_error(_("No subtour, too many cities"));
+
+ child->cycle[cncyc++] = parent1->cycle[i];
+ tcused[parent1->cycle[i]] = 1;
+ G_debug(3, "cncyc: %d", cncyc);
+ }
+ }
+#endif
+
+ child->cost = 0;
+ child->cycle[ncities] = child->cycle[0];
+ for (i = 0; i < ncities; i++)
+ child->cost +=cost_cache[child->cycle[i]][child->cycle[i + 1]];
+
+ child->used = 1;
+
+ return child->used;
+}
+
+/* Greedy Crossover as in java-tsp
+ * http://http://code.google.com/p/java-traveling-salesman/
+ * children are too expensive */
+int ga_recombine2(struct tsp_tour *tour, int p1, int p2, int c)
+{
+ int i;
+ int nxt1, nxt2, picked, last;
+ int cncyc;
+ struct tsp_tour *parent1, *parent2, *child;
+
+ parent1 = &(tour[p1]);
+ parent2 = &(tour[p2]);
+ child = &(tour[c]);
+
+ cncyc = 0;
+ nxt1 = 0;
+ nxt2 = 0;
+
+ for (i = 0; i < ncities; i++)
+ tcused[i] = 0;
+
+ picked = parent1->cycle[nxt1];
+ child->cycle[cncyc] = picked;
+ tcused[picked] = 1;
+ last = picked;
+ cncyc++;
+
+ while (cncyc < ncities) {
+
+ while (parent1->cycle[nxt1] != last)
+ nxt1 = wrap_into(nxt1 + 1, ncities);
+ nxt1 = wrap_into(nxt1 + 1, ncities);
+
+ while (parent2->cycle[nxt2] != last)
+ nxt2 = wrap_into(nxt2 + 1, ncities);
+ nxt2 = wrap_into(nxt2 + 1, ncities);
+
+ if (tcused[parent1->cycle[nxt1]] && tcused[parent2->cycle[nxt2]]) {
+ while (tcused[parent1->cycle[nxt1]])
+ nxt1 = wrap_into(nxt1 + 1, ncities);
+
+ picked = parent1->cycle[nxt1];
+ }
+ else if (tcused[parent1->cycle[nxt1]]) {
+ picked = parent2->cycle[nxt2];
+ }
+ else if (tcused[parent2->cycle[nxt2]]) {
+ picked = parent1->cycle[nxt1];
+ }
+ else {
+ if (cost_cache[last][parent1->cycle[nxt1]] <
+ cost_cache[last][parent2->cycle[nxt2]]) {
+
+ picked = parent1->cycle[nxt1];
+ }
+ else {
+ picked = parent2->cycle[nxt2];
+ }
+ }
+ child->cycle[cncyc] = picked;
+ tcused[picked] = 1;
+ last = picked;
+ cncyc++;
+ }
+
+ child->cycle[ncities] = child->cycle[0];
+ child->cost = 0;
+ for (i = 0; i < ncities; i++)
+ child->cost += cost_cache[child->cycle[i]][child->cycle[i + 1]];
+
+ if (child->cost < parent1->cost && child->cost < parent2->cost)
+ child->used = 1;
+ else
+ child->used = 1;
+
+ child->opt_done = 0;
+
+ return child->used;
+}
+
+
+int ga_recombine(struct tsp_tour *tour, int p1, int p2, int c)
+{
+ int i, start_i, best_i;
+ double best_cost;
+
+ G_debug(1, "ga_recombine()");
+
+ /* select subtour resulting in shortest total tour:
+ * genetic engineering */
+
+ best_cost = best_i = start_i = -1;
+
+ i = start_i + 1;
+ while (i < ncities) {
+ ga_subtour(tour, p1, p2, c, i, &start_i);
+ if (tour[c].used) {
+ if (best_cost < 0 || best_cost > tour[c].cost) {
+ best_cost = tour[c].cost;
+ best_i = start_i;
+ }
+ }
+ i = start_i + 1;
+ }
+
+ if (best_cost < 0 || (best_cost > tour[p1].cost - 0.0001 &&
+ best_cost > tour[p2].cost - 0.0001)) {
+ tour[c].used = 0;
+ }
+ else {
+ ga_subtour(tour, p1, p2, c, best_i, &start_i);
+ tour[c].used = 1;
+ }
+
+ return tour[c].used;
+}
+
+
+int cmp_tour_cost(const void *pa, const void *pb)
+{
+ struct tour_cost *p1 = (struct tour_cost *) pa;
+ struct tour_cost *p2 = (struct tour_cost *) pb;
+
+ if (p1->cost < p2->cost)
+ return -1;
+
+ if (p1->cost > p2->cost)
+ return 1;
+
+ return 0;
+}
Property changes on: grass-addons/grass7/vector/v.net.salesman.opt/evolution.c
___________________________________________________________________
Added: svn:mime-type
+ text/x-csrc
Added: svn:eol-style
+ native
Added: grass-addons/grass7/vector/v.net.salesman.opt/local_proto.h
===================================================================
--- grass-addons/grass7/vector/v.net.salesman.opt/local_proto.h (rev 0)
+++ grass-addons/grass7/vector/v.net.salesman.opt/local_proto.h 2012-05-14 14:37:56 UTC (rev 51622)
@@ -0,0 +1,34 @@
+
+#define MAX_CHAIN_LENGTH 11
+
+
+typedef struct
+{
+ int city;
+ double cost;
+} COST;
+
+
+extern int ncities; /* number of cities */
+extern int nnodes; /* number of nodes */
+extern int *cities; /* array of cities */
+extern COST **costs; /* pointer to array of pointers to arrays of sorted forward costs */
+extern COST **bcosts; /* pointer to array of pointers to arrays of sorted backward costs */
+extern double **cost_cache; /* pointer to array of pointers to arrays of cached costs */
+extern int debug_level;
+
+
+/* tour.c */
+void add_city(int, int, int *, int *, int *);
+int build_tour(int *cycle, int *cused, int *tncyc, int, int);
+
+/* optimize.c */
+int wrap_into(int i, int n);
+int optimize_nbrs(int, int, int *);
+int opt_2opt(int *tcycle, int *tcused, int tncyc);
+int optimize_tour(int *, int *, int, int, int, int, int);
+int optimize_tour_chains(int, int, int, int*, int *, int, int, int, int);
+
+/* ga.c */
+int ga_opt(int ntours, int nelim, int nopt, int ngen, int *);
+void ga_add_city(int city, int after, int *tnc, int *tcycle, int *tused);
Property changes on: grass-addons/grass7/vector/v.net.salesman.opt/local_proto.h
___________________________________________________________________
Added: svn:mime-type
+ text/x-chdr
Added: svn:eol-style
+ native
Modified: grass-addons/grass7/vector/v.net.salesman.opt/main.c
===================================================================
--- grass-addons/grass7/vector/v.net.salesman.opt/main.c 2012-05-14 14:31:54 UTC (rev 51621)
+++ grass-addons/grass7/vector/v.net.salesman.opt/main.c 2012-05-14 14:37:56 UTC (rev 51622)
@@ -5,7 +5,7 @@
*
* AUTHOR(S): Radim Blazek, Markus Metz
*
- * PURPOSE: Create a tour connecting given nodes.
+ * PURPOSE: Create a cycle connecting given nodes.
*
* COPYRIGHT: (C) 2001-2011 by the GRASS Development Team
*
@@ -24,6 +24,9 @@
#include <grass/glocale.h>
#include "local_proto.h"
+/* use EUC_2D distances for TSPLIB test data */
+#define TSP_TEST 0
+
/* TODO: Use some better algorithm */
int ncities; /* number of cities */
@@ -44,32 +47,7 @@
return (cities[city]);
}
-void add_city(int city, int after)
-{ /* index !!! to cycle, after which to put it */
- int i, j;
- if (after == -1) {
- cycle[0] = city;
- }
- else {
- /* for a large number of cities this will become slow */
- for (j = ncyc - 1; j > after; j--)
- cycle[j + 1] = cycle[j];
-
- cycle[after + 1] = city;
- }
- ncyc++;
- cused[city] = 1;
-
- if (debug_level >= 2) {
- G_debug(2, "Cycle:");
- for (i = 0; i < ncyc; i++) {
- G_debug(2, "%d: %d: %d", i, cycle[i], cities[cycle[i]]);
- }
- }
-
-}
-
/* like Vect_list_append(), but allows duplicates */
int tsp_list_append(struct ilist *list, int val)
{
@@ -97,14 +75,14 @@
int node, node1, node2, line;
int last_opt = 0, ostep, optimize;
struct Option *map, *output, *afield_opt, *tfield_opt, *afcol, *abcol,
- *seq, *type_opt, *term_opt;
- struct Flag *geo_f, *opt_f;
+ *seq, *type_opt, *term_opt, *opt_opt;
+ struct Flag *geo_f;
struct GModule *module;
struct Map_info Map, Out;
struct ilist *TList; /* list of terminal nodes */
struct ilist *List;
- struct ilist *StArcs; /* list of arcs on Steiner tree */
- struct ilist *StNodes; /* list of nodes on Steiner tree */
+ struct ilist *StArcs; /* list of arcs on tour */
+ struct ilist *StNodes; /* list of nodes on tour */
double cost, tmpcost, tcost;
struct cat_list *Clist;
struct line_cats *Cats;
@@ -112,6 +90,7 @@
const char *dstr;
const char *seqname;
int seq2stdout;
+ char *desc;
FILE *fp;
/* Initialize the GIS calls */
@@ -122,10 +101,10 @@
G_add_keyword(_("network"));
G_add_keyword(_("salesman"));
module->label =
- _("Creates a tour connecting given nodes (Traveling salesman problem).");
+ _("Creates a cycle connecting given nodes (Traveling salesman problem).");
module->description =
_("Note that TSP is NP-hard, heuristic algorithm is used by "
- "this module and created tour may be suboptimal");
+ "this module and created cycle may be suboptimal");
map = G_define_standard_option(G_OPT_V_INPUT);
output = G_define_standard_option(G_OPT_V_OUTPUT);
@@ -169,15 +148,25 @@
term_opt->description = _("Categories of points ('cities') on nodes "
"(layer is specified by nlayer)");
+ opt_opt = G_define_option();
+ opt_opt->type = TYPE_STRING;
+ opt_opt->key = "method";
+ opt_opt->required = NO;
+ opt_opt->options = "bs,ga";
+ opt_opt->description = _("Optimization method");
+ desc = NULL;
+ G_asprintf(&desc,
+ "bs;%s;"
+ "ga;%s;",
+ _("bootstrapping"),
+ _("genetic algorithm"));
+ opt_opt->descriptions = desc;
+
geo_f = G_define_flag();
geo_f->key = 'g';
geo_f->description =
_("Use geodesic calculation for longitude-latitude locations");
- opt_f = G_define_flag();
- opt_f->key = 'o';
- opt_f->description = _("Optimize the tour");
-
if (G_parser(argc, argv))
exit(EXIT_FAILURE);
@@ -211,7 +200,17 @@
}
geo = geo_f->answer;
- optimize = opt_f->answer;
+
+ if (opt_opt->answer) {
+ if (opt_opt->answer[0] == 'b')
+ optimize = 1;
+ else if (opt_opt->answer[0] == 'g')
+ optimize = 2;
+ else
+ G_fatal_error(_("Unknown method '%s'"), opt_opt->answer);
+ }
+ else
+ optimize = 0;
Vect_check_input_output_name(map->answer, output->answer, G_FATAL_EXIT);
@@ -290,13 +289,30 @@
if (i == j)
continue;
- ret =
- Vect_net_shortest_path(&Map, cities[i], cities[j], NULL,
- &cost);
+ if (!TSP_TEST) {
+ ret =
+ Vect_net_shortest_path(&Map, cities[i], cities[j], NULL,
+ &cost);
- if (ret == -1)
- G_fatal_error(_("Destination node [%d] is unreachable "
- "from node [%d]"), cities[i], cities[j]);
+ if (ret == -1)
+ G_fatal_error(_("Destination node [%d] is unreachable "
+ "from node [%d]"), cities[i], cities[j]);
+ }
+ else {
+ double x1, y1, z1, x2, y2, z2, dx, dy;
+
+ Vect_get_node_coor(&Map, cities[i], &x1, &y1, &z1);
+ Vect_get_node_coor(&Map, cities[j], &x2, &y2, &z2);
+
+ if (geo) {
+ cost = G_distance(x1, y1, x2, y2);
+ }
+ else {
+ dx = x1 - x2;
+ dy = y1 - y2;
+ cost = sqrt(dx * dx + dy * dy);
+ }
+ }
/* add to directional cost cache: from, to, cost */
costs[i][k].city = j;
@@ -337,188 +353,223 @@
}
}
- G_message(_("Searching for the shortest tour..."));
- /* find 2 cities with largest distance */
- cost = city = -1;
- for (i = 0; i < ncities; i++) {
- tmpcost = costs[i][ncities - 2].cost;
- if (tmpcost > cost) {
- cost = tmpcost;
- city = i;
- }
+ if (ncities < 5 && optimize) {
+ G_message(_("Optimization is not necessary for less than 5 cities"));
+ optimize = 0;
}
- G_debug(2, "biggest costs %d - %d", city,
- costs[city][ncities - 2].city);
- /* add these 2 cities to array */
- add_city(city, -1);
- add_city(costs[city][ncities - 2].city, 0);
+ if (optimize < 2) {
+ G_message(_("Searching for the shortest tour..."));
+ /* find 2 cities with largest distance */
+ cost = city = -1;
+ for (i = 0; i < ncities; i++) {
+ tmpcost = costs[i][ncities - 2].cost;
+ if (tmpcost > cost) {
+ cost = tmpcost;
+ city = i;
+ }
+ }
- /* In each step, find not used city, with biggest cost to any used city, and insert
- * into cycle between 2 nearest nodes */
- /* for a large number of cities this will become very slow, can be fixed */
+ G_debug(2, "biggest costs %d - %d", city,
+ costs[city][ncities - 2].city);
- /* for (i = 0; i < ncities - 2; i++) { */
- while (ncyc < ncities) {
- G_percent(ncyc, ncities, 1);
- cost = -1;
- G_debug(2, "---- city %d ----", i);
- for (j = 0; j < ncities; j++) {
- if (cused[j])
- continue;
- tmpcost = 0;
- for (k = 0; k < ncities - 1; k++) {
- G_debug(2, "forward? %d (%d) - %d (%d)", j, cnode(j),
- costs[j][k].city, cnode(costs[j][k].city));
- if (!cused[costs[j][k].city])
- continue; /* only used */
- /* directional costs j -> k */
- tmpcost += costs[j][k].cost;
- break; /* first nearest */
- }
- /* forward/backward: tmpcost = min(fcost) + min(bcost) */
- if (bcosts) {
+ /* add these 2 cities to array */
+ add_city(city, -1, &ncyc, cycle, cused);
+ add_city(costs[city][ncities - 2].city, 0, &ncyc, cycle, cused);
+
+ /* In each step, find not used city with biggest cost to any used city,
+ * and insert into cycle between 2 nearest nodes */
+ /* for a large number of cities this will become very slow, can be fixed */
+
+ /* for (i = 0; i < ncities - 2; i++) { */
+ while (ncyc < ncities) {
+ G_percent(ncyc, ncities, 1);
+ cost = -1;
+ G_debug(2, "---- city %d ----", i);
+ for (j = 0; j < ncities; j++) {
+ if (cused[j])
+ continue;
+ tmpcost = 0;
for (k = 0; k < ncities - 1; k++) {
- G_debug(2, "backward? %d (%d) - %d (%d)", j, cnode(j),
- bcosts[j][k].city, cnode(bcosts[j][k].city));
- if (!cused[bcosts[j][k].city])
+ G_debug(2, "forward? %d (%d) - %d (%d)", j, cnode(j),
+ costs[j][k].city, cnode(costs[j][k].city));
+ if (!cused[costs[j][k].city])
continue; /* only used */
- /* directional costs k -> j */
- tmpcost += bcosts[j][k].cost;
+ /* directional costs j -> k */
+ tmpcost += costs[j][k].cost;
break; /* first nearest */
}
- }
+ /* forward/backward: tmpcost = min(fcost) + min(bcost) */
+ if (bcosts) {
+ for (k = 0; k < ncities - 1; k++) {
+ G_debug(2, "backward? %d (%d) - %d (%d)", j, cnode(j),
+ bcosts[j][k].city, cnode(bcosts[j][k].city));
+ if (!cused[bcosts[j][k].city])
+ continue; /* only used */
+ /* directional costs k -> j */
+ tmpcost += bcosts[j][k].cost;
+ break; /* first nearest */
+ }
+ }
- G_debug(2, " cost = %f x %f", tmpcost, cost);
- if (tmpcost > cost) {
- cost = tmpcost;
- city = j;
+ G_debug(2, " cost = %f x %f", tmpcost, cost);
+ if (tmpcost > cost) {
+ cost = tmpcost;
+ city = j;
+ }
}
- }
- G_debug(2, "add city %d", city);
+ G_debug(2, "add city %d", city);
- /* add to cycle on lowest costs */
- cycle[ncyc] = cycle[0]; /* temporarily close the cycle */
- cost = PORT_DOUBLE_MAX;
- city1 = 0;
- for (j = 0; j < ncyc; j++) {
- /* cost from j to j + 1 (directional) */
- /* get cost from directional cost cache */
- tcost = cost_cache[cycle[j]][cycle[j + 1]];
- tmpcost = -tcost;
+ /* add to cycle on lowest costs */
+ cycle[ncyc] = cycle[0]; /* temporarily close the cycle */
+ cost = PORT_DOUBLE_MAX;
+ city1 = 0;
+ for (j = 0; j < ncyc; j++) {
+ /* cost from j to j + 1 (directional) */
+ /* get cost from directional cost cache */
+ tcost = cost_cache[cycle[j]][cycle[j + 1]];
+ tmpcost = -tcost;
- /* check insertion of city between j and j + 1 */
+ /* check insertion of city between j and j + 1 */
- /* cost from j to city (directional) */
- /* get cost from directional cost cache */
- tcost = cost_cache[cycle[j]][city];
- tmpcost += tcost;
- /* cost from city to j + 1 (directional) */
- /* get cost from directional cost cache */
- tcost = cost_cache[city][cycle[j + 1]];
- tmpcost += tcost;
-
- /* tmpcost must always be > 0 */
+ /* cost from j to city (directional) */
+ /* get cost from directional cost cache */
+ tcost = cost_cache[cycle[j]][city];
+ tmpcost += tcost;
+ /* cost from city to j + 1 (directional) */
+ /* get cost from directional cost cache */
+ tcost = cost_cache[city][cycle[j + 1]];
+ tmpcost += tcost;
+
+ /* tmpcost must always be > 0 */
- /* always true for j = 0 */
- if (tmpcost < cost) {
- city1 = j;
- cost = tmpcost;
+ /* always true for j = 0 */
+ if (tmpcost < cost) {
+ city1 = j;
+ cost = tmpcost;
+ }
}
- }
- add_city(city, city1);
+ add_city(city, city1, &ncyc, cycle, cused);
- /* keep ostep sufficiently large */
- ostep = sqrt(ncyc) * 2 / 3;
+ if (optimize && ncyc > 4) {
+ int chain_length;
+ int success;
+
+ success = optimize_nbrs(city1, ncyc, cycle);
+ if (success)
+ G_debug(3, "Neighbours optimized? %s",
+ success ? "Yes" : "No");
- if (ostep < 2)
- ostep = 2;
- if (ostep > MAX_CHAIN_LENGTH && ncities < 3000) /* max 3000 */
- ostep = MAX_CHAIN_LENGTH;
- if (optimize && ncyc > 5) {
- int chain_length;
- int success;
-
- success = optimize_nbrs(city1);
- if (success)
- G_debug(3, "Neighbours optimized? %s",
- success ? "Yes" : "No");
+ /* keep ostep sufficiently large */
+ ostep = sqrt(ncyc) * 2 / 3;
- if (last_opt < ncyc - ostep) {
-
- /* use brute force optimization
- * don't overdo it here, it's costly in terms of time
- * and the benefits are small */
+ if (ostep < 2)
+ ostep = 2;
+ if (ostep > MAX_CHAIN_LENGTH && ncities < 3000) /* max 3000 */
+ ostep = MAX_CHAIN_LENGTH;
- chain_length = ostep;
- if (chain_length > ncyc / 2)
- chain_length = ncyc / 2;
- if (chain_length > MAX_CHAIN_LENGTH)
- chain_length = MAX_CHAIN_LENGTH;
- if (chain_length == 0)
- chain_length = 1;
+ if (last_opt < ncyc - ostep) {
+
+ /* use brute force optimization
+ * don't overdo it here, it's costly in terms of time
+ * and the benefits are small */
- success = optimize_tour_chains(chain_length, chain_length,
- ostep, 0, 0);
- if (success)
- G_debug(3, "Preliminary tour optimized? %s",
- success ? "Yes" : "No");
+ chain_length = ostep;
+ if (chain_length > ncyc / 2)
+ chain_length = ncyc / 2;
+ if (chain_length > MAX_CHAIN_LENGTH)
+ chain_length = MAX_CHAIN_LENGTH;
+ if (chain_length == 0)
+ chain_length = 1;
- last_opt = ncyc;
- }
-
- } /* optimize done */
- }
-
- /* TODO: optimize tour using some refined Lin–Kernighan method */
- if (optimize && ncyc > 3) {
- double ocost, ocost1;
- int max_chain_length, chain_length;
- int success = 0;
- int optiter = 1;
+ success = optimize_tour_chains(chain_length, chain_length,
+ ostep, cycle, cused, ncyc, 0, 0, 1);
+ if (success)
+ G_debug(3, "Preliminary tour optimized? %s",
+ success ? "Yes" : "No");
- /* brute force bootstrapping */
- ocost = 0.;
- cycle[ncyc] = cycle[0];
- for (j = 0; j < ncyc; j++) {
- ocost += cost_cache[cycle[j]][cycle[j + 1]];
+ last_opt = ncyc;
+ }
+ } /* optimize done */
}
- G_verbose_message(_("Current total cost: %.3f"), ocost);
+
+ /* tour is complete */
- max_chain_length = MAX_CHAIN_LENGTH < ncyc / 2 ? MAX_CHAIN_LENGTH : ncyc / 2;
+ if (optimize && ncyc > 3) {
+ double ocost, ocost1;
+ int max_chain_length, chain_length;
+ int success = 0;
+ int optiter = 1;
- G_message(_("Optimizing..."));
+ /* brute force bootstrapping */
+ ocost = 0.;
+ cycle[ncyc] = cycle[0];
+ for (j = 0; j < ncyc; j++) {
+ ocost += cost_cache[cycle[j]][cycle[j + 1]];
+ }
+ G_verbose_message(_("Current total cost: %.3f"), ocost);
- for (chain_length = max_chain_length; chain_length > 0; chain_length--) {
- success = 0;
+ max_chain_length = MAX_CHAIN_LENGTH < ncyc / 2 ? MAX_CHAIN_LENGTH : ncyc / 2;
- G_message("%d. iteration...", optiter++);
+ G_message(_("Optimizing..."));
- success = optimize_tour_chains(max_chain_length, 1, 1, 1, 1);
-
- if (!success)
- break;
- }
- success = 1;
- while (success) {
- success = 0;
- G_message("%d. iteration...", optiter++);
- for (i = 0; i < ncyc; i++) {
- if (optimize_nbrs(i))
- success = 1;
+ for (chain_length = max_chain_length; chain_length > 0; chain_length--) {
+ success = 0;
+
+ G_message("%d. iteration...", optiter++);
+
+ if (max_chain_length >= 5)
+ optimize_tour_chains(max_chain_length, max_chain_length, 1,
+ cycle, cused, ncyc, 1, 1, 0);
+ else
+ optimize_tour_chains(max_chain_length, max_chain_length, 1,
+ cycle, cused, ncyc, 1, 1, 0);
+
+ success = optimize_tour_chains(max_chain_length, 1, 1,
+ cycle, cused, ncyc, 1, 1, 1);
+
+ if (!success)
+ break;
}
- }
+ success = 1;
+ while (success) {
+ success = 0;
+ G_message("%d. iteration...", optiter++);
+ for (i = 0; i < ncyc; i++) {
+ if (optimize_nbrs(i, ncyc, cycle))
+ success = 1;
+ }
+ }
- ocost1 = 0.;
- cycle[ncyc] = cycle[0];
- for (j = 0; j < ncyc; j++) {
- ocost1 += cost_cache[cycle[j]][cycle[j + 1]];
- }
- G_verbose_message(_("Optimized total cost: %.3f, gain %.3f%%"),
- ocost1, ocost / ocost1 * 100 - 100);
- } /* optimize done */
+ ocost1 = 0.;
+ cycle[ncyc] = cycle[0];
+ for (j = 0; j < ncyc; j++) {
+ ocost1 += cost_cache[cycle[j]][cycle[j + 1]];
+ }
+ G_verbose_message(_("Optimized total cost: %.3f, gain %.3f%%"),
+ ocost1, ocost / ocost1 * 100 - 100);
+ } /* optimize done */
+ }
+ else {
+ int ntours, nelim, nopt, ngen;
+
+ /* number of tours to work with */
+ if (ncities > 20)
+ ntours = 20;
+ else
+ ntours = ncities - 1;
+ /* number of tours to eliminate (< ntours / 2) */
+ nelim = ntours * 0.4;
+ /* number of tours to optimize (<= ntours) */
+ nopt = ntours * 1;
+ /* max number of generations */
+ ngen = 200;
+
+ ga_opt(ntours, nelim, nopt, ngen, cycle);
+ }
+
if (debug_level >= 2) {
/* debug print */
G_debug(2, "Tour:");
@@ -551,7 +602,6 @@
Vect_open_new(&Out, output->answer, Vect_is_3d(&Map));
Vect_hist_command(&Out);
- G_verbose_message(_("Tour with total cost %.3f"), cost);
G_debug(2, "Arcs' categories (layer %d, %d arcs):", afield,
StArcs->n_values);
@@ -620,6 +670,8 @@
Vect_close(&Map);
Vect_close(&Out);
+ G_message(_("Tour with total cost %.3f"), cost);
+
if (fp) {
fclose(fp);
if (seq2stdout) {
Added: grass-addons/grass7/vector/v.net.salesman.opt/optimize.c
===================================================================
--- grass-addons/grass7/vector/v.net.salesman.opt/optimize.c (rev 0)
+++ grass-addons/grass7/vector/v.net.salesman.opt/optimize.c 2012-05-14 14:37:56 UTC (rev 51622)
@@ -0,0 +1,411 @@
+#include <stdlib.h>
+#include <string.h>
+#include <grass/gis.h>
+#include <grass/vector.h>
+#include <grass/glocale.h>
+#include "local_proto.h"
+
+
+int wrap_into(int i, int n)
+{
+ /* wrap to range 0, n - 1 */
+ while (i >= n)
+ i -= n;
+ while (i < 0)
+ i += n;
+
+ return i;
+}
+
+int wrap_towards(int i, int n)
+{
+ /* wrap to range -n, 2n - 1 */
+ while (i - n >= n)
+ i -= n;
+ while (i + n < 0)
+ i += n;
+
+ return i;
+}
+
+/* swap nodes if result is better */
+int optimize_nbrs(int city1, int tncyc, int *tcycle)
+{
+ int j, k, city, city2;
+ int prev, nxt, prevk, nxtk;
+ double ocost1, ocost2, tmpcost, gain;
+ int success = 0;
+
+ /* check alternatives for city1 and city1 + 2 */
+ for (j = city1; j < city1 + 3; j += 2) {
+
+ /* city += j; */
+ city = wrap_into(j, tncyc);
+ prev = wrap_into(city - 1, tncyc);
+ nxt = wrap_into(city + 1, tncyc);
+ ocost1 = cost_cache[tcycle[prev]][tcycle[city]] +
+ cost_cache[tcycle[city]][tcycle[nxt]];
+ gain = 0;
+ tmpcost = 0;
+ city2 = tcycle[city];
+ for (k = 0; k < tncyc; k++) {
+
+ prevk = wrap_into(k - 1, tncyc);
+ nxtk = wrap_into(k + 1, tncyc);
+
+ /* original costs */
+ ocost2 = cost_cache[tcycle[prevk]][tcycle[k]] +
+ cost_cache[tcycle[k]][tcycle[nxtk]];
+
+ /* new costs:
+ * (city1 - 1) -> k -> (city1 + 1)
+ * (k - 1) -> city1 -> (k + 1) */
+ tmpcost = cost_cache[tcycle[prev]][tcycle[k]] +
+ cost_cache[tcycle[k]][tcycle[nxt]] +
+ cost_cache[tcycle[prevk]][tcycle[city]] +
+ cost_cache[tcycle[city]][tcycle[nxtk]];
+
+ if (ocost1 + ocost2 - tmpcost > gain) {
+ gain = ocost1 + ocost2 - tmpcost;
+ city2 = k;
+ }
+ }
+ if (gain > 0 && tcycle[city2] != tcycle[city] &&
+ tcycle[city2] != tcycle[prev] && tcycle[city2] != tcycle[nxt]) {
+ /* swap cities */
+ G_debug(3, "swap cities");
+ tcycle[tncyc] = tcycle[city2];
+ tcycle[city2] = tcycle[city];
+ tcycle[city] = tcycle[tncyc];
+ tcycle[tncyc] = tcycle[0];
+ success = 1;
+ }
+ } /* neighbor check done */
+
+ return success;
+}
+
+/* standard 2opt method
+ * disadvantage: requires reversing a subtour, not good for
+ * asymmetric problems, therefore not used */
+int opt_2opt(int *tcycle, int *tcused, int tncyc)
+{
+ int i, j, k, nxt1, nxt2;
+ int pivot1;
+ double cost, new_cost;
+ static int *tmpcycle = NULL;
+
+ if (tmpcycle == NULL)
+ tmpcycle = G_malloc((ncities + 1) * sizeof(int));
+
+ tcycle[tncyc] = tcycle[0];
+
+ for (i = 0; i < tncyc; i++) {
+ for (j = 0; j < tncyc; j++) {
+ if (j == i || j == i + 1)
+ continue;
+
+ nxt2 = wrap_into(j + 1, tncyc);
+ if (nxt2 == i)
+ continue;
+
+ cost = cost_cache[tcycle[i]][tcycle[i + 1]] +
+ cost_cache[tcycle[j]][tcycle[j + 1]];
+
+ new_cost = cost_cache[tcycle[i]][tcycle[j + 1]] +
+ cost_cache[tcycle[j]][tcycle[i + 1]];
+
+ if (new_cost > cost)
+ continue;
+
+ /* go from i to j,
+ * backwards until we meet i + 1, go from i + 1 to j + 1,
+ * forwards until we meet i again */
+ nxt1 = i;
+ tmpcycle[nxt1] = tcycle[i];
+ nxt1 = wrap_into(nxt1 + 1, tncyc);
+ tmpcycle[nxt1] = tcycle[j];
+
+ nxt2 = j;
+ nxt1 = wrap_into(nxt1 + 1, tncyc);
+ nxt2 = wrap_into(nxt2 - 1, tncyc);
+
+ pivot1 = wrap_into(i + 1, tncyc);
+ while (nxt2 != pivot1) {
+ tmpcycle[nxt1] = tcycle[nxt2];
+ nxt1 = wrap_into(nxt1 + 1, tncyc);
+ nxt2 = wrap_into(nxt2 - 1, tncyc);
+ }
+ /* nxt2 = i + 1 */
+ tmpcycle[nxt1] = tcycle[nxt2];
+ nxt1 = wrap_into(nxt1 + 1, tncyc);
+ nxt2 = wrap_into(j + 1, tncyc);
+ tmpcycle[nxt1] = tcycle[nxt2];
+
+ nxt1 = wrap_into(nxt1 + 1, tncyc);
+ nxt2 = wrap_into(nxt2 + 1, tncyc);
+ while (nxt2 != i) {
+ tmpcycle[nxt1] = tcycle[nxt2];
+ nxt1 = wrap_into(nxt1 + 1, tncyc);
+ nxt2 = wrap_into(nxt2 + 1, tncyc);
+ }
+ for (k = 0; k < tncyc; k++)
+ tcycle[k] = tmpcycle[k];
+
+ tcycle[tncyc] = tcycle[0];
+ }
+ }
+
+ return 1;
+}
+
+/* bootstrapping:
+ * cut out subtour of length chain_length
+ * reinsert the removed nodes */
+
+int optimize_tour(int *tcycle, int *tcused, int tncyc, int chain_length,
+ int show_perc, int opt_high, int subtour)
+{
+ int i, j, k;
+ int city, city1;
+ int cl1, cl2, cl3, cstep;
+ double cost, tcost, tmpcost, ocost1, ocost2;
+ int nl[MAX_CHAIN_LENGTH];
+ static int *ccycle = NULL;
+ int success = 0, found;
+ int ncyc_noopt, n_removed;
+
+ if (!ccycle)
+ ccycle = (int *)G_malloc((ncities + 1) * sizeof(int)); /* + 1 is for output cycle */
+
+ cl1 = chain_length;
+ cstep = cl1 * 2 / 3; /* adjust */
+
+ if (cstep == 0)
+ cstep = 1;
+
+ ncyc_noopt = tncyc;
+
+ if (show_perc)
+ G_message("Chain length %d", cl1);
+
+ for (i = 0; i < ncyc_noopt; i += cstep) {
+ if (show_perc)
+ G_percent(i, ncyc_noopt, 2);
+
+ /* copy of current cycle */
+ for (j = 0; j < tncyc; j++)
+ ccycle[j] = tcycle[j];
+
+ /* remove <cl1> nodes */
+ n_removed = 0;
+
+ if (subtour) {
+ /* cut out subtour */
+ for (cl2 = 0; cl2 < tncyc; cl2++) {
+ j = wrap_into(i, tncyc);
+ nl[n_removed] = tcycle[j];
+
+ if (!tcused[tcycle[j]])
+ continue;
+
+ G_debug(3, "removing city %d, cycle %d", nl[n_removed], j);
+ tcused[tcycle[j]] = 0;
+ for (k = j; k < tncyc; k++)
+ tcycle[k] = tcycle[k + 1];
+ tncyc--;
+ n_removed++;
+ if (n_removed == cl1)
+ break;
+ }
+ }
+ else {
+ /* cut our nearest */
+ j = wrap_into(i, tncyc);
+ city = tcycle[j];
+ nl[n_removed] = city;
+ tcused[city] = 0;
+ found = 0;
+ for (k = 0; k < tncyc; k++) {
+ if (tcycle[k] == city)
+ found = 1;
+ if (found)
+ tcycle[k] = tcycle[k + 1];
+ }
+ tncyc--;
+
+ n_removed++;
+
+ for (cl2 = 0; cl2 < tncyc - 1; cl2++) {
+
+ if (n_removed == cl1)
+ break;
+
+ city = costs[tcycle[j]][cl2].city;
+
+ if (!tcused[city])
+ continue;
+
+ nl[n_removed] = city;
+ G_debug(3, "removing city %d", nl[n_removed]);
+ tcused[city] = 0;
+ found = 0;
+ for (k = 0; k < tncyc; k++) {
+ if (tcycle[k] == city)
+ found = 1;
+ if (found)
+ tcycle[k] = tcycle[k + 1];
+ }
+ tncyc--;
+
+ n_removed++;
+ if (n_removed == cl1)
+ break;
+ }
+ }
+
+ if (!n_removed)
+ return 1;
+
+ /* reinsert <cl1> nodes */
+ for (cl2 = 0; cl2 < n_removed; cl2++) {
+ if (opt_high)
+ cost = -1;
+ else
+ cost = PORT_DOUBLE_MAX;
+ city = -1;
+ G_debug(3, "---- city %d ----", i);
+ for (cl3 = 0; cl3 < n_removed; cl3++) {
+ j = nl[cl3];
+ if (tcused[j])
+ continue;
+ tmpcost = 0;
+ for (k = 0; k < ncities - 1; k++) {
+ G_debug(2, "forward? %d (%d) - %d (%d)", j, cities[j],
+ costs[j][k].city, cities[costs[j][k].city]);
+ if (!tcused[costs[j][k].city])
+ continue; /* only used */
+ /* directional costs j -> k */
+ tmpcost += costs[j][k].cost;
+ break; /* first nearest */
+ }
+ /* forward/backward: tmpcost = min(fcost) + min(bcost) */
+ if (bcosts) {
+ for (k = 0; k < ncities - 1; k++) {
+ G_debug(2, "backward? %d (%d) - %d (%d)", j, cities[j],
+ bcosts[j][k].city, cities[bcosts[j][k].city]);
+ if (!tcused[bcosts[j][k].city])
+ continue; /* only used */
+ /* directional costs k -> j */
+ tmpcost += bcosts[j][k].cost;
+ break; /* first nearest */
+ }
+ }
+
+ G_debug(2, " cost = %f x %f", tmpcost, cost);
+ if (opt_high) {
+ /* pick farthest */
+ if (tmpcost > cost) {
+ cost = tmpcost;
+ city = j;
+ }
+ }
+ else {
+ /* pick nearest */
+ if (tmpcost < cost) {
+ cost = tmpcost;
+ city = j;
+ }
+ }
+ }
+ G_debug(3, "add city %d", city);
+
+ /* add to cycle on lowest costs */
+ tcycle[tncyc] = tcycle[0]; /* temporarily close the cycle */
+ cost = PORT_DOUBLE_MAX;
+ city1 = 0;
+ for (j = 0; j < tncyc; j++) {
+ /* cost from j to j + 1 (directional) */
+ /* get cost from directional cost cache */
+ tcost = cost_cache[tcycle[j]][tcycle[j + 1]];
+ tmpcost = -tcost;
+
+ /* check insertion of city between j and j + 1 */
+
+ /* cost from j to city (directional) */
+ /* get cost from directional cost cache */
+ tcost = cost_cache[tcycle[j]][city];
+ tmpcost += tcost;
+ /* cost from city to j + 1 (directional) */
+ /* get cost from directional cost cache */
+ tcost = cost_cache[city][tcycle[j + 1]];
+ tmpcost += tcost;
+
+ /* tmpcost must always be > 0 */
+
+ G_debug(2, "? %d - %d cost = %f x %f", j, city1, tmpcost,
+ cost);
+ /* always true for j = 0 */
+ if (tmpcost < cost) {
+ city1 = j;
+ cost = tmpcost;
+ }
+ }
+ add_city(city, city1, &tncyc, tcycle, tcused);
+
+ tcycle[tncyc] = tcycle[0];
+
+ if (tncyc > 3)
+ optimize_nbrs(city1, tncyc, tcycle);
+
+ } /* all nodes reinserted */
+
+ /* compare cycles */
+ ocost1 = ocost2 = 0.;
+ tcycle[tncyc] = tcycle[0];
+ ccycle[tncyc] = ccycle[0];
+ for (j = 0; j < tncyc; j++) {
+ ocost1 += cost_cache[ccycle[j]][ccycle[j + 1]];
+ ocost2 += cost_cache[tcycle[j]][tcycle[j + 1]];
+ }
+ if (ocost1 - ocost2 > 0.001) {
+ success = 1;
+ G_debug(3, "success for chain length %d, gain %.4f", cl1, ocost1 - ocost2);
+ }
+ else {
+ for (j = 0; j < tncyc; j++)
+ tcycle[j] = ccycle[j];
+ }
+ } /* chain length for all cities removed */
+
+ return success;
+}
+
+int optimize_tour_chains(int hi, int lo, int cstep, int *tcycle,
+ int *tused, int tncyc, int show_perc,
+ int opt_high, int subtour)
+{
+ int cl1, success = 0;
+
+ if (lo <= 0)
+ lo = 1;
+
+ if (hi > tncyc / 2)
+ hi = tncyc / 2;
+ if (hi > MAX_CHAIN_LENGTH)
+ hi = MAX_CHAIN_LENGTH;
+ if (hi <= 0)
+ hi = 1;
+ if (hi < lo)
+ lo = hi;
+
+ for (cl1 = hi; cl1 >= lo; cl1 -= cstep) {
+
+ if (!optimize_tour(tcycle, tused, tncyc, cl1, show_perc, opt_high, subtour))
+ break;
+ success = 1;
+ }
+
+ return success;
+}
+
Property changes on: grass-addons/grass7/vector/v.net.salesman.opt/optimize.c
___________________________________________________________________
Added: svn:mime-type
+ text/x-csrc
Added: svn:eol-style
+ native
Added: grass-addons/grass7/vector/v.net.salesman.opt/tour.c
===================================================================
--- grass-addons/grass7/vector/v.net.salesman.opt/tour.c (rev 0)
+++ grass-addons/grass7/vector/v.net.salesman.opt/tour.c 2012-05-14 14:37:56 UTC (rev 51622)
@@ -0,0 +1,136 @@
+#include <stdlib.h>
+#include <string.h>
+#include <grass/gis.h>
+#include <grass/vector.h>
+#include <grass/glocale.h>
+#include "local_proto.h"
+
+void add_city(int city, int after, int *tncyc, int *tcycle, int *tcused)
+{
+ int i, j;
+
+ /* after: index !!! to cycle, after which to put it */
+ if (after == -1) {
+ tcycle[0] = city;
+ }
+ else {
+ /* for a large number of cities this will become slow */
+ for (j = *tncyc - 1; j > after; j--)
+ tcycle[j + 1] = tcycle[j];
+
+ tcycle[after + 1] = city;
+ }
+ (*tncyc)++;
+ tcused[city] = 1;
+
+ if (debug_level >= 2) {
+ G_debug(2, "Cycle:");
+ for (i = 0; i < *tncyc; i++) {
+ G_debug(2, "%d: %d: %d", i, tcycle[i], cities[tcycle[i]]);
+ }
+ }
+
+}
+
+
+int build_tour(int *cycle, int *cused, int *tncyc, int optimize, int opt_high)
+{
+ int j, k;
+ int city, city1;
+ int ncyc = *tncyc;
+ double cost, tmpcost, tcost;
+
+ /* complete the tour */
+
+ /* In each step, find not used city with biggest cost to any used city,
+ * and insert into cycle between 2 nearest nodes */
+ /* for a large number of cities this will become very slow, can be fixed */
+
+ while (ncyc < ncities) {
+ if (opt_high)
+ cost = -1;
+ else
+ cost = PORT_DOUBLE_MAX;
+ city = -1;
+ for (j = 0; j < ncities; j++) {
+ if (cused[j])
+ continue;
+ tmpcost = 0;
+ for (k = 0; k < ncities - 1; k++) {
+ G_debug(2, "forward? %d (%d) - %d (%d)", j, cities[j],
+ costs[j][k].city, cities[costs[j][k].city]);
+ if (!cused[costs[j][k].city])
+ continue; /* only used */
+ /* directional costs j -> k */
+ tmpcost += costs[j][k].cost;
+ break; /* first nearest */
+ }
+ /* forward/backward: tmpcost = min(fcost) + min(bcost) */
+ if (bcosts) {
+ for (k = 0; k < ncities - 1; k++) {
+ G_debug(2, "backward? %d (%d) - %d (%d)", j, cities[j],
+ bcosts[j][k].city, cities[bcosts[j][k].city]);
+ if (!cused[bcosts[j][k].city])
+ continue; /* only used */
+ /* directional costs k -> j */
+ tmpcost += bcosts[j][k].cost;
+ break; /* first nearest */
+ }
+ }
+
+ G_debug(2, " cost = %f x %f", tmpcost, cost);
+ if (opt_high) {
+ if (tmpcost > cost) {
+ cost = tmpcost;
+ city = j;
+ }
+ }
+ else {
+ if (tmpcost < cost) {
+ cost = tmpcost;
+ city = j;
+ }
+ }
+ }
+ G_debug(2, "add city %d", city);
+
+ /* add to cycle on lowest costs */
+ cycle[ncyc] = cycle[0]; /* temporarily close the cycle */
+ cost = PORT_DOUBLE_MAX;
+ city1 = 0;
+ for (j = 0; j < ncyc; j++) {
+ /* cost from j to j + 1 (directional) */
+ /* get cost from directional cost cache */
+ tcost = cost_cache[cycle[j]][cycle[j + 1]];
+ tmpcost = -tcost;
+
+ /* check insertion of city between j and j + 1 */
+
+ /* cost from j to city (directional) */
+ /* get cost from directional cost cache */
+ tcost = cost_cache[cycle[j]][city];
+ tmpcost += tcost;
+ /* cost from city to j + 1 (directional) */
+ /* get cost from directional cost cache */
+ tcost = cost_cache[city][cycle[j + 1]];
+ tmpcost += tcost;
+
+ /* tmpcost must always be > 0 */
+
+ /* always true for j = 0 */
+ if (tmpcost < cost) {
+ city1 = j;
+ cost = tmpcost;
+ }
+ }
+ add_city(city, city1, &ncyc, cycle, cused);
+
+ if (optimize && ncyc > 3)
+ optimize_nbrs(city1, ncyc, cycle);
+ }
+
+ *tncyc = ncyc;
+
+ return 1;
+}
+
Property changes on: grass-addons/grass7/vector/v.net.salesman.opt/tour.c
___________________________________________________________________
Added: svn:mime-type
+ text/x-csrc
Added: svn:eol-style
+ native
Deleted: grass-addons/grass7/vector/v.net.salesman.opt/v.net.salesman.html
===================================================================
--- grass-addons/grass7/vector/v.net.salesman.opt/v.net.salesman.html 2012-05-14 14:31:54 UTC (rev 51621)
+++ grass-addons/grass7/vector/v.net.salesman.opt/v.net.salesman.html 2012-05-14 14:37:56 UTC (rev 51622)
@@ -1,145 +0,0 @@
-<h2>DESCRIPTION</h2>
-
-<em>v.net.salesman</em> calculates the optimal route to visit nodes on a
-vector network.
-
-<p>Costs may be either line lengths, or attributes saved in a database
-table. These attribute values are taken as costs of whole segments, not
-as costs to traverse a length unit (e.g. meter) of the segment.
-For example, if the speed limit is 100 km / h, the cost to traverse a
-10 km long road segment must be calculated as
-length / speed = 10 km / (100 km/h) = 0.1 h.
-Supported are cost assignments for arcs,
-and also different costs for both directions of a vector line.
-For areas, costs will be calculated along boundary lines.
-<p>The input vector needs to be prepared with <em>v.net operation=connect</em>
-in order to connect points representing center nodes to the network.
-
-<p>Points specified by category must be exactly on network nodes, and the
-input vector map needs to be prepared with <em>v.net operation=connect</em>.
-
-<h2>NOTES</h2>
-Arcs can be closed using cost = -1.
-
-<h2>EXAMPLE</h2>
-
-Traveling salesman for 6 digitized nodes (Spearfish):
-
-<p>Shortest path, along unimproved roads:
-<p><img src="vnetsalesman.png" alt="v.net.salesman example using distance" border="1">
-
-<p>Fastest path, along highways:
-<p><img src="vnetsalesmantime.png" alt="v.net.salesman example using time" border="1">
-
-<p>Searching for the shortest path using distance and the fastest path using
-traveling time according to the speed limits of different road types:
-
-<div class="code"><pre>
-# Spearfish
-
-g.copy vect=roads,myroads
-
-# we have 6 locations to visit on our trip
-echo "1|601653.5|4922869.2|a
-2|608284|4923776.6|b
-3|601845|4914981.9|c
-4|596270|4917456.3|d
-5|593330.8|4924096.6|e
-6|598005.5|4921439.2|f" | v.in.ascii in=- cat=1 x=2 y=3 out=centers col="cat integer, \
- east double precision, north double precision, label varchar(43)"
-
-# verify data preparation
-v.db.select centers
-v.category centers op=report
-# type count min max
-# point 6 1 6
-
-
-# create lines map connecting points to network (on layer 2)
-v.net myroads points=centers out=myroads_net op=connect thresh=500
-v.category myroads_net op=report
-# Layer / table: 1 / myroads_net
-# type count min max
-# line 837 1 5
-#
-# Layer: 2
-# type count min max
-# point 6 1 5
-
-# find the shortest path
-v.net.salesman myroads_net ccats=1-6 out=mysalesman_distance
-
-# set up costs as traveling time
-
-# create unique categories for each road in layer 3
-v.category in=myroads_net out=myroads_net_time opt=add cat=1 layer=3 type=line
-
-# add new table for layer 3
-v.db.addtable myroads_net_time layer=3 col="cat integer,label varchar(43),length double precision,speed double precision,cost double precision,bcost double precision"
-
-# copy road type to layer 3
-v.to.db myroads_net_time layer=3 qlayer=1 opt=query qcolumn=label columns=label
-
-# upload road length in miles
-v.to.db myroads_net_time layer=3 type=line option=length col=length unit=miles
-
-# set speed limits in miles / hour
-v.db.update myroads_net_time layer=3 col=speed val="5.0"
-v.db.update myroads_net_time layer=3 col=speed val="75.0" where="label='interstate'"
-v.db.update myroads_net_time layer=3 col=speed val="75.0" where="label='primary highway, hard surface'"
-v.db.update myroads_net_time layer=3 col=speed val="50.0" where="label='secondary highway, hard surface'"
-v.db.update myroads_net_time layer=3 col=speed val="25.0" where="label='light-duty road, improved surface'"
-v.db.update myroads_net_time layer=3 col=speed val="5.0" where="label='unimproved road'"
-
-# define traveling costs as traveling time in minutes:
-
-# set forward costs
-v.db.update myroads_net_time layer=3 col=cost val="length / speed * 60"
-# set backward costs
-v.db.update myroads_net_time layer=3 col=bcost val="length / speed * 60"
-
-# find the fastest path
-v.net.salesman myroads_net_time alayer=3 nlayer=2 afcol=cost abcol=bcost ccats=1-6 out=mysalesman_time
-</pre></div>
-
-To display the result, run for example:
-
-<div class="code"><pre>
-# Display the results
-g.region vect=myroads_net
-
-# shortest path
-d.mon x0
-d.vect myroads_net
-d.vect centers -c icon=basic/triangle
-d.vect mysalesman_distance col=green width=2
-d.font Vera
-d.vect centers col=red disp=attr attrcol=label lsize=12
-
-# fastest path
-d.mon x1
-d.vect myroads_net
-d.vect centers -c icon=basic/triangle
-d.vect mysalesman_time col=green width=2
-d.font Vera
-d.vect centers col=red disp=attr attrcol=label lsize=12
-</pre></div>
-
-
-<h2>SEE ALSO</h2>
-
-<em><a href="d.path.html">d.path</a></em>,
-<em><a href="v.net.html">v.net</a></em>,
-<em><a href="v.net.alloc.html">v.net.alloc</a></em>,
-<em><a href="v.net.iso.html">v.net.iso</a></em>,
-<em><a href="v.net.path.html">v.net.path</a></em>,
-<em><a href="v.net.steiner.html">v.net.steiner</a></em>
-
-<h2>AUTHOR</h2>
-
-Radim Blazek, ITC-Irst, Trento, Italy<br>
-Markus Metz<br>
-Documentation: Markus Neteler, Markus Metz
-
-
-<p><i>Last changed: $Date$</i>
Copied: grass-addons/grass7/vector/v.net.salesman.opt/v.net.salesman.opt.html (from rev 51621, grass-addons/grass7/vector/v.net.salesman.opt/v.net.salesman.html)
===================================================================
--- grass-addons/grass7/vector/v.net.salesman.opt/v.net.salesman.opt.html (rev 0)
+++ grass-addons/grass7/vector/v.net.salesman.opt/v.net.salesman.opt.html 2012-05-14 14:37:56 UTC (rev 51622)
@@ -0,0 +1,94 @@
+<h2>DESCRIPTION</h2>
+
+<em>v.net.salesman.opt</em> estimates the optimal route to visit nodes
+on a vector network and optionally tries to improve the result.
+
+<p>Costs may be either line lengths, or attributes saved in a database
+table. These attribute values are taken as costs of whole segments, not
+as costs to traverse a length unit (e.g. meter) of the segment.
+For example, if the speed limit is 100 km / h, the cost to traverse a
+10 km long road segment must be calculated as
+length / speed = 10 km / (100 km/h) = 0.1 h.
+Supported are cost assignments for arcs,
+and also different costs for both directions of a vector line.
+For areas, costs will be calculated along boundary lines.
+<p>The input vector needs to be prepared with <em>v.net operation=connect</em>
+in order to connect points representing center nodes to the network.
+
+<p>Points specified by category must be exactly on network nodes, and the
+input vector map needs to be prepared with <em>v.net operation=connect</em>.
+
+<h3>Optimization</h3>
+For less than 10 nodes, the initial result is usually very close to the
+shortest possible tour and further optimization will have no effect.
+<p>
+<em>v.net.salesman.opt</em> uses the same heuristics like
+<em>v.net.salesman.opt</em> to find a good tour which is not necessarily
+the optimal tour. The tour can be optimized with two methods:
+bootstrapping and genetic algorithm (GA). The bootstrapping method
+removes a subtour from the tour and reinserts the nodes resulting in a
+shorter tour. This is applied on intermediate tours and the final tour.
+<p>
+The genetic algorithm first creates several initial tours. From these
+tours nearly identical toursare eliminated. The remaining tours are
+recombined to create new, better tours. All tours are now mutated. The
+sequence of Selection, Recombination, Mutation is repeated until
+there is no better solution. Finally, the best tour is optimized with
+bootstrapping.
+
+<h2>NOTES</h2>
+Arcs can be closed using cost = -1.
+
+<h2>EXAMPLE</h2>
+
+Visiting all 167 schools (North Carolina):
+
+<div class="code"><pre>
+# North Carolina
+
+# prepare network by connecting schools to streets
+v.net -c input=streets_wake points=schools_wake output=streets_schools \
+ operation=connect alayer=1 nlayer=2 thresh=1000
+
+# verify data preparation
+v.category in=streets_schools op=report
+# type count min max
+# point 6 1 6
+
+# find the shortest path
+v.net.salesman.opt in=streets_schools ccats=1-167 out=schools_tour
+
+# Resulting tour length: 551551.662
+
+# find the shortest path, optimize with genetic algorithm
+v.net.salesman.opt in=streets_schools ccats=1-167 out=schools_tour_ga method=ga
+
+# Resulting tour length: 521191.083
+# The original tour was 5.8% longer
+
+# find the shortest path, optimize with bootstrapping
+v.net.salesman.opt in=streets_schools ccats=1-167 out=schools_tour_bs method=bs
+
+# Resulting tour length: 510693.965
+# The original tour was 8.0% longer
+
+</pre></div>
+
+
+<h2>SEE ALSO</h2>
+
+<em><a href="d.path.html">d.path</a></em>,
+<em><a href="v.net.html">v.net</a></em>,
+<em><a href="v.net.alloc.html">v.net.alloc</a></em>,
+<em><a href="v.net.iso.html">v.net.iso</a></em>,
+<em><a href="v.net.path.html">v.net.path</a></em>,
+<em><a href="v.net.steiner.html">v.net.steiner</a></em>
+
+<h2>AUTHOR</h2>
+
+Radim Blazek, ITC-Irst, Trento, Italy<br>
+Markus Metz<br>
+Documentation: Markus Neteler, Markus Metz
+
+
+<p><i>Last changed: $Date$</i>
More information about the grass-commit
mailing list