Coverage for moptipyapps / tsp / instance.py: 84%

1""" 

2An instance of the Traveling Salesperson Problem (TSP) as distance matrix. 

3 

4An instance of the Traveling Salesperson Problem (TSP) is defined as a 

5fully-connected graph with :attr:`Instance.n_cities` nodes. Each edge in the 

6graph has a weight, which identifies the distance between the nodes. The goal 

7is to find the *shortest* tour that visits every single node in the graph 

8exactly once and then returns back to its starting node. Then nodes are 

9usually called cities. In this file, we present methods for loading instances 

10of the TSP as distance matrices `A`. In other words, the value at `A[i, j]` 

11identifies the travel distance from `i` to `j`. 

12 

13We can load files in a subset of the TSPLib format [1, 2] and also include the 

14instances of TSPLib with no more than 2000 cities. We load instances as full 

15distance matrices, which makes writing algorithms to solve them easier but 

16handling instances with more than 2000 cities problematic due to the high 

17memory consumption. Of course, you could still load them, but it would be 

18ill-advised to do so on a normal PC (at the time of this writing). 

19 

20The TSPLib contains many instances of the symmetric TSP, where the distance 

21`A[i, j]` from city `i` to city `j` is the same as the distance `A[j, i]` from 

22`j` to `i`. Here, we provide 111 of them as resource. The cities of some of 

23these instances are points in the Euclidean plain and the distances are the 

24(approximate) Euclidean distances. Others use geographical coordinates 

25(longitude/latitude), and yet others are provides as distances matrices 

26directly. We also provide 19 of the asymmetric TSPLib instances, where the 

27distance `A[i, j]` from `i` to `j` may be different from the distance 

28`A[j, i]` from `j` to `i`. 

29 

30You can obtain lists of either all, only the symmetric, or only the asymmetric 

31instance resources via 

32 

33>>> print(Instance.list_resources()[0:10]) # get all instances (1..10) 

34('a280', 'ali535', 'att48', 'att532', 'bayg29', 'bays29', 'berlin52', \ 

35'bier127', 'br17', 'brazil58') 

36 

37>>> print(Instance.list_resources(asymmetric=False)[0:10]) # only symmetric 

38('a280', 'ali535', 'att48', 'att532', 'bayg29', 'bays29', 'berlin52', \ 

39'bier127', 'brazil58', 'brg180') 

40 

41>>> print(Instance.list_resources(symmetric=False)[0:10]) # only asymmetric 

42('br17', 'ft53', 'ft70', 'ftv170', 'ftv33', 'ftv35', 'ftv38', 'ftv44', \ 

43'ftv47', 'ftv55') 

44 

45You can load any of these instances from the resources via 

46:meth:`Instance.from_resource` as follows: 

47 

48>>> inst = Instance.from_resource("a280") 

49>>> print(inst.n_cities) 

50280 

51 

52If you want to read an instance from an external TSPLib file, you can use 

53:meth:`Instance.from_file`. Be aware that not the whole TSPLib standard is 

54supported right now, but only a reasonably large subset. 

55 

56Every TSP instance automatically provides a lower bound 

57:attr:`Instance.tour_length_lower_bound` and an upper bound 

58:attr:`Instance.tour_length_upper_bound` of the lengths of valid tours. 

59For the TSPLib instances, the globally optimal solutions and their tour 

60lengths are known, so we can use them as lower bounds directly. Otherwise, 

61we currently use a very crude approximation: We assume that, for each city 

62`i`, the next city `j` to be visited would be the nearest neighbor of `i`. 

63Of course, in reality, such a tour may not be feasible, as it could contain 

64disjoint sets of loops. But no tour can be shorter than this. 

65As upper bound, we do the same but assume that `j` would be the cities 

66farthest away from `i`. 

67 

68>>> print(inst.tour_length_lower_bound) 

692579 

70>>> print(inst.tour_length_upper_bound) 

7165406 

72 

73It should be noted that all TSPLib instances by definition have integer 

74distances. This means that there are never floating point distances and, for 

75example, Euclidean distances are rounded and are only approximate Euclidean. 

76Then again, since floating point numbers could also only represent things such 

77as `sqrt(2)` approximately, using integers instead of floating point numbers 

78does not really change anything - distances would be approximately Euclidean 

79or approximately geographical either way. 

80 

81TSPLib also provides some known optimal solutions in path representation, 

82i.e., as permutations of the numbers `0` to `n_cities-1`. The optimal 

83solutions corresponding to the instances provides as resources can be obtained 

84via :mod:`~moptipyapps.tsp.known_optima`. 

85 

86>>> inst = Instance.from_resource("si175") 

87>>> print(inst.n_cities) 

88175 

89>>> int(inst[0, 1]) 

90113 

91>>> int(inst[173, 174]) 

92337 

93>>> int(inst[1, 174]) 

94384 

95>>> int(inst[2, 174]) 

96384 

97>>> int(inst[2, 3]) 

98248 

99>>> int(inst[3, 5]) 

100335 

101>>> int(inst[4, 6]) 

102134 

103 

104The original data of TSPLib can be found at 

105<http://comopt.ifi.uni-heidelberg.de/software/TSPLIB95/>. Before doing 

106anything with these data directly, you should make sure to read the FAQ 

107<http://comopt.ifi.uni-heidelberg.de/software/TSPLIB95/TSPFAQ.html> and the 

108documentation 

109<http://comopt.ifi.uni-heidelberg.de/software/TSPLIB95/tsp95.pdf>. 

110 

111On 2024-10-18, we added a small eleven city instance based on the car travel 

112distances in km between the Chinese cities 

113Shanghai (上海), Beijing (北京), Nanjing (南京), Hefei (合肥), 

114Harbin (哈尔滨), Kunming (昆明), Wuhan (武汉), 

115Xi'an (西安), Chongqing (重庆), Changsha (长沙), 

116and Hong Kong (香港), determined using BaiDu Maps on 2024-10-18. 

117The optimal solution for this instance has length 9547 (km) and visits the 

118cities in the following order: Shanghai (上海), 

119Nanjing (南京), Hefei (合肥), Wuhan (武汉), Changsha (长沙), 

120Hong Kong (香港), Kunming (昆明), Chongqing (重庆), 

121Xi'an (西安), Beijing (北京), Harbin (哈尔滨). 

122This instance can be used to validate and santity-test algorithms, as its 

123solution can easily be determined using exhaustive enumeration. 

124 

125Interestingly, it holds that the calculated travel distance from 

126Harbin (哈尔滨) to Kunming (昆明) is longer than the travel distance 

127from Harbin (哈尔滨) to Xi'an (西安) and then from 

128Xi'an (西安) to Kunming (昆明). 

129 

130- Harbin (哈尔滨) - Kunming (昆明) = 3781 km 

131- Harbin (哈尔滨) - Xi'an (西安) = 2291 km 

132- Xi'an (西安) - Kunming (昆明) = 1483 km 

133- 2291 km + 1483 km = 3774 km, which is less than 3781 km 

134 

135This may be because BaiduMaps ranked the paths by travel time and not travel 

136distance, or by any other strange effect I cannot account for. Maybe it is 

137related to which lane of a road one would naturally drive on or lengths of 

138intersections depending on the direction one is coming from. 

139Either way, this shows that this simple instance `cn11` may already have 

140interesting properties. It violates the triangle equation and it is 

141non-Euclidean. It is also not based on Geo-coordinates, but on actual street 

142distances and times. 

143 

144>>> inst = Instance.from_resource("cn11") 

145>>> int(inst[0, 0]) 

1460 

147>>> int(inst[1, 2]) 

1481007 

149>>> int(inst[1, 3]) 

1501017 

151>>> int(inst[9, 10]) 

152830 

153>>> int(inst[5, 6]) 

1541560 

155 

156Important initial work on this code has been contributed by Mr. Tianyu LIANG 

157(梁天宇), <liangty@stu.hfuu.edu.cn> a Master's student at the Institute of 

158Applied Optimization (应用优化研究所) of the School 

159of Artificial Intelligence and Big Data (人工智能与大数据学院) at Hefei 

160University (合肥大学) in Hefei, Anhui, China (中国安徽省合肥市) under the 

161supervision of Prof. Dr. Thomas Weise (汤卫思教授). 

162 

1631. Gerhard Reinelt. TSPLIB - A Traveling Salesman Problem Library. 

164 *ORSA Journal on Computing* 3(4):376-384. November 1991. 

165 https://doi.org/10.1287/ijoc.3.4.376. 

166 http://comopt.ifi.uni-heidelberg.de/software/TSPLIB95/ 

1672. Gerhard Reinelt. *TSPLIB95.* Heidelberg, Germany: Universität 

168 Heidelberg, Institut für Angewandte Mathematik. 1995. 

169 http://comopt.ifi.uni-heidelberg.de/software/TSPLIB95/tsp95.pdf 

1703. David Lee Applegate, Robert E. Bixby, Vašek Chvátal, and William John Cook. 

171 *The Traveling Salesman Problem: A Computational Study.* Second Edition, 

172 2007. Princeton, NJ, USA: Princeton University Press. Volume 17 of 

173 Princeton Series in Applied Mathematics. ISBN: 0-691-12993-2. 

1744. Gregory Z. Gutin and Abraham P. Punnen. *The Traveling Salesman Problem and 

175 its Variations.* 2002. Kluwer Academic Publishers. Volume 12 of 

176 Combinatorial Optimization. https://doi.org/10.1007/b101971. 

1775. Pedro Larrañaga, Cindy M. H. Kuijpers, Roberto H. Murga, Iñaki Inza, and 

178 Sejla Dizdarevic. Genetic Algorithms for the Travelling Salesman Problem: A 

179 Review of Representations and Operators. *Artificial Intelligence Review* 

180 13(2):129-170. April 1999. https://doi.org/10.1023/A:1006529012972. 

1816. Eugene Leighton Lawler, Jan Karel Lenstra, Alexander Hendrik George Rinnooy 

182 Kan, and David B. Shmoys. *The Traveling Salesman Problem: A Guided Tour of 

183 Combinatorial Optimization.* September 1985. Chichester, West Sussex, UK: 

184 Wiley Interscience. In Estimation, Simulation, and 

185 Control - Wiley-Interscience Series in Discrete Mathematics and 

186 Optimization. ISBN: 0-471-90413-9 

1877. Tianyu Liang, Zhize Wu, Jörg Lässig, Daan van den Berg, and Thomas Weise. 

188 Solving the Traveling Salesperson Problem using Frequency Fitness 

189 Assignment. In *Proceedings of the IEEE Symposium on Foundations of 

190 Computational Intelligence (IEEE FOCI'22),* part of the *IEEE Symposium 

191 Series on Computational Intelligence (SSCI'22),* December 4-7, 2022, 

192 Singapore. Pages 360-367. IEEE. 

193 https://doi.org/10.1109/SSCI51031.2022.10022296. 

1948. Thomas Weise, Raymond Chiong, Ke Tang, Jörg Lässig, Shigeyoshi Tsutsui, 

195 Wenxiang Chen, Zbigniew Michalewicz, and Xin Yao. Benchmarking Optimization 

196 Algorithms: An Open Source Framework for the Traveling Salesman Problem. 

197 *IEEE Computational Intelligence Magazine.* 9(3):40-52. August 2014. 

198 https://doi.org/10.1109/MCI.2014.2326101. 

199""" 

200 

201from math import acos, cos, isfinite, sqrt 

202from string import digits 

203 

204# pylint: disable=W0611 

205from typing import Any, Callable, Final, Iterable, TextIO, cast 

206 

207import moptipy.utils.nputils as npu 

208import numpy as np 

209from moptipy.api.component import Component 

210from moptipy.utils.logger import KeyValueLogSection 

211from moptipy.utils.nputils import int_range_to_dtype 

212from moptipy.utils.strings import sanitize_name 

213from pycommons.io.path import Path, file_path 

214from pycommons.types import check_int_range, check_to_int_range, type_error 

215 

216from moptipyapps.tsp.tsplib import open_resource_stream 

217 

218#: the known optimal tour lengths and lower bounds of TSP Lib 

219_LOWER_BOUNDS: dict[str, int] = { 

220 "a280": 2579, "ali535": 202339, "att48": 10628, "att532": 27686, 

221 "bayg29": 1610, "bays29": 2020, "berlin52": 7542, "bier127": 118282, 

222 "br17": 39, "brazil58": 25395, "brd14051": 469385, "brg180": 1950, 

223 "burma14": 3323, "ch130": 6110, "ch150": 6528, "cn11": 9547, 

224 "d1291": 50801, "d15112": 1573084, "d1655": 62128, "d18512": 645238, 

225 "d198": 15780, "d2103": 80450, "d493": 35002, "d657": 48912, 

226 "dantzig42": 699, "dsj1000": 18660188, "eil101": 629, "eil51": 426, 

227 "eil76": 538, "fl1400": 20127, "fl1577": 22249, "fl3795": 28772, 

228 "fl417": 11861, "fnl4461": 182566, "fri26": 937, "ft53": 6905, 

229 "ft70": 38673, "ftv100": 1788, "ftv110": 1958, "ftv120": 2166, 

230 "ftv130": 2307, "ftv140": 2420, "ftv150": 2611, "ftv160": 2683, 

231 "ftv170": 2755, "ftv33": 1286, "ftv35": 1473, "ftv38": 1530, 

232 "ftv44": 1613, "ftv47": 1776, "ftv55": 1608, "ftv64": 1839, 

233 "ftv70": 1950, "ftv90": 1579, "gil262": 2378, "gr120": 6942, 

234 "gr137": 69853, "gr17": 2085, "gr202": 40160, "gr21": 2707, 

235 "gr229": 134602, "gr24": 1272, "gr431": 171414, "gr48": 5046, 

236 "gr666": 294358, "gr96": 55209, "hk48": 11461, "kro124p": 36230, 

237 "kroA100": 21282, "kroA150": 26524, "kroA200": 29368, "kroB100": 22141, 

238 "kroB150": 26130, "kroB200": 29437, "kroC100": 20749, "kroD100": 21294, 

239 "kroE100": 22068, "lin105": 14379, "lin318": 42029, "nrw1379": 56638, 

240 "p43": 5620, "p654": 34643, "pa561": 2763, "pcb1173": 56892, 

241 "pcb3038": 137694, "pcb442": 50778, "pla33810": 66048945, 

242 "pla7397": 23260728, "pla85900": 142382641, "pr1002": 259045, 

243 "pr107": 44303, "pr124": 59030, "pr136": 96772, "pr144": 58537, 

244 "pr152": 73682, "pr226": 80369, "pr2392": 378032, "pr264": 49135, 

245 "pr299": 48191, "pr439": 107217, "pr76": 108159, "rat195": 2323, 

246 "rat575": 6773, "rat783": 8806, "rat99": 1211, "rbg323": 1326, 

247 "rbg358": 1163, "rbg403": 2465, "rbg443": 2720, "rd100": 7910, 

248 "rd400": 15281, "rl11849": 923288, "rl1304": 252948, "rl1323": 270199, 

249 "rl1889": 316536, "rl5915": 565530, "rl5934": 556045, "ry48p": 14422, 

250 "si1032": 92650, "si175": 21407, "si535": 48450, "st70": 675, 

251 "swiss42": 1273, "ts225": 126643, "tsp225": 3916, "u1060": 224094, 

252 "u1432": 152970, "u159": 42080, "u1817": 57201, "u2152": 64253, 

253 "u2319": 234256, "u574": 36905, "u724": 41910, "ulysses16": 6859, 

254 "ulysses22": 7013, "usa13509": 19982859, "vm1084": 239297, 

255 "vm1748": 336556} 

256 

257#: the TSPLib instances 

258_SYMMETRIC_INSTANCES: Final[tuple[str, ...]] = ( 

259 "a280", "ali535", "att48", "att532", "bayg29", "bays29", "berlin52", 

260 "bier127", "brazil58", "brg180", "burma14", "ch130", "ch150", "cn11", 

261 "d1291", "d1655", "d198", "d493", "d657", "dantzig42", "eil101", "eil51", 

262 "eil76", "fl1400", "fl1577", "fl417", "fri26", "gil262", "gr120", "gr137", 

263 "gr17", "gr202", "gr21", "gr229", "gr24", "gr431", "gr48", "gr666", 

264 "gr96", "hk48", "kroA100", "kroA150", "kroA200", "kroB100", "kroB150", 

265 "kroB200", "kroC100", "kroD100", "kroE100", "lin105", "lin318", "nrw1379", 

266 "p654", "pcb1173", "pcb442", "pr1002", "pr107", "pr124", "pr136", "pr144", 

267 "pr152", "pr226", "pr264", "pr299", "pr439", "pr76", "rat195", "rat575", 

268 "rat783", "rat99", "rd100", "rd400", "rl1304", "rl1323", "rl1889", 

269 "si1032", "si175", "si535", "st70", "swiss42", "ts225", "tsp225", "u1060", 

270 "u1432", "u159", "u1817", "u574", "u724", "ulysses16", "ulysses22", 

271 "vm1084", "vm1748") 

272 

273#: The set of asymmetric instances 

274_ASYMMETRIC_INSTANCES: Final[tuple[str, ...]] = ( 

275 "br17", "ft53", "ft70", "ftv170", "ftv33", "ftv35", "ftv38", "ftv44", 

276 "ftv47", "ftv55", "ftv64", "ftv70", "kro124p", "p43", "rbg323", "rbg358", 

277 "rbg403", "rbg443", "ry48p") 

278 

279#: The set of all TSP instances 

280_INSTANCES: Final[tuple[str, ...]] = tuple(sorted( 

281 _SYMMETRIC_INSTANCES + _ASYMMETRIC_INSTANCES)) 

282 

283#: the set of asymmetric resources 

284_ASYMMETRIC_RESOURCES: Final[set[str]] = set(_ASYMMETRIC_INSTANCES) 

285 

286#: the problem is a symmetric tsp 

287_TYPE_SYMMETRIC_TSP: Final[str] = "TSP" 

288#: the problem is an asymmetric tsp 

289_TYPE_ASYMMETRIC_TSP: Final[str] = "ATSP" 

290#: the permitted types 

291_TYPES: Final[set[str]] = {_TYPE_SYMMETRIC_TSP, _TYPE_ASYMMETRIC_TSP} 

292#: the name start 

293_KEY_NAME: Final[str] = "NAME" 

294#: the type start 

295_KEY_TYPE: Final[str] = "TYPE" 

296#: the dimension start 

297_KEY_DIMENSION: Final[str] = "DIMENSION" 

298#: the comment start 

299_KEY_COMMENT: Final[str] = "COMMENT" 

300#: EUC_2D coordinates 

301__EWT_EUC_2D: Final[str] = "EUC_2D" 

302#: geographical coordinates 

303__EWT_GEO: Final[str] = "GEO" 

304#: ATT coordinates 

305__EWT_ATT: Final[str] = "ATT" 

306#: ceiling 2D coordinates 

307__EWT_CEIL2D: Final[str] = "CEIL_2D" 

308#: the explicit edge weight type 

309_EWT_EXPLICIT: Final[str] = "EXPLICIT" 

310#: the edge weight type start 

311_KEY_EDGE_WEIGHT_TYPE: Final[str] = "EDGE_WEIGHT_TYPE" 

312#: the permitted edge weight types 

313_EDGE_WEIGHT_TYPES: Final[set[str]] = { 

314 __EWT_EUC_2D, __EWT_GEO, __EWT_ATT, __EWT_CEIL2D, _EWT_EXPLICIT} 

315#: the edge weight format "function" 

316__EWF_FUNCTION: Final[str] = "FUNCTION" 

317#: the full matrix edge weight format 

318_EWF_FULL_MATRIX: Final[str] = "FULL_MATRIX" 

319#: the upper row edge weight format 

320_EWF_UPPER_ROW: Final[str] = "UPPER_ROW" 

321#: the lower diagonal row 

322__EWF_LOWER_DIAG_ROW: Final[str] = "LOWER_DIAG_ROW" 

323#: the upper diagonal row 

324__EWF_UPPER_DIAG_ROW: Final[str] = "UPPER_DIAG_ROW" 

325#: the edge weight format start 

326_KEY_EDGE_WEIGHT_FORMAT: Final[str] = "EDGE_WEIGHT_FORMAT" 

327#: the permitted edge weight formats 

328_EDGE_WEIGHT_FORMATS: Final[set[str]] = { 

329 __EWF_FUNCTION, _EWF_FULL_MATRIX, _EWF_UPPER_ROW, 

330 __EWF_LOWER_DIAG_ROW, __EWF_UPPER_DIAG_ROW} 

331#: the start of the node coord type 

332_KEY_NODE_COORD_TYPE: Final[str] = "NODE_COORD_TYPE" 

333#: 2d coordinates 

334__NODE_COORD_TYPE_2D: Final[str] = "TWOD_COORDS" 

335#: no coordinates 

336__NODE_COORD_TYPE_NONE: Final[str] = "NO_COORDS" 

337#: the permitted node coordinate types 

338_NODE_COORD_TYPES: Final[set[str]] = { 

339 __NODE_COORD_TYPE_2D, __NODE_COORD_TYPE_NONE} 

340#: the node coordinate section starts 

341_START_NODE_COORD_SECTION: Final[str] = "NODE_COORD_SECTION" 

342#: start the edge weight section 

343_START_EDGE_WEIGHT_SECTION: Final[str] = "EDGE_WEIGHT_SECTION" 

344#: the end of the file 

345_EOF: Final[str] = "EOF" 

346#: the fixed edges section 

347_FIXED_EDGES: Final[str] = "FIXED_EDGES_SECTION" 

348 

349 

350def __line_to_nums(line: str, 

351 collector: Callable[[int | float], Any]) -> None: 

352 """ 

353 Split a line to a list of integers. 

354 

355 :param line: the line string 

356 :param collector: the collector 

357 """ 

358 if not isinstance(line, str): 

359 raise type_error(line, "line", str) 

360 idx: int = 0 

361 line = line.strip() 

362 str_len: Final[int] = len(line) 

363 

364 while idx < str_len: 

365 while (idx < str_len) and (line[idx] == " "): 

366 idx += 1 

367 if idx >= str_len: 

368 return 

369 

370 next_space = line.find(" ", idx) 

371 if next_space < 0: 

372 next_space = str_len 

373 

374 part: str = line[idx:next_space] 

375 if ("." in part) or ("E" in part) or ("e" in part): 

376 f: float = float(part) 

377 if not isfinite(f): 

378 raise ValueError( 

379 f"{part!r} translates to non-finite float {f}.") 

380 collector(f) 

381 else: 

382 collector(check_to_int_range( 

383 part, "line fragment", -1_000_000_000_000, 1_000_000_000_000)) 

384 idx = next_space 

385 

386 

387def __nint(v: int | float) -> int: 

388 """ 

389 Get the nearest integer in the TSPLIB format. 

390 

391 :param v: the value 

392 :return: the integer 

393 """ 

394 if isinstance(v, int): 

395 return v 

396 if isinstance(v, float): 

397 return int(0.5 + v) 

398 raise type_error(v, "value", (int, float)) 

399 

400 

401def __dist_2deuc(a: list[int | float], b: list[int | float]) -> int: 

402 """ 

403 Compute the two-dimensional Euclidean distance function. 

404 

405 :param a: the first point 

406 :param b: the second point 

407 :return: the distance 

408 """ 

409 return __nint(sqrt(((a[0] - b[0]) ** 2) + ((a[1] - b[1]) ** 2))) 

410 

411 

412def __matrix_from_points( 

413 n_cities: int, coord_dim: int, stream: Iterable[str], 

414 dist_func: Callable[[list[int | float], list[int | float]], int]) \ 

415 -> np.ndarray: 

416 """ 

417 Load a distance matrix from 2D coordinates. 

418 

419 :param n_cities: the dimension 

420 :param coord_dim: the coordinate dimension 

421 :param stream: the stream 

422 :param dist_func: the distance function 

423 :return: the matrix 

424 """ 

425 index: int = 0 

426 coordinates: Final[list[list[int | float]]] = [] 

427 row: Final[list[int | float]] = [] 

428 for the_line in stream: 

429 if not isinstance(the_line, str): 

430 raise type_error(the_line, "line", str) 

431 line = the_line.strip() 

432 if len(line) <= 0: 

433 continue 

434 if line == _EOF: 

435 break 

436 __line_to_nums(line, row.append) 

437 index += 1 

438 if (len(row) != (coord_dim + 1)) or (not isinstance(row[0], int)) \ 

439 or (row[0] != index): 

440 raise ValueError(f"invalid row {line!r} at index {index}, " 

441 f"gives values {row}.") 

442 coordinates.append(row[1:]) 

443 row.clear() 

444 if index != n_cities: 

445 raise ValueError(f"only found {index} rows, but expected {n_cities}") 

446 

447 # now construct the matrix 

448 matrix: Final[np.ndarray] = np.zeros((n_cities, n_cities), 

449 npu.DEFAULT_INT) 

450 for i in range(n_cities): 

451 a = coordinates[i] 

452 for j in range(i): 

453 b = coordinates[j] 

454 dist = dist_func(a, b) 

455 if not isinstance(dist, int): 

456 raise type_error(dist, f"distance[{i},{j}]", int) 

457 if dist < 0: 

458 raise ValueError( 

459 f"Cannot have node distance {dist}: ({a}) at index=" 

460 f"{i + 1} and ({b} at index={j + 1}.") 

461 matrix[i, j] = dist 

462 matrix[j, i] = dist 

463 return matrix 

464 

465 

466def __coord_to_rad(x: int | float) -> float: 

467 """ 

468 Convert coordinates to longitude/latitude. 

469 

470 :param x: the coordinate 

471 """ 

472 degrees: int = int(x) 

473 return (3.141592 * (degrees + (5.0 * ( # noqa: FURB152 

474 x - degrees)) / 3.0)) / 180.0 

475 

476 

477def __dist_loglat(a: list[int | float], b: list[int | float]) -> int: 

478 """ 

479 Convert a longitude-latitude pair to a distance. 

480 

481 :param a: the first point 

482 :param b: the second point 

483 :return: the distance 

484 """ 

485 lat1: float = __coord_to_rad(a[0]) 

486 long1: float = __coord_to_rad(a[1]) 

487 lat2: float = __coord_to_rad(b[0]) 

488 long2: float = __coord_to_rad(b[1]) 

489 q1: Final[float] = cos(long1 - long2) 

490 q2: Final[float] = cos(lat1 - lat2) 

491 q3: Final[float] = cos(lat1 + lat2) 

492 return int(6378.388 * acos(0.5 * ( 

493 (1.0 + q1) * q2 - (1.0 - q1) * q3)) + 1.0) 

494 

495 

496def __dist_att(a: list[int | float], b: list[int | float]) -> int: 

497 """ 

498 Compute the ATT Pseudo-Euclidean distance function. 

499 

500 :param a: the first point 

501 :param b: the second point 

502 :return: the distance 

503 """ 

504 xd: Final[int | float] = a[0] - b[0] 

505 yd: Final[int | float] = a[1] - b[1] 

506 rij: Final[float] = sqrt((xd * xd + yd * yd) / 10.0) 

507 tij: Final[int] = __nint(rij) 

508 return (tij + 1) if tij < rij else tij 

509 

510 

511def __dist_2dceil(a: list[int | float], b: list[int | float]) -> int: 

512 """ 

513 Compute the ceiling of the two-dimensional Euclidean distance function. 

514 

515 :param a: the first point 

516 :param b: the second point 

517 :return: the distance 

518 """ 

519 dist: Final[float] = sqrt(((a[0] - b[0]) ** 2) + ((a[1] - b[1]) ** 2)) 

520 disti: Final[int] = int(dist) 

521 return disti if dist == disti else (disti + 1) 

522 

523 

524def _matrix_from_node_coord_section( 

525 n_cities: int | None, edge_weight_type: str | None, 

526 node_coord_type: str | None, stream: Iterable[str]) -> np.ndarray: 

527 """ 

528 Get a data matrix from a node coordinate section. 

529 

530 :param n_cities: the dimension 

531 :param edge_weight_type: the edge weight type 

532 :param node_coord_type: the node coordinate type 

533 :param stream: the node coordinate stream 

534 :return: the data matrix 

535 """ 

536 check_int_range(n_cities, "n_cities", 2, 1_000_000_000_000) 

537 if (node_coord_type is not None) and \ 

538 (not isinstance(node_coord_type, str)): 

539 raise type_error(node_coord_type, "node_coord_type", str) 

540 if not isinstance(edge_weight_type, str): 

541 raise type_error(edge_weight_type, "edge_weight_type", str) 

542 

543 dist_fun = None 

544 coord_dim: int | None = None 

545 if (edge_weight_type == __EWT_EUC_2D) \ 

546 and (node_coord_type in {None, __NODE_COORD_TYPE_2D}): 

547 coord_dim = 2 

548 dist_fun = __dist_2deuc 

549 elif (edge_weight_type == __EWT_GEO) \ 

550 and (node_coord_type in {None, __NODE_COORD_TYPE_2D}): 

551 coord_dim = 2 

552 dist_fun = __dist_loglat 

553 elif (edge_weight_type == __EWT_ATT) \ 

554 and (node_coord_type in {None, __NODE_COORD_TYPE_2D}): 

555 coord_dim = 2 

556 dist_fun = __dist_att 

557 elif (edge_weight_type == __EWT_CEIL2D) \ 

558 and (node_coord_type in {None, __NODE_COORD_TYPE_2D}): 

559 coord_dim = 2 

560 dist_fun = __dist_2dceil 

561 

562 if (coord_dim is not None) and (dist_fun is not None): 

563 return __matrix_from_points(n_cities, coord_dim, stream, dist_fun) 

564 

565 raise ValueError(f"invalid combination of {_KEY_EDGE_WEIGHT_TYPE} " 

566 f"and {_KEY_NODE_COORD_TYPE}") 

567 

568 

569def __read_n_ints(n: int, stream: Iterable[str]) -> list[int]: 

570 """ 

571 Read exactly `n` integers from a stream. 

572 

573 :param n: the number of integers to read 

574 :param stream: the stream to read from 

575 :return: the list of integers 

576 """ 

577 res: list[int] = [] 

578 

579 def __append(i: int | float, fwd=res.append) -> None: 

580 if isinstance(i, int): 

581 fwd(i) 

582 else: 

583 i2 = int(i) 

584 if i2 != i: 

585 raise ValueError(f"{i} is not integer") 

586 fwd(i2) 

587 

588 for line in stream: 

589 __line_to_nums(line, __append) 

590 if len(res) == n: 

591 break 

592 

593 if len(res) != n: 

594 raise ValueError(f"expected {n} integers, got {len(res)}.") 

595 return res 

596 

597 

598def _matrix_from_edge_weights( 

599 n_cities: int | None, edge_weight_type: str | None, 

600 edge_weight_format: str | None, stream: Iterable[str]) -> np.ndarray: 

601 """ 

602 Get a data matrix from a edge weights section. 

603 

604 :param n_cities: the dimension 

605 :param edge_weight_type: the edge weight type 

606 :param edge_weight_format: the edge weight format 

607 :param stream: the node coordinate stream 

608 :return: the data matrix 

609 """ 

610 check_int_range(n_cities, "n_cities", 2, 1_000_000_000_000) 

611 if not isinstance(edge_weight_type, str): 

612 raise type_error(edge_weight_type, "node_coord_type", str) 

613 if not isinstance(edge_weight_type, str): 

614 raise type_error(edge_weight_type, "edge_weight_type", str) 

615 if edge_weight_type == _EWT_EXPLICIT: 

616 if edge_weight_format == _EWF_FULL_MATRIX: 

617 res = np.array(__read_n_ints(n_cities * n_cities, stream), 

618 dtype=npu.DEFAULT_INT).reshape( 

619 (n_cities, n_cities)) 

620 np.fill_diagonal(res, 0) 

621 return res 

622 if edge_weight_format == _EWF_UPPER_ROW: 

623 ints = __read_n_ints((n_cities * (n_cities - 1)) // 2, stream) 

624 res = np.zeros((n_cities, n_cities), dtype=npu.DEFAULT_INT) 

625 i: int = 1 

626 j: int = 0 

627 for v in ints: 

628 res[j, i] = res[i, j] = v 

629 i += 1 

630 if i >= n_cities: 

631 j += 1 

632 i = j + 1 

633 return res 

634 if edge_weight_format == __EWF_LOWER_DIAG_ROW: 

635 ints = __read_n_ints( 

636 n_cities + ((n_cities * (n_cities - 1)) // 2), stream) 

637 res = np.zeros((n_cities, n_cities), dtype=npu.DEFAULT_INT) 

638 i = 0 

639 j = 0 

640 for v in ints: 

641 if i != j: 

642 res[j, i] = res[i, j] = v 

643 i += 1 

644 if i > j: 

645 j += 1 

646 i = 0 

647 return res 

648 if edge_weight_format == __EWF_UPPER_DIAG_ROW: 

649 ints = __read_n_ints( 

650 n_cities + ((n_cities * (n_cities - 1)) // 2), stream) 

651 res = np.zeros((n_cities, n_cities), dtype=npu.DEFAULT_INT) 

652 i = 0 

653 j = 0 

654 for v in ints: 

655 if i != j: 

656 res[j, i] = res[i, j] = v 

657 i += 1 

658 if i >= n_cities: 

659 j += 1 

660 i = j 

661 return res 

662 raise ValueError( 

663 f"unsupported combination of {_KEY_EDGE_WEIGHT_TYPE}=" 

664 f"{edge_weight_type!r} and {_KEY_EDGE_WEIGHT_FORMAT}=" 

665 f"{edge_weight_format!r}") 

666 

667 

668def _from_stream( 

669 stream: Iterable[str], 

670 lower_bound_getter: Callable[[str], int] | None = 

671 _LOWER_BOUNDS.__getitem__) -> "Instance": 

672 """ 

673 Read a TSP Lib instance from a TSP-lib formatted stream. 

674 

675 :param stream: the text stream 

676 :param lower_bound_getter: a function returning a lower bound for an 

677 instance name, or `None` to use a simple lower bound approximation 

678 :return: the instance 

679 """ 

680 if (lower_bound_getter is not None) \ 

681 and (not callable(lower_bound_getter)): 

682 raise type_error( 

683 lower_bound_getter, "lower_bound_getter", None, call=True) 

684 

685 the_name: str | None = None 

686 the_type: str | None = None 

687 the_n_cities: int | None = None 

688 the_ewt: str | None = None 

689 the_ewf: str | None = None 

690 the_nct: str | None = None 

691 the_matrix: np.ndarray | None = None 

692 

693 for the_line in stream: 

694 if not isinstance(the_line, str): 

695 raise type_error(the_line, "line", str) 

696 line = the_line.strip() 

697 if len(line) <= 0: 

698 continue 

699 

700 sep_idx: int = line.find(":") 

701 if sep_idx > 0: 

702 key: str = line[:sep_idx].strip() 

703 value: str = line[sep_idx + 1:].strip() 

704 if len(value) <= 0: 

705 raise ValueError(f"{line!r} has empty value " 

706 f"{value!r} for key {key!r}.") 

707 

708 if key == _KEY_NAME: 

709 if the_name is not None: 

710 raise ValueError( 

711 f"{line!r} cannot come at this position," 

712 f" already got {_KEY_NAME}={the_name!r}.") 

713 if value.endswith(".tsp"): 

714 value = value[0:-4] 

715 the_name = value 

716 continue 

717 

718 if key == _KEY_TYPE: 

719 if the_type is not None: 

720 raise ValueError( 

721 f"{line!r} cannot come at this position, already " 

722 f"got {_KEY_TYPE}={the_type!r}.") 

723 the_type = _TYPE_SYMMETRIC_TSP \ 

724 if value == "TSP (M.~Hofmeister)" else value 

725 if the_type not in _TYPES: 

726 raise ValueError( 

727 f"only {_TYPES!r} are permitted as {_KEY_TYPE}, " 

728 f"but got {the_type!r} from {line!r}.") 

729 continue 

730 

731 if key == _KEY_DIMENSION: 

732 if the_n_cities is not None: 

733 raise ValueError( 

734 f"{line!r} cannot come at this position," 

735 f" already got {_KEY_DIMENSION}={the_n_cities}.") 

736 the_n_cities = check_to_int_range(value, "dimension", 

737 2, 1_000_000_000) 

738 

739 if key == _KEY_EDGE_WEIGHT_TYPE: 

740 if the_ewt is not None: 

741 raise ValueError( 

742 f"{line!r} cannot come at this position, already " 

743 f"got {_KEY_EDGE_WEIGHT_TYPE}={the_ewt!r}.") 

744 the_ewt = value 

745 if the_ewt not in _EDGE_WEIGHT_TYPES: 

746 raise ValueError( 

747 f"only {_EDGE_WEIGHT_TYPES!r} are permitted as " 

748 f"{_KEY_EDGE_WEIGHT_TYPE}, but got {the_ewt!r} " 

749 f"in {line!r}") 

750 continue 

751 

752 if key == _KEY_EDGE_WEIGHT_FORMAT: 

753 if the_ewf is not None: 

754 raise ValueError( 

755 f"{line!r} cannot come at this position, already " 

756 f"got {_KEY_EDGE_WEIGHT_FORMAT}={the_ewf!r}.") 

757 the_ewf = value 

758 if the_ewf not in _EDGE_WEIGHT_FORMATS: 

759 raise ValueError( 

760 f"only {_EDGE_WEIGHT_FORMATS!r} are permitted as " 

761 f"{_KEY_EDGE_WEIGHT_FORMAT}, but got {the_ewf} " 

762 f"in {line!r}") 

763 continue 

764 if key == _KEY_NODE_COORD_TYPE: 

765 if the_nct is not None: 

766 raise ValueError( 

767 f"{line!r} cannot come at this position, already " 

768 f"got {_KEY_NODE_COORD_TYPE}={the_nct!r}.") 

769 the_nct = value 

770 if the_nct not in _NODE_COORD_TYPES: 

771 raise ValueError( 

772 f"only {_NODE_COORD_TYPES!r} are permitted as node " 

773 f"{_KEY_NODE_COORD_TYPE}, but got {the_nct!r} " 

774 f"in {line!r}") 

775 continue 

776 elif line == _START_NODE_COORD_SECTION: 

777 if the_matrix is not None: 

778 raise ValueError( 

779 f"already got matrix, cannot have {line!r} here!") 

780 the_matrix = _matrix_from_node_coord_section( 

781 the_n_cities, the_ewt, the_nct, stream) 

782 continue 

783 elif line == _START_EDGE_WEIGHT_SECTION: 

784 if the_matrix is not None: 

785 raise ValueError( 

786 f"already got matrix, cannot have {line!r} here!")