﻿

# tsp

a guest
Dec 11th, 2019
96
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
1. NS = 80;
2. xmax = 20; ymax = 20; delta = 1.0;
3. xs = zeros(NS, 1);
4. ys = zeros(NS, 1);
5. xs(1) = xmax*rand;
6. ys(1) = ymax*rand;
7. n = 1;
8. while n < NS
9.     x = xmax * rand;
10.     y = ymax * rand;
11.     if all(hypot(x-xs,y-ys)> delta)
12.         n = n+1;
13.         xs(n) =x; ys(n) = y;
14.     end
15. end
16.
17.
18.
19. dist = zeros(NS, NS);
20. for i= 1:NS
21.     for j = 1:NS
22.         dist(i, j) = hypot(xs(i)-xs(j), ys(i)-ys(j));
23.     end
24. end
25.
26. f = dist(:);
27. intcon = 1:length(f);
28. lb = zeros(length(f), 1);
29. ub = ones(length(f), 1);
30. i = 1:NS;
31. ub(i + (i-1)*NS) = 0.
32.
33. A = zeros(NS*NS, length(f));
34. b = zeros(size(A, 1), 1);
35. n = 1;
36. for i = 1:NS
37.     for j = 1:NS
38.         A(n, i+(j-1)*NS) = 1;
39.         A(n, j+(i-1)*NS) = 1;
40.         b(n) = 1;
41.         n = n+1;
42.     end
43. end
44. Aeq = zeros(2*NS, length(f));
45. beq = zeros(2*NS, 1);
46. n = 1;
47. for i = 1:NS
48.     for j = 1:NS
49.         Aeq(n, i+(j-1)*NS) = 1;
50.         Aeq(n+1, j+(i-1)*NS) = 1;
51.     end
52.     beq(n) = 1;
53.     beq(n+1) = 1;
54.     n = n + 2;
55. end
56.
57. while 1
58.     [xy, fval, exitflag] = intlinprog(f, intcon, A, b, Aeq, beq, lb, ub);
59.     xy = round(xy);
60.
61.     % draw
62.     plot(xs, ys, 'r*')
63.     hold on
64.     for j = 1:NS
65.         for i = 1:NS
66.             if xy(i+(j-1)*NS) == 1
67.                 plot([xs(i) xs(j)], [ys(i), ys(j)], 'b-')
68.             end
69.         end
70.     end
71.     hold off
72.     pause
73.     % draw end
74.
75.     xy = reshape(xy, NS, NS);
76.
77.     stops = zeros(NS, 1);
78.     for n = 1:NS
79.         subtour = find(stops==0, 1);
80.         while 1
81.            i = subtour(end);
82.            stops(i) = 1;
83.            j = find(xy(i, :), 1);
84.            if j == subtour(1), break; end
85.            subtour = [subtour j];
86.         end
87.         if n == 1 && all(stops==1)
88.             break;
89.         end
90.
91.         A = [A; zeros(1, length(f))];
92.         for i = subtour
93.             for j = subtour
94.                 A(end, i+(j-1)*NS) = 1;
95.             end
96.         end
97.         b = [b; length(subtour)-1];
98.         if all(stops==1) break; end
99.     end
100.
101.     if n == 1 break; end;
102. end
103.
104. plot(xs, ys, 'r*')
105. hold on
106. for j = 1:NS
107.     for i = 1:NS
108.         if xy(i+(j-1)*NS) == 1
109.             plot([xs(i) xs(j)], [ys(i), ys(j)], 'b-')
110.         end
111.     end
112. end
113. hold off
114. pause
RAW Paste Data