SHARE
TWEET

Untitled

a guest May 20th, 2019 59 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. #=
  2. main:
  3. - Julia version:
  4. - Author: quangio
  5. - Date: 2019-05-06
  6. =#
  7.  
  8.  
  9.  
  10. using JuMP
  11. using SCIP
  12. # include("input.jl")
  13.  
  14. model = Model(with_optimizer(SCIP.Optimizer))
  15.  
  16.  
  17. dist_mat_gr17 = [
  18.   0 633 257  91 412 150  80 134 259 505 353 324  70 211 268 246 121
  19. 633   0 390 661 227 488 572 530 555 289 282 638 567 466 420 745 518
  20. 257 390   0 228 169 112 196 154 372 262 110 437 191  74  53 472 142
  21.  91 661 228   0 383 120  77 105 175 476 324 240  27 182 239 237  84
  22. 412 227 169 383   0 267 351 309 338 196  61 421 346 243 199 528 297
  23. 150 488 112 120 267   0  63  34 264 360 208 329  83 105 123 364  35
  24.  80 572 196  77 351  63   0  29 232 444 292 297  47 150 207 332  29
  25. 134 530 154 105 309  34  29   0 249 402 250 314  68 108 165 349  36
  26. 259 555 372 175 338 264 232 249   0 495 352  95 189 326 383 202 236
  27. 505 289 262 476 196 360 444 402 495   0 154 578 439 336 240 685 390
  28. 353 282 110 324  61 208 292 250 352 154   0 435 287 184 140 542 238
  29. 324 638 437 240 421 329 297 314  95 578 435   0 254 391 448 157 301
  30.  70 567 191  27 346  83  47  68 189 439 287 254   0 145 202 289  55
  31. 211 466  74 182 243 105 150 108 326 336 184 391 145   0  57 426  96
  32. 268 420  53 239 199 123 207 165 383 240 140 448 202  57   0 483 153
  33. 246 745 472 237 528 364 332 349 202 685 542 157 289 426 483   0 336
  34. 121 518 142  84 297  35  29  36 236 390 238 301  55  96 153 336   0] # shoud return 2085
  35. dist_mat_5 = [
  36. 0.0  3.0  4.0  2.0  7.0
  37. 3.0  0.0  4.0  6.0  3.0
  38. 4.0  4.0  0.0  5.0  8.0
  39. 2.0  6.0  5.0  0.0  6.0
  40. 7.0  3.0  8.0  6.0  0.0
  41. ] # shoud return 19.0
  42.  
  43. n = 17 # number of vertices
  44. c = zeros(n)   # cost assigned to each vertex (weight needed to be picked up)
  45. d = dist_mat_gr17 # TestSet.dist_mat_48 # distance matrix of the graph
  46.  
  47. w0 = 1
  48. no_salesman = 1
  49. cap = 100
  50.  
  51. range = collect(1:n)
  52.  
  53. for i = 1:n
  54.     d[i, i] = 0
  55. end
  56.  
  57.  
  58. @variable(model, y[1:n, 1:n])
  59. @variable(model, x[1:n, 1:n], Bin)
  60.  
  61. @objective(model, Min, sum(y[i, j] * d[i, j] for i = 1:n, j = range[1:n .!= i]))
  62.  
  63. @constraint(model, sum(x[1, i] for i = 2:n) == no_salesman)
  64. @constraint(model, sum(x[i, 1] for i = 2:n) == no_salesman)
  65.  
  66. for i = 2:n
  67.     @constraint(model, y[1, i] == w0 * x[1, i])
  68.     @constraint(model, sum(x[i, j] for j = range[1:n .!= i]) == 1)
  69.     @constraint(model, sum(x[j, i] for j = range[1:n .!= i]) == 1)
  70.     @constraint(model, c[i] == sum(y[i, j] - y[j, i] for j = range[1:n .!= i]))
  71. end
  72.  
  73. for i = 1:n
  74.     for j = 1:n
  75.         if i == j continue end
  76.         @constraint(model, y[i, j] >= (w0 + c[i]) * x[i, j])
  77.         @constraint(model, y[i, j] <= (w0 + cap - c[j]) * x[i, j])
  78.     end
  79. end
  80.  
  81.  
  82. function solved(m)
  83.     x_val = JuMP.value.(x)
  84.     cycle_idx = [1] # we can adapt this to support mTSP, just keep this for simplicity for now
  85.     while true
  86.         v, idx = findmax(x_val[cycle_idx[length(cycle_idx)], :])
  87.         if idx == 1
  88.             break
  89.         else
  90.             push!(cycle_idx, idx)
  91.         end
  92.     end
  93.     if length(cycle_idx) < n
  94.         @constraint(m, sum(x[i, j] for i = cycle_idx, j = cycle_idx[cycle_idx .!= i]) <= length(cycle_idx) - 1)
  95.         return false
  96.     end
  97.     return true
  98. end
  99.  
  100. if no_salesman == 0
  101.     JuMP.optimize!(model)
  102.     while !solved(model)
  103.         JuMP.optimize!(model)
  104.     end
  105. end
  106.  
  107. if no_salesman >= 1
  108.     @variable(model, 0 <= u[1:n] <= n - 1)
  109.     for i = 2:n
  110.         for j = collect(2:n)[2:n .!= i]
  111.             @constraint(model, u[i] - u[j] + n * x[i, j] <= n - 1)
  112.         end
  113.     end
  114.     @time JuMP.optimize!(model)
  115. end
  116.  
  117. println(string("Objective value:", JuMP.objective_value(model)))
  118.  
  119. println("Arcs: ")
  120.  
  121. for i = range
  122.     for j = range
  123.         if JuMP.value(x[i, j]) > 0.0
  124.             println(string(i, "->", j, "[label=", d[i, j], "]")) # you can use graphviz to visualize this
  125.         end
  126.     end
  127. end
  128.  
  129. open("model2.lp", "w") do f
  130.     print(f, model)
  131. end
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
 
Top