diff options
Diffstat (limited to 'glpk-5.0/examples/tsp')
-rw-r--r-- | glpk-5.0/examples/tsp/README | 37 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/bench.txt | 56 | ||||
-rwxr-xr-x | glpk-5.0/examples/tsp/build.sh | 8 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/dantzig42.tsp | 103 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/gr120.tsp | 534 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/main.c | 535 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/maxflow.c | 170 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/maxflow.h | 20 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/mincut.c | 392 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/mincut.h | 23 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/misc.c | 159 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/misc.h | 24 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/moscow.tsp | 200 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/sample.tsp | 16 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/tsplib.c | 730 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/tsplib.h | 80 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/ulysses16.tsp | 24 | ||||
-rw-r--r-- | glpk-5.0/examples/tsp/ulysses22.tsp | 30 |
18 files changed, 3141 insertions, 0 deletions
diff --git a/glpk-5.0/examples/tsp/README b/glpk-5.0/examples/tsp/README new file mode 100644 index 0000000..7d497e5 --- /dev/null +++ b/glpk-5.0/examples/tsp/README @@ -0,0 +1,37 @@ +This subdirectory contains an example application program, TSPSOL, +which is a stand-alone solver intended to solve the Symmetric Traveling +Salesman Problem (TSP) with the GLPK integer optimizer. + +Please note that this program is only an illustrative example that +illustrates generating "lazy" constraints during the branch-and-bound +search. It is *not* a state-of-the-art code, so only small-sized TSP +instances (perhaps, having up to 150-200 nodes) can be solved with this +program in a reasonable time. For more details see comments in the +source code. + +To build TSPSOL executable you need to run 'build.sh' script. Note that +you should have the GLPK library properly installed. + +To run the TSPSOL program use the following command: + + tspsol tsp-file + +where tsp-file specifies an input text file containing TSP data in +TSPLIB 95 format. + +Detailed description of the input format recognized by TSPSOL is given +in the report: Gerhard Reinelt, "TSPLIB 95". This report as well as +TSPLIB, a library of sample TSP instances, are freely available for +research purposes; see: +<http://www.iwr.uni-heidelberg.de/groups/comopt/software/TSPLIB95/>. + +This subdirectory also includes the following example TSP instances: + +dantzig42.tsp 42 cities (Dantzig) [from TSPLIB] +gr120.tsp 120 cities in Germany (Groetschel) [from TSPLIB] +moscow.tsp 68 cities in Moscow region (Makhorin) +sample.tsp small example from D.Phillips and A.Garcia-Diaz +ulysses16.tsp Odyssey of Ulysses (Groetschel/Padberg) [from TSPLIB] +ulysses22.tsp Odyssey of Ulysses (Groetschel/Padberg) [from TSPLIB] + +Please send comments to the <help-glpk@gnu.org> mailing list. diff --git a/glpk-5.0/examples/tsp/bench.txt b/glpk-5.0/examples/tsp/bench.txt new file mode 100644 index 0000000..4596b2e --- /dev/null +++ b/glpk-5.0/examples/tsp/bench.txt @@ -0,0 +1,56 @@ +Solver: TSPSOL for GLPK 4.56 +Computer: Intel Celeron J1800 2.41 GHz +OS: Debian GNU/Linux 8.1.0 "Jessie" +Compiler: GCC 4.7.2 (options used: -O2) +Test set: TSPLIB 95 <http://www.iwr.uni-heidelberg.de/groups/comopt/software/TSPLIB95/> + +Instance N Solution Lower Bound Nodes Iters Time,s Mem,MB +------------ --- ------------ ------------ -------- ------- ------ ------ +att48 48 10628 opt 1 336 < 1 1.2 +bayg29 29 1610 opt 1 173 < 1 0.3 +bays29 29 2020 opt 1 166 < 1 0.3 +berlin52 52 7542 opt 1 253 < 1 0.7 +bier127 127 118282 opt 29 1330 14 19.1 +brazil58 58 25395 opt 1 458 1 2.0 +brg180 180 1950 opt 131 20012 83 52.0 +burma14 14 3323 opt 1 55 < 1 0.1 +ch130 130 6110 opt 45 2212 38 24.2 +ch150 150 6528 opt 271 4967 138 27.5 +d198 98 15780 opt 259 11371 719 92.3 +dantzig42 42 699 opt 1 171 < 1 0.8 +eil51 51 426 opt 115 1368 2 2.4 +eil76 76 538 opt 1 517 < 1 2.2 +eil101 101 629 opt 1 838 4 9.9 +fri26 26 937 opt 1 125 < 1 0.2 +gr17 17 2085 opt 1 93 < 1 0.1 +gr21 21 2707 opt 1 82 < 1 0.1 +gr24 24 1272 opt 1 137 < 1 0.2 +gr48 48 5046 opt 3 407 1 2.3 +gr96 96 55209 opt 367 5564 63 12.4 +gr120 120 6942 opt 121 2940 46 14.7 +gr137 137 69853 opt 97 1934 27 16.2 +gr202 202 40160 opt 183 4176 287 88.4 +hk48 48 11461 opt 1 322 < 1 1.1 +kroA100 100 21282 opt 57 2227 23 13.2 +kroB100 100 22141 opt 71 1891 27 15.6 +kroC100 100 20749 opt 9 1035 5 9.4 +kroD100 100 21294 opt 9 1203 10 12.3 +kroE100 100 22068 opt 323 5055 100 13.4 +kroA150 150 26524 opt 431 8069 487 50.2 +kroB150 150 26130 opt 309 12599 610 43.1 +lin105 105 14379 opt 5 910 4 7.3 +lin318 318 42029 opt 1527 25885 4972 286.4 +pr76 76 108159 opt 20537 636616 9218 128.7 * +pr107 107 44303 opt 1 3692 38 33.2 +pr124 124 59030 opt 19 1306 23 25.6 +pr136 136 96772 opt 93 2799 60 27.7 +pr144 144 58537 opt 11 1948 37 32.2 +pr152 152 73682 opt 125 3269 99 46.0 +pr226 226 80369 opt 1 6699 240 163.8 +rat99 99 1211 opt 11 467 2 5,6 +rd100 100 7910 opt 1 868 3 7.2 +st70 70 675 opt 1 688 1 3.2 +swiss42 42 1273 opt 1 231 < 1 0.8 +u159 159 42080 opt 9 1356 15 30.0 +ulysses16 16 6859 opt 1 72 < 1 0.1 +ulysses22 22 7013 opt 1 110 < 1 0.2 diff --git a/glpk-5.0/examples/tsp/build.sh b/glpk-5.0/examples/tsp/build.sh new file mode 100755 index 0000000..ab6abc1 --- /dev/null +++ b/glpk-5.0/examples/tsp/build.sh @@ -0,0 +1,8 @@ +#!/bin/sh +# +# Run this script to build TSPSOL executable. +# +# NOTE: you need to have GLPK properly installed. +# +gcc -O2 -otspsol main.c maxflow.c mincut.c misc.c tsplib.c -lglpk -lm +./tspsol sample.tsp diff --git a/glpk-5.0/examples/tsp/dantzig42.tsp b/glpk-5.0/examples/tsp/dantzig42.tsp new file mode 100644 index 0000000..15567bf --- /dev/null +++ b/glpk-5.0/examples/tsp/dantzig42.tsp @@ -0,0 +1,103 @@ +NAME : dantzig42 +TYPE : TSP +COMMENT : 42 cities (Dantzig) +DIMENSION : 42 +EDGE_WEIGHT_TYPE : EXPLICIT +EDGE_WEIGHT_FORMAT : LOWER_DIAG_ROW +DISPLAY_DATA_TYPE : TWOD_DISPLAY +EDGE_WEIGHT_SECTION + 0 8 0 39 45 0 37 47 9 0 50 49 21 15 0 61 62 21 + 20 17 0 58 60 16 17 18 6 0 59 60 15 20 26 17 10 0 + 62 66 20 25 31 22 15 5 0 81 81 40 44 50 41 35 24 20 + 0 103 107 62 67 72 63 57 46 41 23 0 108 117 66 71 77 68 + 61 51 46 26 11 0 145 149 104 108 114 106 99 88 84 63 49 40 + 0 181 185 140 144 150 142 135 124 120 99 85 76 35 0 187 191 146 + 150 156 142 137 130 125 105 90 81 41 10 0 161 170 120 124 130 115 + 110 104 105 90 72 62 34 31 27 0 142 146 101 104 111 97 91 85 + 86 75 51 59 29 53 48 21 0 174 178 133 138 143 129 123 117 118 + 107 83 84 54 46 35 26 31 0 185 186 142 143 140 130 126 124 128 + 118 93 101 72 69 58 58 43 26 0 164 165 120 123 124 106 106 105 + 110 104 86 97 71 93 82 62 42 45 22 0 137 139 94 96 94 80 + 78 77 84 77 56 64 65 90 87 58 36 68 50 30 0 117 122 77 + 80 83 68 62 60 61 50 34 42 49 82 77 60 30 62 70 49 21 + 0 114 118 73 78 84 69 63 57 59 48 28 36 43 77 72 45 27 + 59 69 55 27 5 0 85 89 44 48 53 41 34 28 29 22 23 35 + 69 105 102 74 56 88 99 81 54 32 29 0 77 80 36 40 46 34 + 27 19 21 14 29 40 77 114 111 84 64 96 107 87 60 40 37 8 + 0 87 89 44 46 46 30 28 29 32 27 36 47 78 116 112 84 66 + 98 95 75 47 36 39 12 11 0 91 93 48 50 48 34 32 33 36 + 30 34 45 77 115 110 83 63 97 91 72 44 32 36 9 15 3 0 + 105 106 62 63 64 47 46 49 54 48 46 59 85 119 115 88 66 98 + 79 59 31 36 42 28 33 21 20 0 111 113 69 71 66 51 53 56 + 61 57 59 71 96 130 126 98 75 98 85 62 38 47 53 39 42 29 + 30 12 0 91 92 50 51 46 30 34 38 43 49 60 71 103 141 136 + 109 90 115 99 81 53 61 62 36 34 24 28 20 20 0 83 85 42 + 43 38 22 26 32 36 51 63 75 106 142 140 112 93 126 108 88 60 + 64 66 39 36 27 31 28 28 8 0 89 91 55 55 50 34 39 44 + 49 63 76 87 120 155 150 123 100 123 109 86 62 71 78 52 49 39 + 44 35 24 15 12 0 95 97 64 63 56 42 49 56 60 75 86 97 + 126 160 155 128 104 128 113 90 67 76 82 62 59 49 53 40 29 25 + 23 11 0 74 81 44 43 35 23 30 39 44 62 78 89 121 159 155 + 127 108 136 124 101 75 79 81 54 50 42 46 43 39 23 14 14 21 + 0 67 69 42 41 31 25 32 41 46 64 83 90 130 164 160 133 114 + 146 134 111 85 84 86 59 52 47 51 53 49 32 24 24 30 9 0 + 74 76 61 60 42 44 51 60 66 83 102 110 147 185 179 155 133 159 + 146 122 98 105 107 79 71 66 70 70 60 48 40 36 33 25 18 0 + 57 59 46 41 25 30 36 47 52 71 93 98 136 172 172 148 126 158 + 147 124 121 97 99 71 65 59 63 67 62 46 38 37 43 23 13 17 + 0 45 46 41 34 20 34 38 48 53 73 96 99 137 176 178 151 131 + 163 159 135 108 102 103 73 67 64 69 75 72 54 46 49 54 34 24 + 29 12 0 35 37 35 26 18 34 36 46 51 70 93 97 134 171 176 + 151 129 161 163 139 118 102 101 71 65 65 70 84 78 58 50 56 62 + 41 32 38 21 9 0 29 33 30 21 18 35 33 40 45 65 87 91 + 117 166 171 144 125 157 156 139 113 95 97 67 60 62 67 79 82 62 + 53 59 66 45 38 45 27 15 6 0 3 11 41 37 47 57 55 58 + 63 83 105 109 147 186 188 164 144 176 182 161 134 119 116 86 78 84 + 88 101 108 88 80 86 92 71 64 71 54 41 32 25 0 5 12 55 + 41 53 64 61 61 66 84 111 113 150 186 192 166 147 180 188 167 140 + 124 119 90 87 90 94 107 114 77 86 92 98 80 74 77 60 48 38 + 32 6 0 +DISPLAY_DATA_SECTION + 1 170.0 85.0 + 2 166.0 88.0 + 3 133.0 73.0 + 4 140.0 70.0 + 5 142.0 55.0 + 6 126.0 53.0 + 7 125.0 60.0 + 8 119.0 68.0 + 9 117.0 74.0 + 10 99.0 83.0 + 11 73.0 79.0 + 12 72.0 91.0 + 13 37.0 94.0 + 14 6.0 106.0 + 15 3.0 97.0 + 16 21.0 82.0 + 17 33.0 67.0 + 18 4.0 66.0 + 19 3.0 42.0 + 20 27.0 33.0 + 21 52.0 41.0 + 22 57.0 59.0 + 23 58.0 66.0 + 24 88.0 65.0 + 25 99.0 67.0 + 26 95.0 55.0 + 27 89.0 55.0 + 28 83.0 38.0 + 29 85.0 25.0 + 30 104.0 35.0 + 31 112.0 37.0 + 32 112.0 24.0 + 33 113.0 13.0 + 34 125.0 30.0 + 35 135.0 32.0 + 36 147.0 18.0 + 37 147.5 36.0 + 38 154.5 45.0 + 39 157.0 54.0 + 40 158.0 61.0 + 41 172.0 82.0 + 42 174.0 87.0 +EOF diff --git a/glpk-5.0/examples/tsp/gr120.tsp b/glpk-5.0/examples/tsp/gr120.tsp new file mode 100644 index 0000000..6322489 --- /dev/null +++ b/glpk-5.0/examples/tsp/gr120.tsp @@ -0,0 +1,534 @@ +NAME: gr120 +TYPE: TSP +COMMENT: 120 cities in Germany (Groetschel) +DIMENSION: 120 +EDGE_WEIGHT_TYPE: EXPLICIT +EDGE_WEIGHT_FORMAT: LOWER_DIAG_ROW +DISPLAY_DATA_TYPE: TWOD_DISPLAY +EDGE_WEIGHT_SECTION + 0 534 0 434 107 0 294 241 148 0 593 190 137 374 0 409 351 240 + 190 258 0 332 320 232 139 494 310 0 232 354 261 113 372 188 208 0 + 464 124 88 171 202 328 188 284 0 566 508 397 347 331 171 467 345 485 + 0 552 80 127 259 234 365 249 372 61 522 0 802 316 336 509 222 470 + 588 584 392 502 386 0 633 432 479 552 586 723 417 621 411 874 354 738 + 0 257 641 541 407 706 522 184 391 372 679 433 915 390 0 187 577 477 + 337 636 452 375 321 507 609 595 845 572 196 0 91 450 357 210 509 325 + 248 141 380 482 468 718 661 228 158 0 412 624 531 384 690 506 210 408 + 398 663 459 892 227 169 351 383 0 400 752 659 512 818 634 338 536 526 + 791 587 1020 524 151 270 371 167 0 472 805 712 565 871 687 391 589 579 + 844 640 1073 413 257 342 443 220 57 0 389 665 572 425 731 547 251 449 + 439 704 500 933 274 146 328 360 53 112 165 0 610 76 183 317 192 442 + 396 430 202 515 141 233 492 723 653 526 700 828 881 741 0 340 730 630 + 490 789 605 394 474 582 762 643 998 444 125 185 311 223 67 139 168 806 + 0 510 152 134 217 248 370 175 330 46 527 72 438 380 359 553 426 385 + 513 566 426 213 569 0 153 447 354 207 470 280 246 113 377 437 465 715 + 659 345 275 98 500 488 560 477 523 428 423 0 511 844 751 604 910 726 + 430 628 618 883 679 1112 407 296 381 482 259 96 39 204 920 178 605 599 + 0 269 283 190 42 332 148 169 63 213 305 301 544 582 382 312 185 365 + 493 546 406 359 465 259 170 585 0 525 157 95 232 42 257 316 355 160 + 330 167 254 521 638 568 441 615 743 796 656 188 721 206 438 835 274 0 + 150 539 446 299 598 414 279 283 469 571 557 807 488 112 96 120 267 255 + 327 244 615 195 515 237 366 274 530 0 80 507 414 267 566 382 305 251 + 437 539 525 775 572 196 88 77 351 339 411 328 583 279 483 205 450 242 + 498 63 0 130 520 427 280 579 395 318 264 450 552 538 788 543 167 59 + 101 322 310 382 299 596 250 496 218 412 255 511 56 25 0 401 791 691 + 551 850 666 474 535 662 823 723 1059 524 238 238 372 303 138 126 248 867 + 101 649 489 165 526 782 256 340 311 0 134 524 431 284 583 399 314 268 + 454 556 542 792 530 154 63 105 309 297 369 286 600 237 500 222 408 259 + 515 34 29 22 298 0 666 942 849 702 1008 824 528 726 716 981 777 1210 + 446 423 605 637 357 280 217 279 1018 336 703 754 194 683 933 521 605 576 + 416 563 0 259 281 188 40 364 180 142 72 211 337 299 549 555 372 302 + 175 338 466 519 379 357 455 257 172 558 27 272 264 232 245 516 249 656 + 0 505 447 336 286 354 110 406 284 424 70 461 566 819 618 548 421 602 + 730 783 643 538 701 466 376 822 244 353 509 478 491 762 494 920 276 0 + 453 358 247 234 265 59 354 232 335 182 372 477 767 566 496 369 550 678 + 731 591 449 649 377 324 770 192 264 458 426 439 710 443 868 224 95 0 + 627 334 251 408 168 239 528 406 300 166 348 331 700 740 670 504 724 852 + 905 765 360 823 346 504 944 366 179 632 600 613 884 617 1042 398 162 177 + 0 339 275 187 94 313 281 45 184 143 438 204 543 458 229 382 250 255 + 350 436 296 351 439 130 248 475 136 271 286 312 325 519 321 573 102 377 + 325 425 0 710 283 254 491 117 375 611 489 319 354 351 202 679 823 753 + 626 807 935 988 848 272 906 365 587 1027 449 159 715 683 696 967 700 1125 + 481 345 337 183 430 0 243 353 260 113 419 235 89 133 283 392 368 621 + 502 209 286 159 285 413 466 326 429 439 226 157 505 94 344 190 216 229 + 500 225 603 67 331 279 453 96 536 0 376 520 427 280 586 402 106 300 + 294 559 355 788 340 146 322 306 112 240 293 153 596 296 281 338 332 261 + 511 220 322 293 376 253 430 234 498 446 620 151 703 181 0 449 594 501 + 354 660 476 180 378 368 633 429 862 256 206 388 420 45 204 257 101 670 + 260 355 537 296 335 585 304 388 359 340 346 394 308 572 520 694 225 777 + 255 82 0 505 781 688 541 847 663 367 565 555 820 616 1049 289 262 444 + 476 196 119 136 118 857 175 542 593 129 522 772 360 444 415 255 402 161 + 495 759 707 881 412 964 442 269 233 0 322 611 518 371 677 493 197 395 + 372 650 446 879 359 69 261 293 94 142 224 84 687 146 372 410 263 352 + 602 177 261 232 247 219 361 325 589 537 711 242 794 272 99 106 200 0 + 185 575 482 335 634 450 231 319 419 607 480 843 462 86 124 156 241 226 + 298 218 651 166 406 273 337 310 566 40 124 95 227 82 495 300 546 494 + 668 276 751 207 212 278 334 151 0 353 638 545 398 704 520 224 422 412 + 677 473 906 282 110 292 324 61 125 178 38 714 181 399 441 217 379 629 + 208 292 263 261 250 315 352 616 564 738 269 821 299 126 82 154 46 182 + 0 324 314 191 106 275 91 225 104 244 248 332 487 638 437 367 240 421 + 549 602 462 390 520 290 238 641 64 227 329 297 310 581 314 739 95 187 + 135 309 196 392 150 317 391 578 408 365 435 0 388 234 127 124 219 113 + 289 168 215 270 254 431 606 501 431 304 485 613 666 526 310 584 232 302 + 705 128 163 393 361 374 645 378 803 159 209 120 219 229 336 214 381 455 + 642 472 429 499 64 0 447 664 571 424 730 546 250 448 438 703 499 932 + 188 204 386 418 36 202 255 88 740 258 425 535 294 405 655 302 386 357 + 338 344 392 378 642 590 764 295 847 325 152 68 231 130 276 96 461 525 + 0 360 606 513 366 672 488 192 390 380 645 441 874 273 117 299 331 46 + 150 203 63 682 206 367 448 242 347 597 215 299 270 286 257 340 320 584 + 532 706 237 789 267 94 58 179 48 189 31 403 467 82 0 605 133 180 + 312 287 418 264 425 112 575 55 439 313 448 648 520 474 602 655 515 193 + 658 89 517 694 354 220 610 577 591 738 595 792 352 514 425 401 219 404 + 421 370 444 631 461 495 488 385 307 514 456 0 656 932 839 692 998 814 + 518 716 706 971 767 1200 444 413 595 627 347 270 200 269 1008 326 693 744 + 177 673 923 511 595 566 406 553 42 646 910 858 1032 563 1115 593 420 384 + 151 351 485 305 729 793 382 330 782 0 573 113 101 280 79 329 359 403 + 163 402 157 235 509 686 616 489 663 791 844 704 131 769 209 486 883 322 + 57 578 546 559 830 563 981 320 425 336 236 314 176 392 559 633 820 650 + 614 677 353 220 603 645 210 971 0 293 384 274 143 319 129 262 60 314 + 286 402 531 675 451 381 201 458 586 639 499 460 534 360 151 678 101 318 + 343 311 324 595 328 776 132 225 173 353 233 436 187 354 428 615 445 379 + 472 83 147 498 440 455 766 375 0 372 283 199 153 229 39 273 151 324 + 196 324 441 686 485 415 288 469 597 650 510 413 568 370 241 689 111 228 + 377 345 358 629 362 787 143 135 83 263 244 346 198 365 439 626 456 413 + 483 54 72 509 451 377 777 300 90 0 330 479 386 239 545 361 65 259 + 253 518 314 747 378 119 276 260 150 278 331 191 555 334 240 297 370 220 + 470 174 276 247 414 207 468 193 457 405 579 110 662 140 46 120 307 126 + 166 164 276 340 190 132 329 458 518 313 324 0 610 297 234 391 107 275 + 511 389 322 248 331 254 693 723 653 526 707 835 888 748 302 806 368 487 + 927 349 149 615 583 596 867 600 1025 381 239 231 77 408 106 436 603 677 + 864 694 651 721 292 236 747 689 384 1015 186 336 246 562 0 598 874 781 + 634 940 756 460 658 648 913 709 1142 370 355 537 569 289 212 197 211 950 + 268 635 686 157 615 865 453 537 508 348 495 84 588 852 800 974 505 1057 + 535 362 326 93 293 427 247 671 735 324 272 724 85 913 708 719 400 957 + 0 214 604 504 364 663 479 402 348 534 636 622 872 599 223 49 185 378 + 319 391 355 680 234 580 302 430 339 595 123 115 86 287 90 632 329 575 + 523 697 409 780 313 349 415 471 288 151 319 394 458 413 326 675 622 643 + 408 442 303 680 564 0 154 401 308 161 460 276 199 78 331 433 419 669 + 612 298 228 63 453 441 513 430 477 381 377 47 552 136 392 190 158 171 + 442 175 707 126 372 320 494 201 577 110 291 490 546 363 226 396 191 255 + 488 401 471 697 440 138 239 250 477 639 255 0 70 464 371 224 523 339 + 262 208 394 496 482 732 567 191 121 27 346 334 406 323 540 274 440 162 + 445 199 455 83 47 64 335 68 600 189 435 383 557 269 640 173 295 383 + 439 256 119 287 254 378 381 294 534 590 503 268 302 249 540 532 148 115 + 0 606 349 266 387 183 216 507 385 354 147 363 357 715 719 649 522 703 + 831 884 744 375 802 400 483 923 345 194 611 579 592 863 596 1021 377 139 + 154 22 447 209 432 599 673 860 690 647 717 288 232 743 685 416 1011 262 + 332 242 558 103 953 676 473 536 0 631 228 175 412 38 296 532 410 240 + 306 272 210 640 744 674 547 728 856 909 769 233 827 286 508 948 370 80 + 636 604 617 888 621 1046 402 293 261 138 351 79 457 624 698 885 715 672 + 742 313 257 768 710 325 1036 117 357 267 583 69 978 701 498 561 153 0 + 642 129 176 349 121 369 428 472 232 428 226 187 578 755 685 558 732 860 + 913 773 98 838 278 555 952 391 132 647 615 628 899 632 1050 389 465 376 + 260 383 161 461 628 702 889 719 683 746 422 330 772 714 279 1040 75 430 + 340 587 191 982 712 509 572 275 122 0 503 779 686 539 845 661 365 563 + 553 818 614 1047 245 260 442 474 141 151 168 116 855 207 540 591 162 520 + 770 358 442 413 287 400 201 493 757 705 879 410 962 440 267 231 44 198 + 332 152 576 640 174 177 629 199 818 613 624 305 862 125 469 544 437 858 + 883 887 0 372 762 662 522 821 637 456 506 644 794 705 1030 506 209 209 + 343 285 120 108 230 838 72 631 460 147 497 753 227 311 282 29 269 398 + 487 733 681 855 501 938 471 354 322 237 218 198 243 552 616 320 268 720 + 388 801 566 600 396 838 330 258 413 306 834 859 870 269 0 641 348 265 + 422 152 262 542 420 314 189 362 313 714 754 684 557 738 866 919 779 344 + 837 360 518 958 380 193 646 614 627 898 631 1056 412 185 200 23 439 165 + 467 634 708 895 725 682 752 323 233 778 720 415 1046 231 367 277 593 59 + 988 711 508 571 45 122 244 893 869 0 561 837 744 597 903 719 423 621 + 611 876 672 1105 311 318 500 532 252 175 192 174 913 231 598 649 186 578 + 828 416 500 471 311 458 154 551 815 763 937 468 1020 498 325 289 65 256 + 390 210 634 698 287 235 687 152 876 671 682 363 920 77 527 602 495 916 + 941 945 66 293 951 0 478 754 661 514 820 636 340 538 528 793 589 1022 + 261 235 417 449 116 132 149 91 830 188 515 566 159 495 745 333 417 388 + 268 375 226 468 732 680 854 385 937 415 242 206 55 173 307 127 551 615 + 149 152 604 224 793 588 599 280 837 150 444 519 412 833 858 862 25 250 + 868 91 0 247 316 223 76 360 176 170 39 246 333 334 572 583 360 290 + 163 366 494 547 407 392 443 292 133 586 37 307 252 220 233 504 237 684 + 34 272 220 394 146 477 95 262 336 523 353 288 380 92 156 406 348 387 + 674 355 83 139 221 377 616 317 92 177 373 398 424 521 475 408 579 496 + 0 317 336 212 95 289 105 218 79 266 262 354 501 631 430 360 233 414 + 542 595 455 412 513 312 213 634 53 248 397 290 378 574 382 732 88 201 + 149 323 189 406 143 310 384 571 401 358 428 21 85 454 396 407 722 375 + 62 68 269 306 664 387 184 247 302 327 444 569 545 337 627 544 70 0 + 272 382 289 142 448 264 84 162 234 421 295 650 497 180 315 188 280 408 + 461 321 458 464 221 186 500 123 373 193 245 258 544 228 598 96 360 308 + 482 91 565 29 142 250 437 267 159 294 179 243 320 262 310 588 421 216 + 227 96 465 530 342 139 209 461 486 490 435 526 496 493 410 124 172 0 + 575 276 199 356 86 240 476 354 287 250 296 266 672 688 618 491 672 800 + 853 713 289 771 333 452 892 314 127 580 548 561 832 565 990 346 237 205 + 82 373 141 401 568 642 829 659 616 686 257 201 712 654 349 980 165 301 + 211 527 35 922 645 442 505 97 56 178 827 803 66 885 802 342 271 430 + 0 219 479 386 239 545 361 167 259 317 518 378 747 474 83 172 149 253 + 235 339 230 555 228 304 263 378 220 470 79 139 134 289 112 507 193 457 + 405 579 174 662 126 154 290 346 136 621 194 276 340 288 201 393 497 518 + 313 324 108 562 439 199 216 153 558 583 587 344 260 593 402 319 221 269 + 97 527 0 293 683 583 443 742 558 238 427 426 715 487 951 352 50 232 + 264 131 101 174 108 759 105 413 381 213 418 674 148 232 203 206 190 385 + 408 654 602 776 283 859 248 140 168 224 41 122 72 473 537 166 89 502 + 375 722 487 521 167 759 317 259 334 227 755 780 791 222 177 790 280 197 + 396 466 219 724 134 0 54 529 429 289 588 404 327 273 459 561 547 797 + 595 219 92 82 374 362 434 351 605 302 505 227 473 264 520 119 31 43 + 363 58 628 254 500 448 622 334 705 238 327 411 467 284 147 315 319 383 + 409 322 600 618 568 333 367 281 605 560 84 180 53 601 626 637 465 334 + 636 523 440 242 312 267 570 170 255 0 648 188 182 355 68 316 434 430 + 238 362 232 154 584 761 691 564 738 866 919 779 177 844 284 561 958 390 + 100 653 621 634 905 638 1056 395 412 323 194 389 95 467 634 708 895 725 + 689 752 333 277 778 720 285 1046 81 377 287 593 125 988 718 515 578 209 + 56 66 893 876 178 951 868 418 347 496 112 593 797 643 0 211 601 501 + 361 660 476 231 345 479 633 480 869 466 74 81 182 243 189 261 220 677 + 129 406 299 300 336 592 105 150 121 190 108 497 326 572 520 694 276 777 + 310 204 280 336 143 37 184 391 455 278 191 495 487 640 405 439 166 677 + 429 160 252 145 673 698 709 334 161 708 392 309 314 384 196 642 99 125 + 173 715 0 568 844 751 604 910 726 430 628 618 883 679 1112 348 325 507 + 539 259 182 150 181 920 238 605 656 127 585 835 423 507 478 318 465 98 + 558 822 770 944 475 1027 505 332 296 63 263 397 217 641 705 294 242 694 + 96 883 678 689 370 927 30 534 609 502 923 948 952 103 300 958 56 120 + 586 634 500 892 409 287 530 958 399 0 497 150 67 204 70 257 288 327 + 155 335 164 282 516 610 540 413 587 715 768 628 216 693 201 410 807 246 + 28 502 470 483 754 487 905 244 353 264 184 243 187 316 483 557 744 574 + 538 601 199 135 627 569 217 895 85 292 228 442 167 837 567 364 427 199 + 108 160 742 725 198 800 717 279 220 345 132 442 646 492 128 564 807 0 + 290 680 580 440 739 555 317 424 505 712 566 948 506 139 98 261 285 155 + 227 229 756 88 492 378 266 415 671 144 176 164 140 136 424 405 651 599 + 773 362 856 389 290 322 263 195 116 226 470 534 320 243 581 414 719 484 + 518 252 756 356 147 331 224 752 777 788 295 111 787 319 276 393 463 275 + 721 178 154 190 794 79 326 643 0 475 65 42 182 137 282 261 295 65 + 439 85 321 437 588 518 391 565 693 746 606 141 671 111 388 785 224 95 + 480 448 461 732 465 883 222 378 289 272 216 254 294 461 535 722 552 516 + 579 255 169 605 547 138 873 92 325 241 420 255 815 545 342 405 287 175 + 161 720 703 286 778 695 257 277 323 220 420 624 470 167 542 785 88 621 + 0 654 341 278 435 151 319 555 433 366 266 375 298 755 767 697 570 751 + 879 932 792 346 850 412 531 971 393 193 659 627 640 911 644 1069 425 262 + 254 100 452 103 480 647 721 908 738 695 765 336 280 791 733 428 1059 230 + 380 290 606 44 1001 724 521 584 122 113 235 906 882 77 964 881 421 350 + 509 79 606 803 649 169 721 971 211 800 299 0 445 387 276 226 294 50 + 346 224 364 125 410 506 759 558 488 361 542 670 723 583 478 641 406 316 + 762 184 293 450 418 431 702 435 860 216 64 57 193 317 411 271 438 512 + 699 529 486 556 127 149 582 524 454 850 365 165 75 397 311 792 515 312 + 375 170 332 405 697 673 216 755 672 212 141 300 276 397 594 440 352 512 + 762 293 591 318 355 0 375 765 665 525 824 640 384 509 572 797 633 1033 + 434 160 245 346 213 48 59 158 841 42 559 463 98 500 756 230 314 285 + 90 272 326 490 736 684 858 429 941 474 286 250 165 169 201 171 555 619 + 248 196 648 316 804 569 603 324 841 258 294 416 309 837 862 873 197 72 + 872 221 178 478 548 454 832 263 128 337 879 164 228 728 130 706 885 676 + 0 268 658 558 418 717 533 231 402 419 690 480 926 420 53 138 239 199 + 132 204 202 734 72 406 356 243 393 649 123 207 178 185 165 401 383 629 + 577 751 276 834 367 204 236 240 109 86 140 448 512 234 157 495 391 697 + 462 496 166 734 333 187 309 202 730 755 766 272 156 765 296 253 371 441 + 227 699 130 68 230 772 57 303 627 86 599 778 569 107 0 261 519 426 + 279 585 401 141 299 329 558 390 787 422 43 200 232 211 183 294 178 595 + 162 316 342 333 260 510 98 200 171 275 131 455 233 497 445 619 186 702 + 166 114 248 294 67 90 142 316 380 246 159 405 445 558 353 364 76 602 + 387 227 295 195 598 623 627 292 246 633 350 267 261 309 137 567 69 82 + 223 633 90 357 482 176 460 646 437 197 90 0 710 184 271 417 239 487 + 496 530 300 546 249 168 616 823 753 626 800 928 981 841 108 906 321 623 +1020 459 241 715 683 696 967 700 1118 457 583 494 378 451 279 529 696 770 + 957 787 751 814 490 448 840 782 310 1108 184 560 458 655 309 1050 780 577 + 640 393 240 118 955 938 362 1013 930 492 512 558 296 655 859 705 179 777 +1020 269 856 229 353 523 941 834 695 0 396 291 180 177 198 53 297 175 + 268 218 305 410 710 509 439 312 493 621 674 534 382 592 310 273 713 135 + 197 401 369 382 653 386 811 167 157 67 215 268 315 222 389 463 650 480 + 437 507 78 53 533 475 358 801 269 122 32 348 215 743 466 263 326 206 + 236 309 648 624 238 706 623 163 92 251 180 348 545 391 256 463 713 197 + 542 222 259 97 627 520 388 427 0 295 424 278 183 308 118 302 100 354 + 275 442 520 715 491 421 241 498 626 679 539 500 574 400 191 718 141 307 + 383 351 364 635 368 816 172 214 162 342 273 425 227 394 468 655 485 419 + 512 114 151 538 480 495 806 364 40 79 353 325 748 448 178 308 321 346 + 419 653 606 356 711 628 123 93 256 290 353 527 373 366 445 718 281 524 + 365 369 154 609 502 393 537 111 0 651 927 834 687 993 809 513 711 701 + 966 762 1195 407 408 590 622 342 265 282 264 1003 321 688 739 276 668 918 + 506 590 561 401 548 174 641 905 853 1027 558 1110 588 415 379 155 346 480 + 300 724 788 377 325 777 173 966 761 772 453 1010 90 617 692 585 1006 1031 +1035 162 383 1041 96 187 669 717 583 975 492 370 613 1041 482 112 890 409 + 868 1054 845 311 386 440 1103 796 801 0 175 565 472 325 624 440 311 309 + 495 597 583 833 504 128 76 146 283 271 343 260 641 211 541 263 382 300 + 556 32 76 47 272 30 537 290 536 484 658 318 741 222 256 320 376 193 + 56 224 355 419 318 231 636 527 604 369 403 210 641 469 103 216 109 637 + 662 673 374 243 672 432 349 278 423 225 606 104 164 99 679 57 439 528 + 112 506 685 476 246 114 134 741 427 409 522 0 585 67 146 292 135 385 + 371 405 175 458 147 249 499 698 628 501 675 803 856 716 57 781 221 498 + 895 334 131 590 558 571 842 575 993 332 481 392 303 326 215 404 571 645 + 832 662 626 689 365 261 715 657 200 983 74 435 356 530 245 925 655 452 + 515 318 176 62 830 813 287 888 805 367 387 433 232 530 734 580 120 652 + 895 159 731 104 289 421 816 709 570 121 325 438 978 616 0 250 640 540 + 400 699 515 277 384 465 672 526 908 466 99 89 221 245 178 250 220 716 + 96 452 338 289 375 631 105 189 160 151 147 447 365 611 559 733 322 816 + 349 250 282 286 155 76 186 430 494 280 203 541 437 679 444 478 212 716 + 379 138 291 184 712 737 748 318 122 747 342 299 353 423 235 681 138 114 + 212 754 39 349 603 40 581 760 551 138 46 136 816 502 484 432 96 691 + 0 717 221 251 424 137 385 503 499 307 417 301 95 653 830 760 633 807 + 935 988 848 190 913 353 631 1027 459 169 722 690 703 974 707 1125 464 481 + 392 246 458 117 536 703 777 964 794 758 821 402 346 847 789 354 1115 150 + 446 356 662 169 1057 787 584 647 272 125 92 962 945 228 1020 937 487 416 + 565 181 662 866 712 69 784 1027 197 863 236 213 421 948 841 702 162 325 + 435 1110 748 154 823 0 246 454 361 213 373 183 332 130 384 340 472 585 + 745 472 402 237 528 656 709 569 530 555 430 167 748 171 372 364 332 345 + 616 349 846 202 279 227 407 303 490 257 424 498 685 515 400 542 157 221 + 568 510 525 836 429 70 144 383 390 778 429 174 289 386 411 484 683 587 + 421 741 658 153 132 286 355 383 508 354 431 426 748 346 505 395 434 219 + 590 483 423 630 176 65 831 390 505 465 500 0 788 302 322 495 208 456 + 574 570 378 488 372 23 724 901 831 704 818 1006 1059 919 203 984 424 701 +1098 530 240 793 761 774 1045 778 1196 535 552 463 317 529 188 607 774 848 +1035 865 829 892 473 417 918 860 425 1186 221 517 427 733 240 1128 858 655 + 718 343 196 173 1033 1016 299 1091 1008 558 487 636 252 733 937 783 154 855 +1098 268 934 307 284 492 1019 912 773 138 396 506 1181 819 235 894 81 571 + 0 426 596 503 356 662 478 182 380 370 635 431 864 253 183 365 397 27 + 181 234 82 672 237 357 514 273 337 587 281 365 336 317 323 371 310 574 + 522 696 227 779 257 84 19 210 83 255 67 393 457 65 35 446 361 635 + 430 441 122 679 303 392 467 360 675 700 704 208 299 710 266 183 338 386 + 252 644 267 145 388 710 257 273 559 299 537 723 514 227 213 225 772 465 + 470 355 297 647 259 779 500 850 0 596 575 330 377 237 179 497 375 552 + 100 589 407 941 709 639 512 693 821 874 734 421 792 460 467 913 335 236 + 601 569 582 853 586 1011 367 91 117 72 468 259 422 589 663 850 680 637 + 707 278 221 733 675 642 1001 308 316 212 548 153 943 666 463 526 47 212 + 334 848 824 95 906 823 363 292 451 156 548 745 591 268 663 913 241 742 + 506 172 133 827 720 588 452 185 305 996 627 364 702 322 370 393 665 0 + 634 910 817 670 976 792 496 694 684 949 745 1178 414 391 573 605 325 248 + 185 247 986 304 671 722 162 651 901 489 573 544 384 531 32 624 888 836 +1010 541 1093 571 398 362 129 329 463 283 707 771 360 308 760 33 949 744 + 755 436 983 52 600 675 568 989 1014 1018 169 366 1024 122 194 652 700 566 + 958 475 353 596 1024 465 66 873 392 851 1037 828 294 369 423 1086 779 784 + 142 505 961 415 1003 814 1164 339 979 0 507 209 111 201 139 172 408 286 + 213 329 223 351 575 620 550 423 604 732 785 645 275 703 259 384 824 202 + 87 512 480 493 764 497 922 278 268 141 173 252 256 333 500 574 761 591 + 548 618 138 74 644 586 276 912 144 233 143 459 156 854 577 374 437 208 + 177 219 759 735 187 817 734 274 159 362 121 459 656 502 197 574 824 59 + 653 147 200 208 738 631 499 328 112 222 907 538 218 613 266 287 337 576 + 208 890 0 463 186 79 155 157 161 251 216 167 318 206 369 558 576 506 + 379 553 681 734 594 262 659 213 376 773 176 115 468 436 449 720 453 871 + 210 257 168 219 206 274 282 449 523 710 540 504 567 112 48 593 535 259 + 861 172 195 120 408 202 803 533 330 393 254 195 247 708 691 233 766 683 + 204 133 311 167 408 612 458 215 530 773 87 609 121 246 197 694 587 448 + 350 101 199 856 494 225 569 284 269 355 525 251 839 46 0 408 169 105 + 116 242 298 131 229 57 455 118 437 468 315 451 325 341 469 522 382 245 + 525 72 322 561 158 200 413 382 394 605 398 659 155 394 305 344 86 359 + 228 237 311 498 328 362 355 188 160 381 323 169 649 208 259 269 196 327 + 591 478 276 339 391 280 277 496 587 358 554 471 191 211 177 292 260 369 + 403 283 362 561 172 448 110 371 334 515 362 272 345 238 299 644 439 220 + 408 352 329 423 313 388 627 183 137 0 529 389 278 310 236 127 430 308 + 366 125 403 447 755 642 572 445 626 754 807 667 420 725 408 400 846 268 + 235 534 502 515 786 519 944 300 80 68 112 401 299 355 522 596 783 613 + 570 640 211 169 666 608 456 934 307 249 225 481 193 876 599 396 459 89 + 231 353 781 757 135 839 756 296 225 384 175 481 678 524 287 596 846 240 + 675 320 212 84 760 653 521 471 133 238 929 560 363 635 362 303 433 598 + 52 912 156 199 336 0 192 412 319 172 477 293 160 154 342 450 430 680 + 573 228 235 108 383 371 443 360 488 311 388 167 482 153 403 119 165 178 + 372 154 637 126 389 337 511 162 594 71 207 420 476 232 136 324 209 273 + 418 331 482 627 451 246 256 161 494 569 262 120 110 490 515 590 474 343 + 525 532 449 127 202 74 459 96 264 187 526 182 539 375 261 353 538 329 + 346 239 165 588 280 286 622 151 463 221 595 294 666 397 480 605 391 341 + 287 413 0 529 286 231 310 177 165 430 308 319 182 328 379 680 642 572 + 445 626 754 807 667 344 725 365 406 846 268 159 534 502 515 786 519 944 + 300 137 106 75 370 231 355 522 596 783 613 570 640 211 155 666 608 381 + 934 231 255 165 481 125 876 599 396 459 77 155 247 781 757 89 839 756 + 296 225 384 99 481 678 524 211 596 846 164 675 252 148 168 760 653 521 + 365 137 244 929 560 287 635 294 309 365 598 80 912 131 177 314 87 413 + 0 434 710 617 470 776 592 296 494 484 749 545 978 346 191 373 405 125 + 82 142 47 786 145 471 522 145 451 701 289 373 344 225 331 238 424 688 + 636 810 341 893 371 198 162 77 129 263 83 507 571 160 108 560 228 749 + 544 555 236 793 170 400 475 368 789 814 818 76 207 824 133 51 452 500 + 366 758 275 153 396 824 265 140 673 233 651 837 628 135 210 223 886 579 + 584 223 305 761 241 893 614 964 139 779 206 690 639 427 712 405 712 0 + 535 811 718 571 877 693 397 595 585 850 646 1079 336 292 474 506 226 94 + 76 129 887 154 572 623 75 552 802 390 474 445 202 432 175 525 789 737 + 911 442 994 472 299 263 58 230 364 184 608 672 261 209 661 160 850 645 + 656 337 894 140 501 576 469 890 915 919 91 184 925 113 88 553 601 467 + 859 376 176 497 925 276 110 774 242 752 938 729 123 219 324 987 680 685 + 205 406 862 265 994 715 1065 240 880 143 791 740 528 813 506 813 82 0 + 630 108 191 337 165 424 416 450 220 483 188 190 540 743 673 546 720 848 + 901 761 43 826 266 543 940 379 161 635 603 616 887 620 1038 877 520 431 + 315 371 216 449 616 690 877 707 671 734 410 306 760 702 241 1028 104 480 + 395 575 246 970 700 497 560 348 177 55 875 858 299 933 850 412 432 478 + 233 575 779 625 121 697 940 189 776 149 290 460 861 754 615 80 364 474 +1023 661 41 736 147 550 160 692 394 1006 248 270 265 393 508 317 806 907 + 0 446 252 145 227 162 111 347 225 233 268 272 374 624 559 489 362 543 + 671 724 584 346 642 279 323 763 185 161 451 419 432 703 436 861 217 207 + 141 162 272 279 272 439 513 700 530 487 557 128 50 583 525 325 851 233 + 172 82 398 179 793 516 313 376 175 200 273 698 674 176 756 673 213 142 + 301 144 398 595 441 220 513 763 108 592 187 223 147 677 570 438 391 51 + 161 846 477 289 552 289 226 360 515 167 829 49 66 203 115 330 98 629 + 730 319 0 166 518 425 277 437 247 336 114 448 404 536 649 749 435 365 + 150 590 578 650 567 594 518 494 90 689 235 436 402 295 383 579 387 844 + 189 343 291 471 295 554 247 428 627 683 500 363 531 221 285 625 538 589 + 834 493 134 208 372 454 776 392 137 252 450 475 548 681 550 485 739 656 + 150 196 276 419 353 471 317 495 389 746 410 468 459 498 283 553 446 432 + 694 240 129 829 353 569 428 564 80 635 604 434 812 351 333 393 367 257 + 373 612 713 614 290 0 471 313 202 252 166 99 372 250 290 210 358 378 + 679 584 514 387 568 696 749 609 350 667 332 348 788 210 165 476 444 457 + 728 461 886 242 156 61 135 311 276 297 464 538 725 555 512 582 153 93 + 608 550 380 876 237 197 107 423 170 818 541 338 401 147 191 277 723 699 + 149 781 698 238 167 326 135 423 620 466 224 538 788 139 617 244 214 118 + 702 595 462 395 67 186 871 502 293 577 293 251 364 540 128 854 80 123 + 260 76 355 70 654 755 323 39 315 0 442 718 625 478 784 600 304 502 + 492 757 553 986 300 199 381 413 81 137 184 53 794 209 479 530 194 459 + 709 297 381 352 285 339 261 432 696 644 818 349 901 379 206 126 90 137 + 271 91 515 579 111 116 568 259 757 552 563 244 801 185 408 483 376 797 + 822 826 60 267 832 126 35 460 508 374 766 283 161 404 832 273 155 681 + 293 659 845 636 195 270 231 894 587 592 222 313 769 301 901 622 972 108 + 787 229 698 647 435 720 413 720 55 123 814 637 620 662 0 523 230 147 + 304 81 188 424 302 235 255 174 293 596 636 566 439 620 748 801 661 265 + 719 281 400 840 262 75 528 496 509 780 513 938 294 284 195 104 321 193 + 349 516 590 777 607 564 634 205 149 660 602 297 928 147 249 159 475 87 + 870 593 390 453 119 108 192 775 751 118 833 750 290 219 378 52 475 672 + 518 139 590 840 80 669 168 131 224 754 647 515 310 128 238 923 554 208 + 629 208 303 279 592 161 906 69 115 240 160 407 84 706 807 238 92 367 + 87 714 0 566 45 139 273 228 383 352 386 121 540 60 314 411 679 609 + 482 656 784 837 697 81 762 132 479 876 315 189 571 539 552 823 556 974 + 313 479 390 366 307 308 385 552 626 813 643 607 670 346 266 696 638 112 + 964 158 416 342 511 335 906 636 433 496 381 266 155 811 794 380 869 786 + 348 368 414 314 511 715 561 213 633 876 182 712 97 379 419 797 690 551 + 189 323 456 959 597 93 672 247 486 284 628 607 942 244 218 178 421 444 + 346 742 843 124 284 550 345 750 262 0 235 313 220 73 371 187 168 43 + 243 344 331 583 581 348 278 151 364 492 545 405 389 431 289 137 584 48 + 304 240 208 221 492 225 682 32 283 231 405 144 488 93 260 334 521 351 + 276 378 103 167 404 346 384 672 352 95 150 219 388 614 305 94 165 384 + 409 421 519 463 419 577 494 12 92 122 353 219 384 230 429 302 584 276 + 381 254 432 223 466 359 259 489 174 135 667 266 364 341 498 165 569 336 + 374 650 285 215 188 307 115 307 450 551 409 224 154 249 458 301 345 0 + 432 822 722 582 881 697 441 566 629 854 690 1090 491 217 302 403 270 105 + 69 215 898 99 616 520 108 557 813 287 371 342 77 329 383 547 793 741 + 915 486 998 531 343 307 222 226 258 228 612 676 305 253 705 373 861 626 + 660 381 898 315 351 473 366 894 919 930 254 66 929 278 235 535 605 511 + 863 320 185 394 936 221 285 785 177 763 942 733 57 164 254 998 684 666 + 368 303 873 195 1005 647 1076 284 884 351 795 751 572 817 403 817 192 145 + 819 734 610 759 252 811 854 523 0 435 653 560 413 719 535 239 437 427 + 692 488 921 220 192 374 406 29 190 243 64 729 246 414 523 282 394 644 + 290 374 345 326 332 380 367 631 579 753 284 836 314 141 78 219 123 264 + 84 450 514 34 75 503 370 692 487 498 179 736 312 401 476 369 732 757 + 761 137 308 767 275 112 395 443 309 701 276 154 397 767 266 282 616 308 + 594 780 571 236 222 234 829 522 527 365 306 704 268 836 557 907 56 722 + 348 633 582 370 655 406 655 148 249 749 572 613 597 77 649 685 393 293 + 0 369 167 79 77 205 259 153 190 97 416 185 435 537 482 412 286 459 + 587 640 500 243 555 111 283 679 119 163 375 343 356 626 360 777 116 355 + 266 305 108 322 189 356 429 616 446 411 473 149 121 499 441 238 767 206 + 220 230 315 288 709 439 237 300 352 243 275 614 587 319 672 589 152 172 + 218 253 315 518 364 281 436 679 135 515 108 332 295 600 493 355 343 199 + 260 762 401 218 475 350 290 421 431 349 745 144 98 39 297 248 275 545 + 646 263 164 354 221 553 201 199 149 657 488 0 121 511 418 271 570 386 + 309 255 441 543 529 779 518 142 99 84 297 285 357 274 587 225 487 209 + 396 246 502 35 29 42 286 36 551 236 482 430 604 316 687 220 268 334 + 390 207 70 238 301 365 332 245 581 541 550 315 349 222 587 483 126 162 + 55 583 608 619 388 257 618 446 363 224 294 249 552 104 178 60 625 96 + 453 474 175 452 631 422 260 153 146 687 373 355 536 47 562 135 694 336 + 765 311 573 519 484 440 386 506 169 506 319 520 607 423 299 448 327 500 + 543 212 317 320 347 0 +DISPLAY_DATA_SECTION + 1 8.0 124.0 + 2 125.0 80.0 + 3 97.0 74.0 + 4 69.0 96.0 + 5 106.0 46.0 + 6 49.0 57.0 + 7 80.0 125.0 + 8 42.0 93.0 + 9 104.0 94.0 + 10 35.0 17.0 + 11 118.0 96.0 + 12 151.0 22.0 + 13 154.0 182.0 + 14 57.0 165.0 + 15 18.0 159.0 + 16 27.0 123.0 + 17 96.0 170.0 + 18 63.0 198.0 + 19 59.0 211.0 + 20 88.0 182.0 + 21 142.0 72.0 + 22 48.0 190.0 + 23 106.0 106.0 + 24 28.0 102.0 + 25 63.0 224.0 + 26 58.0 93.0 + 27 103.0 56.0 + 28 38.0 149.0 + 29 23.0 138.0 + 30 22.0 146.0 + 31 32.0 208.0 + 32 27.0 144.0 + 33 75.0 258.0 + 34 59.0 101.0 + 35 41.0 32.0 + 36 53.0 46.0 + 37 76.0 19.0 + 38 79.0 115.0 + 39 109.0 13.0 + 40 59.0 118.0 + 41 84.0 147.0 + 42 95.0 160.0 + 43 87.0 213.0 + 44 73.0 166.0 + 45 43.0 153.0 + 46 81.0 175.0 + 47 59.0 77.0 + 48 70.0 68.0 + 49 106.0 169.0 + 50 86.0 168.0 + 51 127.0 109.0 + 52 68.0 243.0 + 53 116.0 57.0 + 54 39.0 78.0 + 55 54.0 65.0 + 56 77.0 141.0 + 57 95.0 24.0 + 58 89.0 238.0 + 59 9.0 158.0 + 60 39.0 109.0 + 61 25.0 129.0 + 62 69.0 20.0 + 63 104.0 34.0 + 64 132.0 51.0 + 65 98.0 207.0 + 66 37.0 203.0 + 67 80.0 16.0 + 68 103.0 224.0 + 69 94.0 202.0 + 70 49.0 96.0 + 71 55.0 80.0 + 72 62.0 123.0 + 73 91.0 31.0 + 74 51.0 142.0 + 75 64.0 172.0 + 76 14.0 138.0 + 77 120.0 37.0 + 78 38.0 160.0 + 79 86.0 230.0 + 80 96.0 59.0 + 81 33.0 177.0 + 82 108.0 78.0 + 83 93.0 12.0 + 84 43.0 47.0 + 85 52.0 200.0 + 86 48.0 171.0 + 87 62.0 153.0 + 88 159.0 53.0 + 89 60.0 60.0 + 90 36.0 73.0 + 91 111.0 243.0 + 92 31.0 150.0 + 93 130.0 67.0 + 94 36.0 172.0 + 95 132.0 28.0 + 96 29.0 73.0 + 97 150.0 28.0 + 98 89.0 166.0 + 99 58.0 21.0 + 100 78.0 244.0 + 101 82.0 58.0 + 102 81.0 68.0 + 103 92.0 98.0 + 104 56.0 33.0 + 105 47.0 126.0 + 106 70.0 34.0 + 107 78.0 195.0 + 108 77.0 215.0 + 109 140.0 62.0 + 110 70.0 57.0 + 111 16.0 89.0 + 112 66.0 50.0 + 113 98.0 194.0 + 114 87.0 45.0 + 115 132.0 87.0 + 116 52.0 99.0 + 117 50.0 212.0 + 118 103.0 176.0 + 119 84.0 91.0 + 120 31.0 140.0 +EOF diff --git a/glpk-5.0/examples/tsp/main.c b/glpk-5.0/examples/tsp/main.c new file mode 100644 index 0000000..0685742 --- /dev/null +++ b/glpk-5.0/examples/tsp/main.c @@ -0,0 +1,535 @@ +/* main.c */ + +/* Written by Andrew Makhorin <mao@gnu.org>, October 2015. */ + +/*********************************************************************** +* This program is a stand-alone solver intended for solving Symmetric +* Traveling Salesman Problem (TSP) with the branch-and-bound method. +* +* Please note that this program is only an illustrative example. It is +* *not* a state-of-the-art code, so only small TSP instances (perhaps, +* having up to 150-200 nodes) can be solved with this code. +* +* To run this program use the following command: +* +* tspsol tsp-file +* +* where tsp-file specifies an input text file containing TSP data in +* TSPLIB 95 format. +* +* Detailed description of the input format recognized by this program +* is given in the report: Gerhard Reinelt, "TSPLIB 95". This report as +* well as TSPLIB, a library of sample TSP instances (and other related +* problems), are freely available for research purposes at the webpage +* <http://www.iwr.uni-heidelberg.de/groups/comopt/software/TSPLIB95/>. +* +* Symmetric Traveling Salesman Problem +* ------------------------------------ +* Let a complete undirected graph be given: +* +* K = (V, E), (1) +* +* where V = {1, ..., n} is a set of nodes, E = V cross V is a set of +* edges. Let also each edge e = (i,j) be assigned a positive number +* c[i,j], which is the length of e. The Symmetric Traveling Salesman +* Problem (TSP) is to find a tour in K of minimal length. +* +* Integer programming formulation of TSP +* -------------------------------------- +* For a set of nodes W within V introduce the following notation: +* +* d(W) = {(i,j):i in W and j not in W or i not in W and j in W}, (2) +* +* i.e. d(W) is the set of edges which have exactly one endnode in W. +* If W = {v}, i.e. W consists of the only node, we write simply d(v). +* +* The integer programming formulation of TSP is the following: +* +* minimize sum c[i,j] * x[i,j] (3) +* i,j +* +* subject to sum x[i,j] = 2 for all v in V (4) +* (i,j) in d(v) +* +* sum x[i,j] >= 2 for all W within V, (5) +* (i,j) in d(W) W != empty, W != V +* +* x[i,j] in {0, 1} for all i, j (6) +* +* The binary variables x[i,j] have conventional meaning: if x[i,j] = 1, +* the edge (i,j) is included in the tour, otherwise, if x[i,j] = 0, the +* edge is not included in the tour. +* +* The constraints (4) are called degree constraints. They require that +* for each node v in V there must be exactly two edges included in the +* tour which are incident to v. +* +* The constraints (5) are called subtour elimination constraints. They +* are used to forbid subtours. Note that the number of the subtour +* elimination constraints grows exponentially on the size of the TSP +* instance, so these constraints are not included explicitly in the +* IP, but generated dynamically during the B&B search. +***********************************************************************/ + +#include <errno.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> +#include <string.h> +#include <time.h> + +#include <glpk.h> +#include "maxflow.h" +#include "mincut.h" +#include "misc.h" +#include "tsplib.h" + +int n; +/* number of nodes in the problem, n >= 2 */ + +int *c; /* int c[1+n*(n-1)/2]; */ +/* upper triangle (without diagonal entries) of the (symmetric) matrix + * C = (c[i,j]) in row-wise format, where c[i,j] specifies a length of + * edge e = (i,j), 1 <= i < j <= n */ + +int *tour; /* int tour[1+n]; */ +/* solution to TSP, which is a tour specified by the list of node + * numbers tour[1] -> ... -> tour[nn] -> tour[1] in the order the nodes + * are visited; note that any tour is a permutation of node numbers */ + +glp_prob *P; +/* integer programming problem object */ + +/*********************************************************************** +* loc - determine reduced index of element of symmetric matrix +* +* Given indices i and j of an element of a symmetric nxn-matrix, +* 1 <= i, j <= n, i != j, this routine returns the index of that +* element in an array, which is the upper triangle (without diagonal +* entries) of the matrix in row-wise format. */ + +int loc(int i, int j) +{ xassert(1 <= i && i <= n); + xassert(1 <= j && j <= n); + xassert(i != j); + if (i < j) + return ((n - 1) + (n - i + 1)) * (i - 1) / 2 + (j - i); + else + return loc(j, i); +} + +/*********************************************************************** +* read_data - read TSP data +* +* This routine reads TSP data from a specified text file in TSPLIB 95 +* format. */ + +void read_data(const char *fname) +{ TSP *tsp; + int i, j; + tsp = tsp_read_data(fname); + if (tsp == NULL) + { xprintf("TSP data file processing error\n"); + exit(EXIT_FAILURE); + } + if (tsp->type != TSP_TSP) + { xprintf("Invalid TSP data type\n"); + exit(EXIT_FAILURE); + } + n = tsp->dimension; + xassert(n >= 2); + if (n > 32768) + { xprintf("TSP instance too large\n"); + exit(EXIT_FAILURE); + } + c = xalloc(1+loc(n-1, n), sizeof(int)); + for (i = 1; i <= n; i++) + { for (j = i+1; j <= n; j++) + c[loc(i, j)] = tsp_distance(tsp, i, j); + } + tsp_free_data(tsp); + return; +} + +/*********************************************************************** +* build_prob - build initial integer programming problem +* +* This routine builds the initial integer programming problem, which +* includes all variables (6), objective (3) and all degree constraints +* (4). Subtour elimination constraints (5) are considered "lazy" and +* not included in the initial problem. */ + +void build_prob(void) +{ int i, j, k, *ind; + double *val; + char name[50]; + /* create problem object */ + P = glp_create_prob(); + /* add all binary variables (6) */ + for (i = 1; i <= n; i++) + { for (j = i+1; j <= n; j++) + { k = glp_add_cols(P, 1); + xassert(k == loc(i,j)); + sprintf(name, "x[%d,%d]", i, j); + glp_set_col_name(P, k, name); + glp_set_col_kind(P, k, GLP_BV); + /* set objective coefficient (3) */ + glp_set_obj_coef(P, k, c[k]); + } + } + /* add all degree constraints (4) */ + ind = xalloc(1+n, sizeof(int)); + val = xalloc(1+n, sizeof(double)); + for (i = 1; i <= n; i++) + { k = glp_add_rows(P, 1); + xassert(k == i); + sprintf(name, "v[%d]", i); + glp_set_row_name(P, i, name); + glp_set_row_bnds(P, i, GLP_FX, 2, 2); + k = 0; + for (j = 1; j <= n; j++) + { if (i != j) + k++, ind[k] = loc(i,j), val[k] = 1; + } + xassert(k == n-1); + glp_set_mat_row(P, i, n-1, ind, val); + } + xfree(ind); + xfree(val); + return; +} + +/*********************************************************************** +* build_tour - build tour for corresponding solution to IP +* +* Given a feasible solution to IP (3)-(6) this routine builds the +* corresponding solution to TSP, which is a tour specified by the list +* of node numbers tour[1] -> ... -> tour[nn] -> tour[1] in the order +* the nodes are to be visited */ + +void build_tour(void) +{ int i, j, k, kk, *beg, *end; + /* solution to MIP should be feasible */ + switch (glp_mip_status(P)) + { case GLP_FEAS: + case GLP_OPT: + break; + default: + xassert(P != P); + } + /* build the list of edges included in the tour */ + beg = xalloc(1+n, sizeof(int)); + end = xalloc(1+n, sizeof(int)); + k = 0; + for (i = 1; i <= n; i++) + { for (j = i+1; j <= n; j++) + { double x; + x = glp_mip_col_val(P, loc(i,j)); + xassert(x == 0 || x == 1); + if (x) + { k++; + xassert(k <= n); + beg[k] = i, end[k] = j; + } + } + } + xassert(k == n); + /* reorder edges in the list as they follow in the tour */ + for (k = 1; k <= n; k++) + { /* find k-th edge of the tour */ + j = (k == 1 ? 1 : end[k-1]); + for (kk = k; kk <= n; kk++) + { if (beg[kk] == j) + break; + if (end[kk] == j) + { end[kk] = beg[kk], beg[kk] = j; + break; + } + } + xassert(kk <= n); + /* put the edge to k-th position in the list */ + i = beg[k], beg[k] = beg[kk], beg[kk] = i; + j = end[k], end[k] = end[kk], end[kk] = j; + } + /* build the tour starting from node 1 */ + xassert(beg[1] == 1); + for (k = 1; k <= n; k++) + { if (k > 1) + xassert(end[k-1] == beg[k]); + tour[k] = beg[k]; + } + xassert(end[n] == 1); + xfree(beg); + xfree(end); + return; +} + +/*********************************************************************** +* tour_length - calculate tour length +* +* This routine calculates the length of the specified tour, which is +* the sum of corresponding edge lengths. */ + +int tour_length(const int tour[/*1+n*/]) +{ int i, j, sum; + sum = 0; + for (i = 1; i <= n; i++) + { j = (i < n ? i+1 : 1); + sum += c[loc(tour[i], tour[j])]; + } + return sum; +} + +/*********************************************************************** +* write_tour - write tour to text file in TSPLIB format +* +* This routine writes the specified tour to a text file in TSPLIB +* format. */ + +void write_tour(const char *fname, const int tour[/*1+n*/]) +{ FILE *fp; + int i; + xprintf("Writing TSP solution to '%s'...\n", fname); + fp = fopen(fname, "w"); + if (fp == NULL) + { xprintf("Unable to create '%s' - %s\n", fname, + strerror(errno)); + return; + } + fprintf(fp, "NAME : %s\n", fname); + fprintf(fp, "COMMENT : Tour length is %d\n", tour_length(tour)); + fprintf(fp, "TYPE : TOUR\n"); + fprintf(fp, "DIMENSION : %d\n", n); + fprintf(fp, "TOUR_SECTION\n"); + for (i = 1; i <= n; i++) + fprintf(fp, "%d\n", tour[i]); + fprintf(fp, "-1\n"); + fprintf(fp, "EOF\n"); + fclose(fp); + return; +} + +/*********************************************************************** +* gen_subt_row - generate violated subtour elimination constraint +* +* This routine is called from the MIP solver to generate a violated +* subtour elimination constraint (5). +* +* Constraints of this class has the form: +* +* sum x[i,j] >= 2, i in W, j in V \ W, +* +* for all W, where W is a proper nonempty subset of V, V is the set of +* nodes of the given graph. +* +* In order to find a violated constraint of this class this routine +* finds a min cut in a capacitated network, which has the same sets of +* nodes and edges as the original graph, and where capacities of edges +* are values of variables x[i,j] in a basic solution to LP relaxation +* of the current subproblem. */ + +void gen_subt(glp_tree *T) +{ int i, j, ne, nz, *beg, *end, *cap, *cut, *ind; + double sum, *val; + /* MIP preprocessor should not be used */ + xassert(glp_ios_get_prob(T) == P); + /* if some variable x[i,j] is zero in basic solution, then the + * capacity of corresponding edge in the associated network is + * zero, so we may not include such edge in the network */ + /* count number of edges having non-zero capacity */ + ne = 0; + for (i = 1; i <= n; i++) + { for (j = i+1; j <= n; j++) + { if (glp_get_col_prim(P, loc(i,j)) >= .001) + ne++; + } + } + /* build the capacitated network */ + beg = xalloc(1+ne, sizeof(int)); + end = xalloc(1+ne, sizeof(int)); + cap = xalloc(1+ne, sizeof(int)); + nz = 0; + for (i = 1; i <= n; i++) + { for (j = i+1; j <= n; j++) + { if (glp_get_col_prim(P, loc(i,j)) >= .001) + { nz++; + xassert(nz <= ne); + beg[nz] = i, end[nz] = j; + /* scale all edge capacities to make them integral */ + cap[nz] = ceil(1000 * glp_get_col_prim(P, loc(i,j))); + } + } + } + xassert(nz == ne); + /* find minimal cut in the capacitated network */ + cut = xalloc(1+n, sizeof(int)); + min_cut(n, ne, beg, end, cap, cut); + /* determine the number of non-zero coefficients in the subtour + * elimination constraint and calculate its left-hand side which + * is the (unscaled) capacity of corresponding min cut */ + ne = 0, sum = 0; + for (i = 1; i <= n; i++) + { for (j = i+1; j <= n; j++) + { if (cut[i] && !cut[j] || !cut[i] && cut[j]) + { ne++; + sum += glp_get_col_prim(P, loc(i,j)); + } + } + } + /* if the (unscaled) capacity of min cut is less than 2, the + * corresponding subtour elimination constraint is violated */ + if (sum <= 1.999) + { /* build the list of non-zero coefficients */ + ind = xalloc(1+ne, sizeof(int)); + val = xalloc(1+ne, sizeof(double)); + nz = 0; + for (i = 1; i <= n; i++) + { for (j = i+1; j <= n; j++) + { if (cut[i] && !cut[j] || !cut[i] && cut[j]) + { nz++; + xassert(nz <= ne); + ind[nz] = loc(i,j); + val[nz] = 1; + } + } + } + xassert(nz == ne); + /* add violated tour elimination constraint to the current + * subproblem */ + i = glp_add_rows(P, 1); + glp_set_row_bnds(P, i, GLP_LO, 2, 0); + glp_set_mat_row(P, i, nz, ind, val); + xfree(ind); + xfree(val); + } + /* free working arrays */ + xfree(beg); + xfree(end); + xfree(cap); + xfree(cut); + return; +} + +/*********************************************************************** +* cb_func - application callback routine +* +* This routine is called from the MIP solver at various points of +* the branch-and-cut algorithm. */ + +void cb_func(glp_tree *T, void *info) +{ xassert(info == info); + switch (glp_ios_reason(T)) + { case GLP_IROWGEN: + /* generate one violated subtour elimination constraint */ + gen_subt(T); + break; + } + return; +} + +/*********************************************************************** +* main - TSP solver main program +* +* This main program parses command-line arguments, reads specified TSP +* instance from a text file, and calls the MIP solver to solve it. */ + +int main(int argc, char *argv[]) +{ int j; + char *in_file = NULL, *out_file = NULL; + time_t start; + glp_iocp iocp; + /* parse command-line arguments */ +# define p(str) (strcmp(argv[j], str) == 0) + for (j = 1; j < argc; j++) + { if (p("--output") || p("-o")) + { j++; + if (j == argc || argv[j][0] == '\0' || argv[j][0] == '-') + { xprintf("No solution output file specified\n"); + exit(EXIT_FAILURE); + } + if (out_file != NULL) + { xprintf("Only one solution output file allowed\n"); + exit(EXIT_FAILURE); + } + out_file = argv[j]; + } + else if (p("--help") || p("-h")) + { xprintf("Usage: %s [options...] tsp-file\n", argv[0]); + xprintf("\n"); + xprintf("Options:\n"); + xprintf(" -o filename, --output filename\n"); + xprintf(" write solution to filename\n") + ; + xprintf(" -h, --help display this help information" + " and exit\n"); + exit(EXIT_SUCCESS); + } + else if (argv[j][0] == '-' || + (argv[j][0] == '-' && argv[j][1] == '-')) + { xprintf("Invalid option '%s'; try %s --help\n", argv[j], + argv[0]); + exit(EXIT_FAILURE); + } + else + { if (in_file != NULL) + { xprintf("Only one input file allowed\n"); + exit(EXIT_FAILURE); + } + in_file = argv[j]; + } + } + if (in_file == NULL) + { xprintf("No input file specified; try %s --help\n", argv[0]); + exit(EXIT_FAILURE); + } +# undef p + /* display program banner */ + xprintf("TSP Solver for GLPK %s\n", glp_version()); + /* remove output solution file specified in command-line */ + if (out_file != NULL) + remove(out_file); + /* read TSP instance from input data file */ + read_data(in_file); + /* build initial IP problem */ + start = time(NULL); + build_prob(); + tour = xalloc(1+n, sizeof(int)); + /* solve LP relaxation of initial IP problem */ + xprintf("Solving initial LP relaxation...\n"); + xassert(glp_simplex(P, NULL) == 0); + xassert(glp_get_status(P) == GLP_OPT); + /* solve IP problem with "lazy" constraints */ + glp_init_iocp(&iocp); + iocp.br_tech = GLP_BR_MFV; /* most fractional variable */ + iocp.bt_tech = GLP_BT_BLB; /* best local bound */ + iocp.sr_heur = GLP_OFF; /* disable simple rounding heuristic */ + iocp.gmi_cuts = GLP_ON; /* enable Gomory cuts */ + iocp.cb_func = cb_func; + glp_intopt(P, &iocp); + build_tour(); + /* display some statistics */ + xprintf("Time used: %.1f secs\n", difftime(time(NULL), start)); + { size_t tpeak; + glp_mem_usage(NULL, NULL, NULL, &tpeak); + xprintf("Memory used: %.1f Mb (%.0f bytes)\n", + (double)tpeak / 1048576.0, (double)tpeak); + } + /* write solution to output file, if required */ + if (out_file != NULL) + write_tour(out_file, tour); + /* deallocate working objects */ + xfree(c); + xfree(tour); + glp_delete_prob(P); + /* check that no memory blocks are still allocated */ + { int count; + size_t total; + glp_mem_usage(&count, NULL, &total, NULL); + if (count != 0) + xerror("Error: %d memory block(s) were lost\n", count); + xassert(total == 0); + } + return 0; +} + +/* eof */ diff --git a/glpk-5.0/examples/tsp/maxflow.c b/glpk-5.0/examples/tsp/maxflow.c new file mode 100644 index 0000000..c28f527 --- /dev/null +++ b/glpk-5.0/examples/tsp/maxflow.c @@ -0,0 +1,170 @@ +/* maxflow.c */ + +/* Written by Andrew Makhorin <mao@gnu.org>, October 2015. */ + +#include <math.h> + +#include <glpk.h> +#include "maxflow.h" +#include "misc.h" + +/*********************************************************************** +* NAME +* +* max_flow - find max flow in undirected capacitated network +* +* SYNOPSIS +* +* #include "maxflow.h" +* int max_flow(int nn, int ne, const int beg[], const int end[], +* const int cap[], int s, int t, int x[]); +* +* DESCRIPTION +* +* This routine finds max flow in a given undirected network. +* +* The undirected capacitated network is specified by the parameters +* nn, ne, beg, end, and cap. The parameter nn specifies the number of +* vertices (nodes), nn >= 2, and the parameter ne specifies the number +* of edges, ne >= 0. The network edges are specified by triplets +* (beg[k], end[k], cap[k]) for k = 1, ..., ne, where beg[k] < end[k] +* are numbers of the first and second nodes of k-th edge, and +* cap[k] > 0 is a capacity of k-th edge. Loops and multiple edges are +* not allowed. +* +* The parameter s is the number of a source node, and the parameter t +* is the number of a sink node, s != t. +* +* On exit the routine computes elementary flows thru edges and stores +* their values to locations x[1], ..., x[ne]. Positive value of x[k] +* means that the elementary flow goes from node beg[k] to node end[k], +* and negative value means that the flow goes in opposite direction. +* +* RETURNS +* +* The routine returns the total maximum flow through the network. */ + +int max_flow(int nn, int ne, const int beg[/*1+ne*/], + const int end[/*1+ne*/], const int cap[/*1+ne*/], int s, int t, + int x[/*1+ne*/]) +{ int k; + /* sanity checks */ + xassert(nn >= 2); + xassert(ne >= 0); + xassert(1 <= s && s <= nn); + xassert(1 <= t && t <= nn); + xassert(s != t); + for (k = 1; k <= ne; k++) + { xassert(1 <= beg[k] && beg[k] < end[k] && end[k] <= nn); + xassert(cap[k] > 0); + } + /* find max flow */ + return max_flow_lp(nn, ne, beg, end, cap, s, t, x); +} + +/*********************************************************************** +* NAME +* +* max_flow_lp - find max flow with simplex method +* +* SYNOPSIS +* +* #include "maxflow.h" +* int max_flow_lp(int nn, int ne, const int beg[], const int end[], +* const int cap[], int s, int t, int x[]); +* +* DESCRIPTION +* +* This routine finds max flow in a given undirected network with the +* simplex method. +* +* Parameters of this routine have the same meaning as for the routine +* max_flow (see above). +* +* RETURNS +* +* The routine returns the total maximum flow through the network. */ + +int max_flow_lp(int nn, int ne, const int beg[/*1+ne*/], + const int end[/*1+ne*/], const int cap[/*1+ne*/], int s, int t, + int x[/*1+ne*/]) +{ glp_prob *lp; + glp_smcp smcp; + int i, k, nz, flow, *rn, *cn; + double temp, *aa; + /* create LP problem instance */ + lp = glp_create_prob(); + /* create LP rows; i-th row is the conservation condition of the + * flow at i-th node, i = 1, ..., nn */ + glp_add_rows(lp, nn); + for (i = 1; i <= nn; i++) + glp_set_row_bnds(lp, i, GLP_FX, 0.0, 0.0); + /* create LP columns; k-th column is the elementary flow thru + * k-th edge, k = 1, ..., ne; the last column with the number + * ne+1 is the total flow through the network, which goes along + * a dummy feedback edge from the sink to the source */ + glp_add_cols(lp, ne+1); + for (k = 1; k <= ne; k++) + { xassert(cap[k] > 0); + glp_set_col_bnds(lp, k, GLP_DB, -cap[k], +cap[k]); + } + glp_set_col_bnds(lp, ne+1, GLP_FR, 0.0, 0.0); + /* build the constraint matrix; structurally this matrix is the + * incidence matrix of the network, so each its column (including + * the last column for the dummy edge) has exactly two non-zero + * entries */ + rn = xalloc(1+2*(ne+1), sizeof(int)); + cn = xalloc(1+2*(ne+1), sizeof(int)); + aa = xalloc(1+2*(ne+1), sizeof(double)); + nz = 0; + for (k = 1; k <= ne; k++) + { /* x[k] > 0 means the elementary flow thru k-th edge goes from + * node beg[k] to node end[k] */ + nz++, rn[nz] = beg[k], cn[nz] = k, aa[nz] = -1.0; + nz++, rn[nz] = end[k], cn[nz] = k, aa[nz] = +1.0; + } + /* total flow thru the network goes from the sink to the source + * along the dummy feedback edge */ + nz++, rn[nz] = t, cn[nz] = ne+1, aa[nz] = -1.0; + nz++, rn[nz] = s, cn[nz] = ne+1, aa[nz] = +1.0; + /* check the number of non-zero entries */ + xassert(nz == 2*(ne+1)); + /* load the constraint matrix into the LP problem object */ + glp_load_matrix(lp, nz, rn, cn, aa); + xfree(rn); + xfree(cn); + xfree(aa); + /* objective function is the total flow through the network to + * be maximized */ + glp_set_obj_dir(lp, GLP_MAX); + glp_set_obj_coef(lp, ne + 1, 1.0); + /* solve LP instance with the (primal) simplex method */ + glp_term_out(0); + glp_adv_basis(lp, 0); + glp_term_out(1); + glp_init_smcp(&smcp); + smcp.msg_lev = GLP_MSG_ON; + smcp.out_dly = 5000; + xassert(glp_simplex(lp, &smcp) == 0); + xassert(glp_get_status(lp) == GLP_OPT); + /* obtain optimal elementary flows thru edges of the network */ + /* (note that the constraint matrix is unimodular and the data + * are integral, so all elementary flows in basic solution should + * also be integral) */ + for (k = 1; k <= ne; k++) + { temp = glp_get_col_prim(lp, k); + x[k] = (int)floor(temp + .5); + xassert(fabs(x[k] - temp) <= 1e-6); + } + /* obtain the maximum flow thru the original network which is the + * flow thru the dummy feedback edge */ + temp = glp_get_col_prim(lp, ne+1); + flow = (int)floor(temp + .5); + xassert(fabs(flow - temp) <= 1e-6); + /* delete LP problem instance */ + glp_delete_prob(lp); + /* return to the calling program */ + return flow; +} + +/* eof */ diff --git a/glpk-5.0/examples/tsp/maxflow.h b/glpk-5.0/examples/tsp/maxflow.h new file mode 100644 index 0000000..245c5ec --- /dev/null +++ b/glpk-5.0/examples/tsp/maxflow.h @@ -0,0 +1,20 @@ +/* maxflow.h */ + +/* Written by Andrew Makhorin <mao@gnu.org>, October 2015. */ + +#ifndef MAXFLOW_H +#define MAXFLOW_H + +int max_flow(int nn, int ne, const int beg[/*1+ne*/], + const int end[/*1+ne*/], const int cap[/*1+ne*/], int s, int t, + int x[/*1+ne*/]); +/* find max flow in undirected capacitated network */ + +int max_flow_lp(int nn, int ne, const int beg[/*1+ne*/], + const int end[/*1+ne*/], const int cap[/*1+ne*/], int s, int t, + int x[/*1+ne*/]); +/* find max flow with simplex method */ + +#endif + +/* eof */ diff --git a/glpk-5.0/examples/tsp/mincut.c b/glpk-5.0/examples/tsp/mincut.c new file mode 100644 index 0000000..225fb7f --- /dev/null +++ b/glpk-5.0/examples/tsp/mincut.c @@ -0,0 +1,392 @@ +/* mincut.c */ + +/* Written by Andrew Makhorin <mao@gnu.org>, October 2015. */ + +#include <limits.h> + +#include "maxflow.h" +#include "mincut.h" +#include "misc.h" + +/*********************************************************************** +* NAME +* +* min_cut - find min cut in undirected capacitated network +* +* SYNOPSIS +* +* #include "mincut.h" +* int min_cut(int nn, int ne, const int beg[], const int end[], +* const cap[], int cut[]); +* +* DESCRIPTION +* +* This routine finds min cut in a given undirected network. +* +* The undirected capacitated network is specified by the parameters +* nn, ne, beg, end, and cap. The parameter nn specifies the number of +* vertices (nodes), nn >= 2, and the parameter ne specifies the number +* of edges, ne >= 0. The network edges are specified by triplets +* (beg[k], end[k], cap[k]) for k = 1, ..., ne, where beg[k] < end[k] +* are numbers of the first and second nodes of k-th edge, and +* cap[k] > 0 is a capacity of k-th edge. Loops and multiple edges are +* not allowed. +* +* Let V be the set of nodes of the network and let W be an arbitrary +* non-empty proper subset of V. A cut associated with the subset W is +* a subset of all the edges, one node of which belongs to W and other +* node belongs to V \ W. The capacity of a cut (W, V \ W) is the sum +* of the capacities of all the edges, which belong to the cut. Minimal +* cut is a cut, whose capacity is minimal. +* +* On exit the routine stores flags of nodes v[i], i = 1, ..., nn, to +* locations cut[i], where cut[i] = 1 means that v[i] belongs to W and +* cut[i] = 0 means that v[i] belongs to V \ W, where W corresponds to +* the minimal cut found. +* +* RETURNS +* +* The routine returns the capacity of the min cut found. */ + +int min_cut(int nn, int ne, const int beg[/*1+ne*/], + const int end[/*1+ne*/], const cap[/*1+ne*/], int cut[/*1+nn*/]) +{ int k; + /* sanity checks */ + xassert(nn >= 2); + xassert(ne >= 0); + for (k = 1; k <= ne; k++) + { xassert(1 <= beg[k] && beg[k] < end[k] && end[k] <= nn); + xassert(cap[k] > 0); + } + /* find min cut */ + return min_cut_sw(nn, ne, beg, end, cap, cut); +} + +/*********************************************************************** +* NAME +* +* min_st_cut - find min (s,t)-cut for known max flow +* +* SYNOPSIS +* +* #include "mincut.h" +* +* DESCRIPTION +* +* This routine finds min (s,t)-cut in a given undirected network that +* corresponds to a known max flow from s to t in the network. +* +* Parameters nn, ne, beg, end, and cap specify the network in the same +* way as for the routine min_cut (see above). +* +* Parameters s and t specify, resp., the number of the source node +* and the number of the sink node, s != t, for which the min (s,t)-cut +* has to be found. +* +* Parameter x specifies the known max flow from s to t in the network, +* where locations x[1], ..., x[ne] contains elementary flow thru edges +* of the network (positive value of x[k] means that the elementary +* flow goes from node beg[k] to node end[k], and negative value means +* that the flow goes in opposite direction). +* +* This routine splits the set of nodes V of the network into two +* non-empty subsets V(s) and V(t) = V \ V(s), where the source node s +* belongs to V(s), the sink node t belongs to V(t), and all edges, one +* node of which belongs to V(s) and other one belongs to V(t), are +* saturated (that is, x[k] = +cap[k] or x[k] = -cap[k]). +* +* On exit the routine stores flags of the nodes v[i], i = 1, ..., nn, +* to locations cut[i], where cut[i] = 1 means that v[i] belongs to V(s) +* and cut[i] = 0 means that v[i] belongs to V(t) = V \ V(s). +* +* RETURNS +* +* The routine returns the capacity of min (s,t)-cut, which is the sum +* of the capacities of all the edges, which belong to the cut. (Note +* that due to theorem by Ford and Fulkerson this value is always equal +* to corresponding max flow.) +* +* ALGORITHM +* +* To determine the set V(s) the routine simply finds all nodes, which +* can be reached from the source node s via non-saturated edges. The +* set V(t) is determined as the complement V \ V(s). */ + +int min_st_cut(int nn, int ne, const int beg[/*1+ne*/], + const int end[/*1+ne*/], const int cap[/*1+ne*/], int s, int t, + const int x[/*1+ne*/], int cut[/*1+nn*/]) +{ int i, j, k, p, q, temp, *head1, *next1, *head2, *next2, *list; + /* head1[i] points to the first edge with beg[k] = i + * next1[k] points to the next edge with the same beg[k] + * head2[i] points to the first edge with end[k] = i + * next2[k] points to the next edge with the same end[k] */ + head1 = xalloc(1+nn, sizeof(int)); + head2 = xalloc(1+nn, sizeof(int)); + next1 = xalloc(1+ne, sizeof(int)); + next2 = xalloc(1+ne, sizeof(int)); + for (i = 1; i <= nn; i++) + head1[i] = head2[i] = 0; + for (k = 1; k <= ne; k++) + { i = beg[k], next1[k] = head1[i], head1[i] = k; + j = end[k], next2[k] = head2[j], head2[j] = k; + } + /* on constructing the set V(s) list[1], ..., list[p-1] contain + * nodes, which can be reached from source node and have been + * visited, and list[p], ..., list[q] contain nodes, which can be + * reached from source node but havn't been visited yet */ + list = xalloc(1+nn, sizeof(int)); + for (i = 1; i <= nn; i++) + cut[i] = 0; + p = q = 1, list[1] = s, cut[s] = 1; + while (p <= q) + { /* pick next node, which is reachable from the source node and + * has not visited yet, and visit it */ + i = list[p++]; + /* walk through edges with beg[k] = i */ + for (k = head1[i]; k != 0; k = next1[k]) + { j = end[k]; + xassert(beg[k] == i); + /* from v[i] we can reach v[j], if the elementary flow from + * v[i] to v[j] is non-saturated */ + if (cut[j] == 0 && x[k] < +cap[k]) + list[++q] = j, cut[j] = 1; + } + /* walk through edges with end[k] = i */ + for (k = head2[i]; k != 0; k = next2[k]) + { j = beg[k]; + xassert(end[k] == i); + /* from v[i] we can reach v[j], if the elementary flow from + * v[i] to v[j] is non-saturated */ + if (cut[j] == 0 && x[k] > -cap[k]) + list[++q] = j, cut[j] = 1; + } + } + /* sink cannot belong to V(s) */ + xassert(!cut[t]); + /* free working arrays */ + xfree(head1); + xfree(head2); + xfree(next1); + xfree(next2); + xfree(list); + /* compute capacity of the minimal (s,t)-cut found */ + temp = 0; + for (k = 1; k <= ne; k++) + { i = beg[k], j = end[k]; + if (cut[i] && !cut[j] || !cut[i] && cut[j]) + temp += cap[k]; + } + /* return to the calling program */ + return temp; +} + +/*********************************************************************** +* NAME +* +* min_cut_sw - find min cut with Stoer and Wagner algorithm +* +* SYNOPSIS +* +* #include "mincut.h" +* int min_cut_sw(int nn, int ne, const int beg[], const int end[], +* const cap[], int cut[]); +* +* DESCRIPTION +* +* This routine find min cut in a given undirected network with the +* algorithm proposed by Stoer and Wagner (see references below). +* +* Parameters of this routine have the same meaning as for the routine +* min_cut (see above). +* +* RETURNS +* +* The routine returns the capacity of the min cut found. +* +* ALGORITHM +* +* The basic idea of Stoer&Wagner algorithm is the following. Let G be +* a capacitated network, and G(s,t) be a network, in which the nodes s +* and t are merged into one new node, loops are deleted, but multiple +* edges are retained. It is obvious that a minimum cut in G is the +* minimum of two quantities: the minimum cut in G(s,t) and a minimum +* cut that separates s and t. This allows to find a minimum cut in the +* original network solving at most nn max flow problems. +* +* REFERENCES +* +* M. Stoer, F. Wagner. A Simple Min Cut Algorithm. Algorithms, ESA'94 +* LNCS 855 (1994), pp. 141-47. +* +* J. Cheriyan, R. Ravi. Approximation Algorithms for Network Problems. +* Univ. of Waterloo (1998), p. 147. */ + +int min_cut_sw(int nn, int ne, const int beg[/*1+ne*/], + const int end[/*1+ne*/], const cap[/*1+ne*/], int cut[/*1+nn*/]) +{ int i, j, k, min_cut, flow, temp, *head1, *next1, *head2, *next2; + int I, J, K, S, T, DEG, NV, NE, *HEAD, *NEXT, *NUMB, *BEG, *END, + *CAP, *X, *ADJ, *SUM, *CUT; + /* head1[i] points to the first edge with beg[k] = i + * next1[k] points to the next edge with the same beg[k] + * head2[i] points to the first edge with end[k] = i + * next2[k] points to the next edge with the same end[k] */ + head1 = xalloc(1+nn, sizeof(int)); + head2 = xalloc(1+nn, sizeof(int)); + next1 = xalloc(1+ne, sizeof(int)); + next2 = xalloc(1+ne, sizeof(int)); + for (i = 1; i <= nn; i++) + head1[i] = head2[i] = 0; + for (k = 1; k <= ne; k++) + { i = beg[k], next1[k] = head1[i], head1[i] = k; + j = end[k], next2[k] = head2[j], head2[j] = k; + } + /* an auxiliary network used in the algorithm is resulted from + * the original network by merging some nodes into one supernode; + * all variables and arrays related to this auxiliary network are + * denoted in CAPS */ + /* HEAD[I] points to the first node of the original network that + * belongs to the I-th supernode + * NEXT[i] points to the next node of the original network that + * belongs to the same supernode as the i-th node + * NUMB[i] is a supernode, which the i-th node belongs to */ + /* initially the auxiliary network is equivalent to the original + * network, i.e. each supernode consists of one node */ + NV = nn; + HEAD = xalloc(1+nn, sizeof(int)); + NEXT = xalloc(1+nn, sizeof(int)); + NUMB = xalloc(1+nn, sizeof(int)); + for (i = 1; i <= nn; i++) + HEAD[i] = i, NEXT[i] = 0, NUMB[i] = i; + /* number of edges in the auxiliary network is never greater than + * in the original one */ + BEG = xalloc(1+ne, sizeof(int)); + END = xalloc(1+ne, sizeof(int)); + CAP = xalloc(1+ne, sizeof(int)); + X = xalloc(1+ne, sizeof(int)); + /* allocate some auxiliary arrays */ + ADJ = xalloc(1+nn, sizeof(int)); + SUM = xalloc(1+nn, sizeof(int)); + CUT = xalloc(1+nn, sizeof(int)); + /* currently no min cut is found so far */ + min_cut = INT_MAX; + /* main loop starts here */ + while (NV > 1) + { /* build the set of edges of the auxiliary network */ + NE = 0; + /* multiple edges are not allowed in the max flow algorithm, + * so we can replace each multiple edge, which is the result + * of merging nodes into supernodes, by a single edge, whose + * capacity is the sum of capacities of particular edges; + * these summary capacities will be stored in the array SUM */ + for (I = 1; I <= NV; I++) + SUM[I] = 0.0; + for (I = 1; I <= NV; I++) + { /* DEG is number of single edges, which connects I-th + * supernode and some J-th supernode, where I < J */ + DEG = 0; + /* walk thru nodes that belong to I-th supernode */ + for (i = HEAD[I]; i != 0; i = NEXT[i]) + { /* i-th node belongs to I-th supernode */ + /* walk through edges with beg[k] = i */ + for (k = head1[i]; k != 0; k = next1[k]) + { j = end[k]; + /* j-th node belongs to J-th supernode */ + J = NUMB[j]; + /* ignore loops and edges with I > J */ + if (I >= J) + continue; + /* add an edge that connects I-th and J-th supernodes + * (if not added yet) */ + if (SUM[J] == 0.0) + ADJ[++DEG] = J; + /* sum up the capacity of the original edge */ + xassert(cap[k] > 0.0); + SUM[J] += cap[k]; + } + /* walk through edges with end[k] = i */ + for (k = head2[i]; k != 0; k = next2[k]) + { j = beg[k]; + /* j-th node belongs to J-th supernode */ + J = NUMB[j]; + /* ignore loops and edges with I > J */ + if (I >= J) + continue; + /* add an edge that connects I-th and J-th supernodes + * (if not added yet) */ + if (SUM[J] == 0.0) + ADJ[++DEG] = J; + /* sum up the capacity of the original edge */ + xassert(cap[k] > 0.0); + SUM[J] += cap[k]; + } + } + /* add single edges connecting I-th supernode with other + * supernodes to the auxiliary network; restore the array + * SUM for subsequent use */ + for (K = 1; K <= DEG; K++) + { NE++; + xassert(NE <= ne); + J = ADJ[K]; + BEG[NE] = I, END[NE] = J, CAP[NE] = SUM[J]; + SUM[J] = 0.0; + } + } + /* choose two arbitrary supernodes of the auxiliary network, + * one of which is the source and other is the sink */ + S = 1, T = NV; + /* determine max flow from S to T */ + flow = max_flow(NV, NE, BEG, END, CAP, S, T, X); + /* if the min cut that separates supernodes S and T is less + * than the currently known, remember it */ + if (min_cut > flow) + { min_cut = flow; + /* find min (s,t)-cut in the auxiliary network */ + temp = min_st_cut(NV, NE, BEG, END, CAP, S, T, X, CUT); + /* (Ford and Fulkerson insist on this) */ + xassert(flow == temp); + /* build corresponding min cut in the original network */ + for (i = 1; i <= nn; i++) cut[i] = CUT[NUMB[i]]; + /* if the min cut capacity is zero (i.e. the network has + * unconnected components), the search can be prematurely + * terminated */ + if (min_cut == 0) + break; + } + /* now merge all nodes of the original network, which belong + * to the supernodes S and T, into one new supernode; this is + * attained by carrying all nodes from T to S (for the sake of + * convenience T should be the last supernode) */ + xassert(T == NV); + /* assign new references to nodes from T */ + for (i = HEAD[T]; i != 0; i = NEXT[i]) + NUMB[i] = S; + /* find last entry in the node list of S */ + i = HEAD[S]; + xassert(i != 0); + while (NEXT[i] != 0) + i = NEXT[i]; + /* and attach to it the node list of T */ + NEXT[i] = HEAD[T]; + /* decrease number of nodes in the auxiliary network */ + NV--; + } + /* free working arrays */ + xfree(HEAD); + xfree(NEXT); + xfree(NUMB); + xfree(BEG); + xfree(END); + xfree(CAP); + xfree(X); + xfree(ADJ); + xfree(SUM); + xfree(CUT); + xfree(head1); + xfree(head2); + xfree(next1); + xfree(next2); + /* return to the calling program */ + return min_cut; +} + +/* eof */ diff --git a/glpk-5.0/examples/tsp/mincut.h b/glpk-5.0/examples/tsp/mincut.h new file mode 100644 index 0000000..aefdbb7 --- /dev/null +++ b/glpk-5.0/examples/tsp/mincut.h @@ -0,0 +1,23 @@ +/* mincut.c */ + +/* Written by Andrew Makhorin <mao@gnu.org>, October 2015. */ + +#ifndef MINCUT_H +#define MINCUT_H + +int min_cut(int nn, int ne, const int beg[/*1+ne*/], + const int end[/*1+ne*/], const cap[/*1+ne*/], int cut[/*1+nn*/]); +/* find min cut in undirected capacitated network */ + +int min_st_cut(int nn, int ne, const int beg[/*1+ne*/], + const int end[/*1+ne*/], const int cap[/*1+ne*/], int s, int t, + const int x[/*1+ne*/], int cut[/*1+nn*/]); +/* find min (s,t)-cut for known max flow */ + +int min_cut_sw(int nn, int ne, const int beg[/*1+ne*/], + const int end[/*1+ne*/], const cap[/*1+ne*/], int cut[/*1+nn*/]); +/* find min cut with Stoer and Wagner algorithm */ + +#endif + +/* eof */ diff --git a/glpk-5.0/examples/tsp/misc.c b/glpk-5.0/examples/tsp/misc.c new file mode 100644 index 0000000..e0b08fb --- /dev/null +++ b/glpk-5.0/examples/tsp/misc.c @@ -0,0 +1,159 @@ +/* misc.c */ + +/* Written by Andrew Makhorin <mao@gnu.org>, October 2015. */ + +#include <ctype.h> +#include <float.h> +#include <limits.h> +#include <stdlib.h> +#include "misc.h" + +/*********************************************************************** +* NAME +* +* str2int - convert character string to value of int type +* +* SYNOPSIS +* +* #include "misc.h" +* int str2int(const char *str, int *val); +* +* DESCRIPTION +* +* The routine str2int converts the character string str to a value of +* integer type and stores the value into location, which the parameter +* val points to (in the case of error content of this location is not +* changed). +* +* RETURNS +* +* The routine returns one of the following error codes: +* +* 0 - no error; +* 1 - value out of range; +* 2 - character string is syntactically incorrect. */ + +int str2int(const char *str, int *val_) +{ int d, k, s, val = 0; + /* scan optional sign */ + if (str[0] == '+') + s = +1, k = 1; + else if (str[0] == '-') + s = -1, k = 1; + else + s = +1, k = 0; + /* check for the first digit */ + if (!isdigit((unsigned char)str[k])) + return 2; + /* scan digits */ + while (isdigit((unsigned char)str[k])) + { d = str[k++] - '0'; + if (s > 0) + { if (val > INT_MAX / 10) + return 1; + val *= 10; + if (val > INT_MAX - d) + return 1; + val += d; + } + else /* s < 0 */ + { if (val < INT_MIN / 10) + return 1; + val *= 10; + if (val < INT_MIN + d) + return 1; + val -= d; + } + } + /* check for terminator */ + if (str[k] != '\0') + return 2; + /* conversion has been done */ + *val_ = val; + return 0; +} + +/*********************************************************************** +* NAME +* +* str2num - convert character string to value of double type +* +* SYNOPSIS +* +* #include "misc.h" +* int str2num(const char *str, double *val); +* +* DESCRIPTION +* +* The routine str2num converts the character string str to a value of +* double type and stores the value into location, which the parameter +* val points to (in the case of error content of this location is not +* changed). +* +* RETURNS +* +* The routine returns one of the following error codes: +* +* 0 - no error; +* 1 - value out of range; +* 2 - character string is syntactically incorrect. */ + +int str2num(const char *str, double *val_) +{ int k; + double val; + /* scan optional sign */ + k = (str[0] == '+' || str[0] == '-' ? 1 : 0); + /* check for decimal point */ + if (str[k] == '.') + { k++; + /* a digit should follow it */ + if (!isdigit((unsigned char)str[k])) + return 2; + k++; + goto frac; + } + /* integer part should start with a digit */ + if (!isdigit((unsigned char)str[k])) + return 2; + /* scan integer part */ + while (isdigit((unsigned char)str[k])) + k++; + /* check for decimal point */ + if (str[k] == '.') k++; +frac: /* scan optional fraction part */ + while (isdigit((unsigned char)str[k])) + k++; + /* check for decimal exponent */ + if (str[k] == 'E' || str[k] == 'e') + { k++; + /* scan optional sign */ + if (str[k] == '+' || str[k] == '-') + k++; + /* a digit should follow E, E+ or E- */ + if (!isdigit((unsigned char)str[k])) + return 2; + } + /* scan optional exponent part */ + while (isdigit((unsigned char)str[k])) + k++; + /* check for terminator */ + if (str[k] != '\0') + return 2; + /* perform conversion */ + { char *endptr; + val = strtod(str, &endptr); + if (*endptr != '\0') + return 2; + } + /* check for overflow */ + if (!(-DBL_MAX <= val && val <= +DBL_MAX)) + return 1; + /* check for underflow */ + if (-DBL_MIN < val && val < +DBL_MIN) + val = 0.0; + /* conversion has been done */ + *val_ = val; + return 0; +} + +/* eof */ diff --git a/glpk-5.0/examples/tsp/misc.h b/glpk-5.0/examples/tsp/misc.h new file mode 100644 index 0000000..0973831 --- /dev/null +++ b/glpk-5.0/examples/tsp/misc.h @@ -0,0 +1,24 @@ +/* misc.h */ + +/* Written by Andrew Makhorin <mao@gnu.org>, October 2015. */ + +#ifndef MISC_H +#define MISC_H + +#include <glpk.h> + +#define xprintf glp_printf +#define xerror glp_error +#define xassert glp_assert +#define xalloc glp_alloc +#define xfree glp_free + +int str2int(const char *str, int *val); +/* convert character string to value of int type */ + +int str2num(const char *str, double *val); +/* convert character string to value of double type */ + +#endif + +/* eof */ diff --git a/glpk-5.0/examples/tsp/moscow.tsp b/glpk-5.0/examples/tsp/moscow.tsp new file mode 100644 index 0000000..a2d4ea2 --- /dev/null +++ b/glpk-5.0/examples/tsp/moscow.tsp @@ -0,0 +1,200 @@ +NAME: MOSCOW +TYPE: TSP +COMMENT: 68 cities in Moscow region (Andrew Makhorin <mao@gnu.org>) +DIMENSION: 68 +EDGE_WEIGHT_TYPE: EXPLICIT +EDGE_WEIGHT_FORMAT: LOWER_DIAG_ROW +EDGE_WEIGHT_SECTION +0 +80 0 +99 66 0 +73 153 156 0 +118 156 205 116 0 +136 103 37 190 242 0 +134 128 180 141 90 217 0 +72 75 127 126 81 164 86 0 +131 94 146 184 133 183 67 91 0 +75 38 90 136 126 127 94 45 56 0 +67 52 46 124 169 83 165 88 132 76 0 +134 66 69 207 219 68 189 138 136 101 106 0 +179 142 194 204 153 231 87 139 52 104 180 188 0 +145 98 63 216 240 26 212 159 178 122 109 53 226 0 +83 7 69 156 159 106 131 78 97 41 55 60 145 95 0 +85 38 68 158 180 105 152 99 118 62 57 74 166 100 41 0 +188 163 97 222 294 77 277 223 243 187 135 145 291 93 166 165 0 +38 101 125 92 80 162 96 34 101 71 93 164 149 176 104 116 214 0 +88 66 118 142 97 155 62 40 67 36 103 129 115 150 69 90 215 50 0 +113 50 102 174 164 139 136 83 92 46 88 94 140 134 53 74 199 109 74 0 +67 94 146 113 62 183 81 19 86 64 107 157 134 178 97 118 242 29 35 102 0 +133 132 96 167 239 91 245 168 212 156 80 159 260 117 135 137 55 159 183 + 168 187 0 +51 64 54 108 157 91 159 82 139 83 25 118 187 117 67 69 143 77 97 100 101 + 88 0 +122 116 168 129 78 205 12 74 55 82 153 177 75 200 119 140 265 84 50 124 + 69 233 147 0 +148 115 49 202 254 29 229 176 195 139 95 97 243 45 118 117 48 174 167 + 151 195 83 103 217 0 +95 32 84 156 146 121 118 65 84 28 70 95 132 116 35 56 181 91 56 18 84 + 150 82 106 133 0 +133 70 122 194 184 159 156 103 112 66 108 99 160 148 73 94 219 129 94 25 + 122 188 120 144 171 38 0 +57 58 110 118 98 147 103 17 84 28 71 121 132 142 61 82 206 43 42 66 36 + 151 65 91 159 48 86 0 +164 101 142 225 188 167 122 134 55 97 139 109 107 152 104 125 239 150 + 116 67 135 219 151 110 191 69 87 117 0 +145 93 78 218 240 41 212 159 153 122 117 27 205 26 87 100 118 176 150 + 121 178 132 129 200 70 116 126 142 126 0 +142 74 77 215 227 60 197 146 134 109 114 8 186 45 68 82 137 172 137 102 + 165 151 126 185 89 103 107 129 107 19 0 +85 48 100 146 119 137 84 55 63 10 86 111 111 132 51 72 197 72 38 56 57 + 166 93 72 149 38 76 38 107 132 119 0 +106 26 73 179 182 110 154 101 120 64 78 48 168 97 30 62 170 127 92 76 + 120 158 90 142 122 58 84 84 102 75 56 74 0 +171 138 72 225 277 52 252 199 218 162 118 120 266 68 141 140 25 197 190 + 174 218 80 126 240 23 156 194 182 214 93 112 172 145 0 +65 18 48 138 160 85 132 79 98 42 37 69 146 80 21 20 145 96 70 54 98 117 + 49 120 97 36 74 62 105 80 77 52 44 120 0 +90 163 173 17 99 207 124 111 167 138 141 217 187 228 166 168 239 77 127 + 176 96 184 125 112 219 158 196 120 222 228 225 148 189 242 148 0 +47 22 67 116 121 104 106 40 80 24 44 85 128 99 25 39 164 66 44 46 59 124 + 46 94 116 28 66 23 97 99 93 34 48 139 19 118 0 +89 26 78 150 140 115 112 59 78 22 64 89 126 110 29 50 175 85 50 24 78 + 144 76 100 127 6 44 42 75 110 97 32 52 150 30 152 22 0 +31 111 118 42 134 155 153 91 157 101 86 165 205 176 114 116 207 57 107 + 139 86 152 70 141 167 121 159 83 190 176 173 111 137 190 96 59 78 115 0 +115 35 56 188 191 93 158 110 111 73 87 31 159 80 42 45 153 136 101 63 + 129 152 99 146 105 65 68 93 86 58 39 83 17 128 53 198 57 61 146 0 +37 70 103 91 121 140 121 44 101 45 62 124 149 135 73 75 197 41 59 83 63 + 142 56 109 152 65 103 27 134 135 132 55 96 175 55 93 25 59 56 105 0 +143 142 106 177 249 101 255 178 222 166 90 169 270 127 145 147 65 169 + 193 178 197 10 98 243 93 160 198 161 229 142 161 176 168 90 127 194 134 + 154 162 162 152 0 +152 151 115 186 258 110 264 187 231 175 99 178 279 136 154 156 100 178 + 202 187 206 45 107 252 98 169 207 170 238 151 170 185 177 121 136 203 + 143 163 171 171 161 55 0 +146 74 81 219 230 70 191 149 124 112 118 12 176 55 72 84 147 175 140 102 + 168 161 130 179 99 104 107 132 97 29 10 122 56 122 81 229 96 100 177 39 + 136 171 180 0 +118 50 53 191 203 84 173 122 126 85 90 16 174 69 44 58 150 148 113 78 + 141 149 102 161 102 79 83 105 101 43 24 95 32 125 53 201 69 73 149 15 + 108 159 168 28 0 +50 56 62 107 155 99 151 74 131 75 17 110 179 121 59 61 151 76 89 92 93 + 96 8 139 111 74 112 57 143 121 118 85 82 134 41 124 38 68 69 91 48 106 + 115 122 94 0 +107 126 110 113 209 117 221 144 201 145 87 179 249 143 129 131 149 133 + 159 162 162 94 62 209 129 144 182 127 213 158 177 155 152 152 111 130 + 108 138 98 161 118 104 113 187 163 70 0 +105 42 94 166 156 131 128 75 84 38 80 86 132 126 45 66 191 101 66 8 94 + 160 92 116 143 10 28 58 59 113 94 48 68 166 46 168 38 16 131 55 75 170 + 179 94 70 84 154 0 +114 133 117 120 216 124 228 151 208 152 94 186 256 150 136 138 156 140 + 166 169 169 101 69 216 136 151 189 134 220 165 184 162 159 159 118 137 + 115 145 105 168 125 111 120 194 170 77 35 161 0 +89 42 72 162 184 109 156 103 122 66 61 70 170 104 45 4 169 120 94 78 122 + 141 73 144 121 60 98 86 127 97 78 76 58 144 24 172 43 54 120 41 79 151 + 160 80 54 65 135 70 142 0 +218 166 151 291 313 114 285 232 226 195 190 100 278 99 160 173 191 249 + 223 194 251 205 202 273 143 189 199 215 199 73 92 205 148 166 153 301 + 172 183 249 131 208 215 224 102 116 194 231 186 238 170 0 +81 147 168 48 68 205 93 80 136 117 136 208 156 219 150 159 250 46 96 155 + 65 195 120 81 217 137 175 89 191 219 216 118 173 240 139 31 109 131 66 + 182 84 205 214 220 192 119 141 147 148 163 292 0 +147 84 125 208 180 160 115 117 48 80 122 100 100 145 87 108 222 133 99 + 50 118 202 134 103 174 52 70 100 17 119 100 90 85 197 88 210 80 58 173 + 69 117 212 221 90 84 126 196 42 203 110 192 179 0 +93 112 96 99 195 103 207 130 187 131 73 165 235 129 115 117 135 119 145 + 148 148 80 48 195 115 130 168 113 199 144 163 141 138 138 97 116 94 124 + 84 147 104 90 99 173 149 56 14 140 21 121 217 127 182 0 +97 91 143 151 103 180 37 49 76 57 128 152 100 175 94 115 240 59 25 99 44 + 208 122 25 192 81 119 66 125 175 160 47 117 215 95 136 69 75 116 121 84 + 218 227 160 136 114 184 91 191 119 248 105 108 170 0 +122 121 85 156 228 80 234 157 201 145 69 148 249 106 124 126 70 148 172 + 157 176 15 77 222 68 139 177 140 208 121 140 155 147 91 106 173 113 133 + 141 141 131 25 30 150 138 85 83 149 90 130 194 184 191 69 197 0 +179 142 194 232 181 222 115 139 48 104 180 164 38 207 145 166 291 149 + 115 126 134 260 187 103 243 128 146 132 83 181 162 111 161 266 146 215 + 128 126 205 145 149 270 279 152 160 179 249 118 256 170 254 184 76 235 + 124 249 0 +38 67 82 99 144 119 146 69 126 70 50 121 174 132 70 72 171 64 84 103 88 + 116 34 134 131 85 123 52 154 132 129 80 93 154 52 116 40 79 57 102 43 + 126 135 133 105 33 90 95 97 76 205 107 137 76 109 105 174 0 +100 34 86 173 168 123 140 87 100 50 72 70 148 118 37 58 183 113 78 52 + 106 152 84 128 135 44 57 70 75 97 78 60 27 158 38 180 39 38 131 39 87 + 162 171 78 54 76 146 44 153 62 170 159 58 132 103 141 134 87 0 +69 46 98 130 117 135 82 39 72 16 83 109 120 130 49 70 195 65 20 54 55 + 163 77 70 147 36 74 22 105 130 117 26 72 170 50 132 24 30 95 81 39 173 + 182 120 93 69 139 46 146 74 203 111 88 125 45 152 120 64 58 0 +160 97 138 219 168 173 102 126 35 91 135 113 87 158 100 121 235 136 102 + 63 121 215 147 90 187 65 83 113 30 132 113 98 98 210 101 202 93 71 186 + 82 130 225 234 103 97 139 209 55 216 123 205 171 13 195 111 204 83 150 + 71 101 0 +70 89 73 110 176 104 184 107 164 108 50 142 212 130 92 94 136 96 122 125 + 125 81 25 172 116 107 145 90 176 145 150 118 115 139 74 127 71 101 89 + 124 81 91 100 154 126 33 37 117 44 98 218 138 159 23 147 70 212 53 109 + 102 172 0 +188 136 121 261 283 84 255 202 196 165 160 70 248 69 130 143 161 219 193 + 164 221 175 172 243 113 159 169 185 169 43 62 175 118 136 123 271 142 + 153 219 101 178 185 194 72 86 164 201 156 208 140 30 262 162 187 218 + 164 224 175 140 173 175 188 0 +94 28 80 167 162 117 134 81 100 44 66 69 148 112 31 52 177 107 72 56 100 + 146 78 122 129 38 63 64 81 96 77 54 21 152 32 174 33 32 125 38 81 156 + 165 77 53 70 140 48 147 56 169 153 64 126 97 135 140 81 6 52 77 103 139 + 0 +45 43 75 118 142 112 138 61 118 62 30 97 166 108 46 48 164 76 76 79 80 + 109 21 126 124 61 99 44 130 108 105 72 69 147 28 128 25 55 76 78 35 119 + 128 109 81 13 83 71 90 52 181 119 113 69 101 98 166 32 63 56 126 46 151 + 57 0 +113 35 48 186 191 85 163 110 119 73 85 37 167 72 39 37 145 136 101 71 + 129 144 97 151 97 67 76 93 94 64 45 83 25 120 48 196 57 61 144 8 103 + 154 163 47 21 89 158 63 165 33 137 182 77 144 126 133 153 100 47 81 90 + 121 107 46 76 0 +97 21 57 170 173 94 145 92 111 55 69 46 159 81 14 55 154 118 83 67 111 + 149 81 133 106 49 87 75 118 73 54 65 18 129 35 180 39 43 128 33 87 159 + 168 58 30 73 143 59 150 58 146 164 101 129 108 138 159 84 45 63 114 106 + 116 39 60 25 0 +116 79 131 177 133 168 98 76 51 41 117 127 63 163 82 103 228 86 52 77 71 + 197 124 86 180 69 97 69 100 154 135 48 105 203 83 163 65 63 142 96 86 + 207 216 135 111 116 186 69 193 107 227 132 83 172 61 186 63 111 85 57 + 86 149 197 85 103 104 96 0 +EOF + +Vertices of the graph: + 1 Aprelevka 35 Lyubertsy + 2 Balashikha 36 Mozhaysk + 3 Bronnitsy 37 Moskva + 4 Vereya 38 Mytischi + 5 Volokolamsk 39 Naro-Fominsk + 6 Voskresensk 40 Noginsk + 7 Vysokovsk 41 Odintsovo + 8 Dedovsk 42 Ozherel'ye + 9 Dmitrov 43 Ozyory +10 Dolgoprudny 44 Orekhovo-Zuevo +11 Domodedovo 45 Pavlovsky Posad +12 Drezna 46 Podol'sk +13 Dubna 47 Protvino +14 Egor'yevsk 48 Pushkino +15 Zheleznodorozhny 49 Puschino +16 Zhukovsky 50 Ramenskoye +17 Zaraysk 51 Roshal +18 Zvenigorod 52 Ruza +19 Zelenograd 53 Sergiyev Posad +20 Ivanteyevka 54 Serpukhov +21 Istra 55 Solnechnogorsk +22 Kashira 56 Stupino +23 Klimovsk 57 Taldom +24 Klin 58 Troitsk +25 Kolomna 59 Fryazino +26 Korolyov (Podlipki) 60 Khimki +27 Krasnoarmeysk 61 Khot'kovo +28 Krasnogorsk 62 Chekhov +29 Krasnozavodsk 63 Shatura +30 Kurovskoye 64 Schyolkovo +31 Likino-Dulyovo 65 Scherbinka +32 Lobnya 66 Elektrostal +33 Losino-Petrovsky 67 Elektrougli +34 Lukhovitsy 68 Yakhroma + +Optimal tour length is 1994 kilometers (obtained by glpk tspsol): +1 39 4 36 52 5 7 24 55 19 60 10 32 68 13 57 9 61 53 29 27 20 48 26 38 +59 64 33 67 15 2 35 16 50 66 40 45 12 44 31 63 51 30 14 6 3 25 34 17 22 +42 56 43 47 49 54 62 23 46 11 65 58 41 37 28 8 21 18 1 diff --git a/glpk-5.0/examples/tsp/sample.tsp b/glpk-5.0/examples/tsp/sample.tsp new file mode 100644 index 0000000..9e1367a --- /dev/null +++ b/glpk-5.0/examples/tsp/sample.tsp @@ -0,0 +1,16 @@ +NAME: SAMPLE +TYPE: TSP +COMMENT: Example from D.Phillips, A.Garcia-Diaz, p.124 +DIMENSION: 8 +EDGE_WEIGHT_TYPE: EXPLICIT +EDGE_WEIGHT_FORMAT: LOWER_DIAG_ROW +EDGE_WEIGHT_SECTION +0 +190 0 +210 380 0 +680 760 890 0 +690 790 900 480 0 +460 610 340 760 890 0 +780 670 410 510 490 720 0 +750 450 600 250 560 600 500 0 +EOF diff --git a/glpk-5.0/examples/tsp/tsplib.c b/glpk-5.0/examples/tsp/tsplib.c new file mode 100644 index 0000000..189ff8a --- /dev/null +++ b/glpk-5.0/examples/tsp/tsplib.c @@ -0,0 +1,730 @@ +/* tsplib.c */ + +/* Written by Andrew Makhorin <mao@gnu.org>, October 2015. */ + +#include <ctype.h> +#include <errno.h> +#include <float.h> +#include <math.h> +#include <stdio.h> +#include <stdlib.h> +#include <string.h> + +#include "misc.h" +#include "tsplib.h" + +/*********************************************************************** +* NAME +* +* tsp_read_data - read TSP instance data +* +* SYNOPSIS +* +* #include "tsplib.h" +* TSP *tsp_read_data(const char *fname); +* +* DESCRIPTION +* +* The routine tsp_read_data reads a TSP (or related problem) instance +* data from a text file, whose name is the character string fname. +* +* For detailed description of the format recognized by the routine see +* the report: G.Reinelt, TSPLIB 95. +* +* RETURNS +* +* If no error occurred, the routine tsp_read_data returns a pointer to +* the TSP instance data block, which contains loaded data. In the case +* of error the routine prints an error message and returns NULL. */ + +struct csa +{ /* common storage area used by the routine tsp_read_data */ + const char *fname; + /* name of the input text file */ + FILE *fp; + /* stream assigned to the input text file */ + int seqn; + /* line sequential number */ + int c; + /* current character */ + char token[255+1]; + /* current token */ +}; + +static int get_char(struct csa *csa) +{ csa->c = fgetc(csa->fp); + if (ferror(csa->fp)) + { xprintf("%s:%d: read error - %s\n", csa->fname, csa->seqn, + strerror(errno)); + return 1; + } + if (feof(csa->fp)) + csa->c = EOF; + else if (csa->c == '\n') + csa->seqn++; + else if (isspace(csa->c)) + csa->c = ' '; + else if (iscntrl(csa->c)) + { xprintf("%s:%d: invalid control character 0x%02X\n", + csa->fname, csa->seqn, csa->c); + return 1; + } + return 0; +} + +static int skip_spaces(struct csa *csa, int across) +{ while (csa->c == ' ' || (across && csa->c == '\n')) + if (get_char(csa)) return 1; + return 0; +} + +static int scan_keyword(struct csa *csa) +{ int len = 0; + if (skip_spaces(csa, 0)) + return 1; + if (csa->c == EOF) + { xprintf("%s:%d: warning: missing EOF inserted\n", csa->fname, + csa->seqn); + strcpy(csa->token, "EOF"); + return 0; + } + csa->token[0] = '\0'; + while (isalnum(csa->c) || csa->c == '_') + { if (len == 31) + { xprintf("%s:%d: keyword '%s...' too long\n", csa->fname, + csa->seqn, csa->token); + return 1; + } + csa->token[len++] = (char)csa->c, csa->token[len] = '\0'; + if (get_char(csa)) + return 1; + } + if (len == 0) + { xprintf("%s:%d: missing keyword\n", csa->fname, csa->seqn); + return 1; + } + return 0; +} + +static int check_colon(struct csa *csa) +{ if (skip_spaces(csa, 0)) + return 1; + if (csa->c != ':') + { xprintf("%s:%d: missing colon after '%s'\n", csa->fname, + csa->seqn, csa->token); + return 1; + } + if (get_char(csa)) + return 1; + return 0; +} + +static int scan_token(struct csa *csa, int across) +{ int len = 0; + if (skip_spaces(csa, across)) + return 1; + csa->token[0] = '\0'; + while (!(csa->c == EOF || csa->c == '\n' || csa->c == ' ')) + { if (len == 255) + { csa->token[31] = '\0'; + xprintf("%s:%d: token '%s...' too long\n", csa->fname, + csa->seqn, csa->token); + return 1; + } + csa->token[len++] = (char)csa->c, csa->token[len] = '\0'; + if (get_char(csa)) + return 1; + } + return 0; +} + +static int check_newline(struct csa *csa) +{ if (skip_spaces(csa, 0)) + return 1; + if (!(csa->c == EOF || csa->c == '\n')) + { xprintf("%s:%d: extra symbols detected\n", csa->fname, + csa->seqn); + return 1; + } + if (get_char(csa)) + return 1; + return 0; +} + +static int scan_comment(struct csa *csa) +{ int len = 0; + if (skip_spaces(csa, 0)) + return 1; + csa->token[0] = '\0'; + while (!(csa->c == EOF || csa->c == '\n')) + { if (len == 255) + { xprintf("%s:%d: comment too long\n", csa->fname, csa->seqn); + return 1; + } + csa->token[len++] = (char)csa->c, csa->token[len] = '\0'; + if (get_char(csa)) + return 1; + } + return 0; +} + +static int scan_integer(struct csa *csa, int across, int *val) +{ if (scan_token(csa, across)) + return 1; + if (strlen(csa->token) == 0) + { xprintf("%s:%d: missing integer\n", csa->fname, csa->seqn); + return 1; + } + if (str2int(csa->token, val)) + { xprintf("%s:%d: integer '%s' invalid\n", csa->fname, csa->seqn, + csa->token); + return 1; + } + return 0; +} + +static int scan_number(struct csa *csa, int across, double *val) +{ if (scan_token(csa, across)) + return 1; + if (strlen(csa->token) == 0) + { xprintf("%s:%d: missing number\n", csa->fname, csa->seqn); + return 1; + } + if (str2num(csa->token, val)) + { xprintf("%s:%d: number '%s' invalid\n", csa->fname, csa->seqn, + csa->token); + return 1; + } + return 0; +} + +TSP *tsp_read_data(const char *fname) +{ struct csa _dsa, *csa = &_dsa; + TSP *tsp = NULL; + csa->fname = fname; + xprintf("Reading TSP data from '%s'...\n", csa->fname); + csa->fp = fopen(csa->fname, "r"); + if (csa->fp == NULL) + { xprintf("Unable to open '%s' - %s\n", csa->fname, + strerror(errno)); + goto fail; + } + tsp = xalloc(1, sizeof(TSP)); + tsp->name = NULL; + tsp->type = TSP_UNDEF; + tsp->comment = NULL; + tsp->dimension = 0; + tsp->edge_weight_type = TSP_UNDEF; + tsp->edge_weight_format = TSP_UNDEF; + tsp->display_data_type = TSP_UNDEF; + tsp->node_x_coord = NULL; + tsp->node_y_coord = NULL; + tsp->dply_x_coord = NULL; + tsp->dply_y_coord = NULL; + tsp->tour = NULL; + tsp->edge_weight = NULL; + csa->seqn = 1; + if (get_char(csa)) + goto fail; +loop: if (scan_keyword(csa)) + goto fail; + if (strcmp(csa->token, "NAME") == 0) + { if (tsp->name != NULL) + { xprintf("%s:%d: NAME entry multiply defined\n", csa->fname, + csa->seqn); + goto fail; + } + if (check_colon(csa)) + goto fail; + if (scan_token(csa, 0)) + goto fail; + if (strlen(csa->token) == 0) + { xprintf("%s:%d: NAME entry incomplete\n", csa->fname, + csa->seqn); + goto fail; + } + tsp->name = xalloc(strlen(csa->token)+1, sizeof(char)); + strcpy(tsp->name, csa->token); + xprintf("NAME: %s\n", tsp->name); + if (check_newline(csa)) + goto fail; + } + else if (strcmp(csa->token, "TYPE") == 0) + { if (tsp->type != TSP_UNDEF) + { xprintf("%s:%d: TYPE entry multiply defined\n", csa->fname, + csa->seqn); + goto fail; + } + if (check_colon(csa)) + goto fail; + if (scan_keyword(csa)) + goto fail; + if (strcmp(csa->token, "TSP") == 0) + tsp->type = TSP_TSP; + else if (strcmp(csa->token, "ATSP") == 0) + tsp->type = TSP_ATSP; + else if (strcmp(csa->token, "TOUR") == 0) + tsp->type = TSP_TOUR; + else + { xprintf("%s:%d: data type '%s' not recognized\n", + csa->fname, csa->seqn, csa->token); + goto fail; + } + xprintf("TYPE: %s\n", csa->token); + if (check_newline(csa)) + goto fail; + } + else if (strcmp(csa->token, "COMMENT") == 0) + { if (check_colon(csa)) + goto fail; + if (scan_comment(csa)) + goto fail; + xprintf("COMMENT: %s\n", csa->token); + if (tsp->comment == NULL) + { tsp->comment = xalloc(strlen(csa->token)+1, sizeof(char)); + strcpy(tsp->comment, csa->token); + } + else + { xprintf("%s:%d: warning: extra COMMENT entry ignored\n", + csa->fname, csa->seqn); + } + if (check_newline(csa)) + goto fail; + } + else if (strcmp(csa->token, "DIMENSION") == 0) + { if (tsp->dimension != 0) + { xprintf("%s:%d: DIMENSION entry multiply defined\n", + csa->fname, csa->seqn); + goto fail; + } + if (check_colon(csa)) + goto fail; + if (scan_integer(csa, 0, &tsp->dimension)) + goto fail; + if (tsp->dimension < 1) + { xprintf("%s:%d: invalid dimension\n", csa->fname, + csa->seqn); + goto fail; + } + xprintf("DIMENSION: %d\n", tsp->dimension); + if (check_newline(csa)) + goto fail; + } + else if (strcmp(csa->token, "EDGE_WEIGHT_TYPE") == 0) + { if (tsp->edge_weight_type != TSP_UNDEF) + { xprintf("%s:%d: EDGE_WEIGHT_TYPE entry multiply defined\n", + csa->fname, csa->seqn); + goto fail; + } + if (check_colon(csa)) + goto fail; + if (scan_keyword(csa)) + goto fail; + if (strcmp(csa->token, "GEO") == 0) + tsp->edge_weight_type = TSP_GEO; + else if (strcmp(csa->token, "EUC_2D") == 0) + tsp->edge_weight_type = TSP_EUC_2D; + else if (strcmp(csa->token, "ATT") == 0) + tsp->edge_weight_type = TSP_ATT; + else if (strcmp(csa->token, "EXPLICIT") == 0) + tsp->edge_weight_type = TSP_EXPLICIT; + else if (strcmp(csa->token, "CEIL_2D") == 0) + tsp->edge_weight_type = TSP_CEIL_2D; + else + { xprintf("%s:%d: edge weight type '%s' not recognized\n", + csa->fname, csa->seqn, csa->token); + goto fail; + } + xprintf("EDGE_WEIGHT_TYPE: %s\n", csa->token); + if (check_newline(csa)) + goto fail; + } + else if (strcmp(csa->token, "EDGE_WEIGHT_FORMAT") == 0) + { if (tsp->edge_weight_format != TSP_UNDEF) + { xprintf("%s:%d: EDGE_WEIGHT_FORMAT entry multiply defined\n" + , csa->fname, csa->seqn); + goto fail; + } + if (check_colon(csa)) + goto fail; + if (scan_keyword(csa)) + goto fail; + if (strcmp(csa->token, "UPPER_ROW") == 0) + tsp->edge_weight_format = TSP_UPPER_ROW; + else if (strcmp(csa->token, "FULL_MATRIX") == 0) + tsp->edge_weight_format = TSP_FULL_MATRIX; + else if (strcmp(csa->token, "FUNCTION") == 0) + tsp->edge_weight_format = TSP_FUNCTION; + else if (strcmp(csa->token, "LOWER_DIAG_ROW") == 0) + tsp->edge_weight_format = TSP_LOWER_DIAG_ROW; + else + { xprintf("%s:%d: edge weight format '%s' not recognized\n", + csa->fname, csa->seqn, csa->token); + goto fail; + } + xprintf("EDGE_WEIGHT_FORMAT: %s\n", csa->token); + if (check_newline(csa)) + goto fail; + } + else if (strcmp(csa->token, "DISPLAY_DATA_TYPE") == 0) + { if (tsp->display_data_type != TSP_UNDEF) + { xprintf("%s:%d: DISPLAY_DATA_TYPE entry multiply defined\n", + csa->fname, csa->seqn); + goto fail; + } + if (check_colon(csa)) + goto fail; + if (scan_keyword(csa)) + goto fail; + if (strcmp(csa->token, "COORD_DISPLAY") == 0) + tsp->display_data_type = TSP_COORD_DISPLAY; + else if (strcmp(csa->token, "TWOD_DISPLAY") == 0) + tsp->display_data_type = TSP_TWOD_DISPLAY; + else + { xprintf("%s:%d: display data type '%s' not recognized\n", + csa->fname, csa->seqn, csa->token); + goto fail; + } + xprintf("DISPLAY_DATA_TYPE: %s\n", csa->token); + if (check_newline(csa)) + goto fail; + } + else if (strcmp(csa->token, "NODE_COORD_SECTION") == 0) + { int n = tsp->dimension, k, node; + if (n == 0) + { xprintf("%s:%d: DIMENSION entry not specified\n", + csa->fname, csa->seqn); + goto fail; + } + if (tsp->node_x_coord != NULL) + { xprintf("%s:%d: NODE_COORD_SECTION multiply specified\n", + csa->fname, csa->seqn); + goto fail; + } + if (check_newline(csa)) + goto fail; + tsp->node_x_coord = xalloc(1+n, sizeof(double)); + tsp->node_y_coord = xalloc(1+n, sizeof(double)); + for (node = 1; node <= n; node++) + tsp->node_x_coord[node] = tsp->node_y_coord[node] = DBL_MAX; + for (k = 1; k <= n; k++) + { if (scan_integer(csa, 0, &node)) + goto fail; + if (!(1 <= node && node <= n)) + { xprintf("%s:%d: invalid node number %d\n", csa->fname, + csa->seqn, node); + goto fail; + } + if (tsp->node_x_coord[node] != DBL_MAX) + { xprintf("%s:%d: node number %d multiply specified\n", + csa->fname, csa->seqn, node); + goto fail; + } + if (scan_number(csa, 0, &tsp->node_x_coord[node])) + goto fail; + if (scan_number(csa, 0, &tsp->node_y_coord[node])) + goto fail; + if (check_newline(csa)) + goto fail; + } + } + else if (strcmp(csa->token, "DISPLAY_DATA_SECTION") == 0) + { int n = tsp->dimension, k, node; + if (n == 0) + { xprintf("%s:%d: DIMENSION entry not specified\n", + csa->fname, csa->seqn); + goto fail; + } + if (tsp->dply_x_coord != NULL) + { xprintf("%s:%d: DISPLAY_DATA_SECTION multiply specified\n", + csa->fname, csa->seqn); + goto fail; + } + if (check_newline(csa)) + goto fail; + tsp->dply_x_coord = xalloc(1+n, sizeof(double)); + tsp->dply_y_coord = xalloc(1+n, sizeof(double)); + for (node = 1; node <= n; node++) + tsp->dply_x_coord[node] = tsp->dply_y_coord[node] = DBL_MAX; + for (k = 1; k <= n; k++) + { if (scan_integer(csa, 0, &node)) + goto fail; + if (!(1 <= node && node <= n)) + { xprintf("%s:%d: invalid node number %d\n", csa->fname, + csa->seqn, node); + goto fail; + } + if (tsp->dply_x_coord[node] != DBL_MAX) + { xprintf("%s:%d: node number %d multiply specified\n", + csa->fname, csa->seqn, node); + goto fail; + } + if (scan_number(csa, 0, &tsp->dply_x_coord[node])) + goto fail; + if (scan_number(csa, 0, &tsp->dply_y_coord[node])) + goto fail; + if (check_newline(csa)) + goto fail; + } + } + else if (strcmp(csa->token, "TOUR_SECTION") == 0) + { int n = tsp->dimension, k, node; + if (n == 0) + { xprintf("%s:%d: DIMENSION entry not specified\n", + csa->fname, csa->seqn); + goto fail; + } + if (tsp->tour != NULL) + { xprintf("%s:%d: TOUR_SECTION multiply specified\n", + csa->fname, csa->seqn); + goto fail; + } + if (check_newline(csa)) + goto fail; + tsp->tour = xalloc(1+n, sizeof(int)); + for (k = 1; k <= n; k++) + { if (scan_integer(csa, 1, &node)) + goto fail; + if (!(1 <= node && node <= n)) + { xprintf("%s:%d: invalid node number %d\n", csa->fname, + csa->seqn, node); + goto fail; + } + tsp->tour[k] = node; + } + if (scan_integer(csa, 1, &node)) + goto fail; + if (node != -1) + { xprintf("%s:%d: extra node(s) detected\n", csa->fname, + csa->seqn); + goto fail; + } + if (check_newline(csa)) + goto fail; + } + else if (strcmp(csa->token, "EDGE_WEIGHT_SECTION") == 0) + { int n = tsp->dimension, i, j, temp; + if (n == 0) + { xprintf("%s:%d: DIMENSION entry not specified\n", + csa->fname, csa->seqn); + goto fail; + } + if (tsp->edge_weight_format == TSP_UNDEF) + { xprintf("%s:%d: EDGE_WEIGHT_FORMAT entry not specified\n", + csa->fname, csa->seqn); + goto fail; + } + if (tsp->edge_weight != NULL) + { xprintf("%s:%d: EDGE_WEIGHT_SECTION multiply specified\n", + csa->fname, csa->seqn); + goto fail; + } + if (check_newline(csa)) + goto fail; + tsp->edge_weight = xalloc(1+n*n, sizeof(int)); + switch (tsp->edge_weight_format) + { case TSP_FULL_MATRIX: + for (i = 1; i <= n; i++) + { for (j = 1; j <= n; j++) + { if (scan_integer(csa, 1, &temp)) + goto fail; + tsp->edge_weight[(i - 1) * n + j] = temp; + } + } + break; + case TSP_UPPER_ROW: + for (i = 1; i <= n; i++) + { tsp->edge_weight[(i - 1) * n + i] = 0; + for (j = i + 1; j <= n; j++) + { if (scan_integer(csa, 1, &temp)) + goto fail; + tsp->edge_weight[(i - 1) * n + j] = temp; + tsp->edge_weight[(j - 1) * n + i] = temp; + } + } + break; + case TSP_LOWER_DIAG_ROW: + for (i = 1; i <= n; i++) + { for (j = 1; j <= i; j++) + { if (scan_integer(csa, 1, &temp)) + goto fail; + tsp->edge_weight[(i - 1) * n + j] = temp; + tsp->edge_weight[(j - 1) * n + i] = temp; + } + } + break; + default: + goto fail; + } + if (check_newline(csa)) + goto fail; + } + else if (strcmp(csa->token, "EOF") == 0) + { if (check_newline(csa)) + goto fail; + goto done; + } + else + { xprintf("%s:%d: keyword '%s' not recognized\n", csa->fname, + csa->seqn, csa->token); + goto fail; + } + goto loop; +done: xprintf("%d lines were read\n", csa->seqn-1); + fclose(csa->fp); + return tsp; +fail: if (tsp != NULL) + { if (tsp->name != NULL) + xfree(tsp->name); + if (tsp->comment != NULL) + xfree(tsp->comment); + if (tsp->node_x_coord != NULL) + xfree(tsp->node_x_coord); + if (tsp->node_y_coord != NULL) + xfree(tsp->node_y_coord); + if (tsp->dply_x_coord != NULL) + xfree(tsp->dply_x_coord); + if (tsp->dply_y_coord != NULL) + xfree(tsp->dply_y_coord); + if (tsp->tour != NULL) + xfree(tsp->tour); + if (tsp->edge_weight != NULL) + xfree(tsp->edge_weight); + xfree(tsp); + } + if (csa->fp != NULL) + fclose(csa->fp); + return NULL; +} + +/*********************************************************************** +* NAME +* +* tsp_distance - compute distance between two nodes +* +* SYNOPSIS +* +* #include "tsplib.h" +* int tsp_distance(TSP *tsp, int i, int j); +* +* DESCRIPTION +* +* The routine tsp_distance computes the distance between i-th and j-th +* nodes for the TSP instance, which tsp points to. +* +* RETURNS +* +* The routine tsp_distance returns the computed distance. */ + +#define nint(x) ((int)((x) + 0.5)) + +static double rad(double x) +{ /* convert input coordinate to longitude/latitude, in radians */ + double pi = 3.141592, deg, min; + deg = (int)x; + min = x - deg; + return pi * (deg + 5.0 * min / 3.0) / 180.0; +} + +int tsp_distance(const TSP *tsp, int i, int j) +{ int n = tsp->dimension, dij; + if (!(tsp->type == TSP_TSP || tsp->type == TSP_ATSP)) + xerror("tsp_distance: invalid TSP instance\n"); + if (!(1 <= i && i <= n && 1 <= j && j <= n)) + xerror("tsp_distance: node number out of range\n"); + switch (tsp->edge_weight_type) + { case TSP_UNDEF: + xerror("tsp_distance: edge weight type not specified\n"); + case TSP_EXPLICIT: + if (tsp->edge_weight == NULL) + xerror("tsp_distance: edge weights not specified\n"); + dij = tsp->edge_weight[(i - 1) * n + j]; + break; + case TSP_EUC_2D: + if (tsp->node_x_coord == NULL || tsp->node_y_coord == NULL) + xerror("tsp_distance: node coordinates not specified\n"); + { double xd, yd; + xd = tsp->node_x_coord[i] - tsp->node_x_coord[j]; + yd = tsp->node_y_coord[i] - tsp->node_y_coord[j]; + dij = nint(sqrt(xd * xd + yd * yd)); + } + break; + case TSP_CEIL_2D: + if (tsp->node_x_coord == NULL || tsp->node_y_coord == NULL) + xerror("tsp_distance: node coordinates not specified\n"); + { double xd, yd; + xd = tsp->node_x_coord[i] - tsp->node_x_coord[j]; + yd = tsp->node_y_coord[i] - tsp->node_y_coord[j]; + dij = (int)ceil(sqrt(xd * xd + yd * yd)); + } + break; + case TSP_GEO: + if (tsp->node_x_coord == NULL || tsp->node_y_coord == NULL) + xerror("tsp_distance: node coordinates not specified\n"); + { double rrr = 6378.388; + double latitude_i = rad(tsp->node_x_coord[i]); + double latitude_j = rad(tsp->node_x_coord[j]); + double longitude_i = rad(tsp->node_y_coord[i]); + double longitude_j = rad(tsp->node_y_coord[j]); + double q1 = cos(longitude_i - longitude_j); + double q2 = cos(latitude_i - latitude_j); + double q3 = cos(latitude_i + latitude_j); + dij = (int)(rrr * acos(0.5 * ((1.0 + q1) * q2 - + (1.0 - q1) *q3)) + 1.0); + } + break; + case TSP_ATT: + if (tsp->node_x_coord == NULL || tsp->node_y_coord == NULL) + xerror("tsp_distance: node coordinates not specified\n"); + { int tij; + double xd, yd, rij; + xd = tsp->node_x_coord[i] - tsp->node_x_coord[j]; + yd = tsp->node_y_coord[i] - tsp->node_y_coord[j]; + rij = sqrt((xd * xd + yd * yd) / 10.0); + tij = nint(rij); + if (tij < rij) dij = tij + 1; else dij = tij; + } + break; + default: + xassert(tsp->edge_weight_type != tsp->edge_weight_type); + } + return dij; +} + +/*********************************************************************** +* NAME +* +* tsp_free_data - free TSP instance data +* +* SYNOPSIS +* +* #include "tsplib.h" +* void tsp_free_data(TSP *tsp); +* +* DESCRIPTION +* +* The routine tsp_free_data frees all the memory allocated to the TSP +* instance data block, which the parameter tsp points to. */ + +void tsp_free_data(TSP *tsp) +{ if (tsp->name != NULL) + xfree(tsp->name); + if (tsp->comment != NULL) + xfree(tsp->comment); + if (tsp->node_x_coord != NULL) + xfree(tsp->node_x_coord); + if (tsp->node_y_coord != NULL) + xfree(tsp->node_y_coord); + if (tsp->dply_x_coord != NULL) + xfree(tsp->dply_x_coord); + if (tsp->dply_y_coord != NULL) + xfree(tsp->dply_y_coord); + if (tsp->tour != NULL) + xfree(tsp->tour); + if (tsp->edge_weight != NULL) + xfree(tsp->edge_weight); + xfree(tsp); + return; +} + +/* eof */ diff --git a/glpk-5.0/examples/tsp/tsplib.h b/glpk-5.0/examples/tsp/tsplib.h new file mode 100644 index 0000000..19936ad --- /dev/null +++ b/glpk-5.0/examples/tsp/tsplib.h @@ -0,0 +1,80 @@ +/* tsplib.h */ + +/* Written by Andrew Makhorin <mao@gnu.org>, October 2015. */ + +#ifndef TSPLIB_H +#define TSPLIB_H + +typedef struct TSP TSP; + +struct TSP +{ /* TSP (or related problem) instance in the format described in + * [G.Reinelt, TSPLIB 95] */ + /*--------------------------------------------------------------*/ + /* specification part */ + char *name; + /* identifies the data file */ + int type; + /* specifies the type of data: */ +#define TSP_UNDEF 0 /* undefined */ +#define TSP_TSP 1 /* symmetric TSP */ +#define TSP_ATSP 2 /* asymmetric TSP */ +#define TSP_TOUR 3 /* collection of tours */ + char *comment; + /* additional comments (usually the name of the contributor or + * creator of the problem instance is given here) */ + int dimension; + /* for a TSP or ATSP, the dimension is the number of its nodes + * for a TOUR it is the dimension of the corresponding problem */ + int edge_weight_type; + /* specifies how the edge weights (or distances) are given: */ +#define TSP_UNDEF 0 /* undefined */ +#define TSP_EXPLICIT 1 /* listed explicitly */ +#define TSP_EUC_2D 2 /* Eucl. distances in 2-D */ +#define TSP_CEIL_2D 3 /* Eucl. distances in 2-D rounded up */ +#define TSP_GEO 4 /* geographical distances */ +#define TSP_ATT 5 /* special distance function */ + int edge_weight_format; + /* describes the format of the edge weights if they are given + * explicitly: */ +#define TSP_UNDEF 0 /* undefined */ +#define TSP_FUNCTION 1 /* given by a function */ +#define TSP_FULL_MATRIX 2 /* given by a full matrix */ +#define TSP_UPPER_ROW 3 /* upper triangular matrix (row-wise + * without diagonal entries) */ +#define TSP_LOWER_DIAG_ROW 4 /* lower triangular matrix (row-wise + * including diagonal entries) */ + int display_data_type; + /* specifies how a graphical display of the nodes can be + * obtained: */ +#define TSP_UNDEF 0 /* undefined */ +#define TSP_COORD_DISPLAY 1 /* display is generated from the node + * coordinates */ +#define TSP_TWOD_DISPLAY 2 /* explicit coordinates in 2-D are + * given */ + /*--------------------------------------------------------------*/ + /* data part */ + /* NODE_COORD_SECTION: */ + double *node_x_coord; /* double node_x_coord[1+dimension]; */ + double *node_y_coord; /* double node_y_coord[1+dimension]; */ + /* DISPLAY_DATA_SECTION: */ + double *dply_x_coord; /* double dply_x_coord[1+dimension]; */ + double *dply_y_coord; /* double dply_y_coord[1+dimension]; */ + /* TOUR_SECTION: */ + int *tour; /* int tour[1+dimension]; */ + /* EDGE_WEIGHT_SECTION: */ + int *edge_weight; /* int edge_weight[1+dimension*dimension]; */ +}; + +TSP *tsp_read_data(const char *fname); +/* read TSP instance data */ + +int tsp_distance(const TSP *tsp, int i, int j); +/* compute distance between two nodes */ + +void tsp_free_data(TSP *tsp); +/* free TSP instance data */ + +#endif + +/* eof */ diff --git a/glpk-5.0/examples/tsp/ulysses16.tsp b/glpk-5.0/examples/tsp/ulysses16.tsp new file mode 100644 index 0000000..4ce24a6 --- /dev/null +++ b/glpk-5.0/examples/tsp/ulysses16.tsp @@ -0,0 +1,24 @@ +NAME: ulysses16.tsp +TYPE: TSP +COMMENT: Odyssey of Ulysses (Groetschel/Padberg) +DIMENSION: 16 +EDGE_WEIGHT_TYPE: GEO +DISPLAY_DATA_TYPE: COORD_DISPLAY +NODE_COORD_SECTION + 1 38.24 20.42 + 2 39.57 26.15 + 3 40.56 25.32 + 4 36.26 23.12 + 5 33.48 10.54 + 6 37.56 12.19 + 7 38.42 13.11 + 8 37.52 20.44 + 9 41.23 9.10 + 10 41.17 13.05 + 11 36.08 -5.21 + 12 38.47 15.13 + 13 38.15 15.35 + 14 37.51 15.17 + 15 35.49 14.32 + 16 39.36 19.56 + EOF diff --git a/glpk-5.0/examples/tsp/ulysses22.tsp b/glpk-5.0/examples/tsp/ulysses22.tsp new file mode 100644 index 0000000..9e95fb8 --- /dev/null +++ b/glpk-5.0/examples/tsp/ulysses22.tsp @@ -0,0 +1,30 @@ +NAME: ulysses22.tsp +TYPE: TSP +COMMENT: Odyssey of Ulysses (Groetschel/Padberg) +DIMENSION: 22 +EDGE_WEIGHT_TYPE: GEO +DISPLAY_DATA_TYPE: COORD_DISPLAY +NODE_COORD_SECTION + 1 38.24 20.42 + 2 39.57 26.15 + 3 40.56 25.32 + 4 36.26 23.12 + 5 33.48 10.54 + 6 37.56 12.19 + 7 38.42 13.11 + 8 37.52 20.44 + 9 41.23 9.10 + 10 41.17 13.05 + 11 36.08 -5.21 + 12 38.47 15.13 + 13 38.15 15.35 + 14 37.51 15.17 + 15 35.49 14.32 + 16 39.36 19.56 + 17 38.09 24.36 + 18 36.09 23.00 + 19 40.44 13.57 + 20 40.33 14.15 + 21 40.37 14.23 + 22 37.57 22.56 +EOF |