Lesson Route Optimization Problem - Artificial Intelligence - ثالث ثانوي

Lesson 3 Route Optimization Problem Link to digital lesson www.ien.edu.sa Mathematical Programming In Optimization Problems Mathematical Programming Mathematical programming is a technique to solve optimization problems by formulating them as mathematical models. The previous two lessons demonstrated how heuristic algorithms could be used to solve different types of optimization problems. While heuristics can be very fast and often produce good solutions, they do not always guarantee the optimal solution and may not be suitable for all types of problems. In this lesson you will focus on a different optimization approach: mathematical programming. Mathematical programming can solve many optimization problems, such as resource allocation, production planning, logistics, and scheduling. The technique has the advantage of providing a guaranteed optimal solution and can handle complex problems with multiple constraints. A mathematical programming solution starts with formulating the given optimization problem as a mathematical model using variables. These variables represent the values that have to be optimized. They are used to define the objective function and constraints, which together describe the problem and enable the use of mathematical programming algorithms. Mathematical programming utilizes decision variables which can be controlled and tuned by the decision-maker to find the solution, or state variables which the decision-maker has no control over and are imposed by the external environment. State variables cannot be tuned. The following lists provide examples of decision and state variables for some popular optimization problems: Table 5.2: Examples of decision and state variables Production Planning Resource Transportation Job Scheduling Personnel Rostering Decision Variables The quantity of each product that has to be produced. The number of goods to be transported from one location to another. The order and time duration of each job to be performed. The assignment and scheduling of workers to different tasks at different times. State Variables The availability of raw materials, production machines' capacity, and production labor availability. The distance between the locations that must be visited and the capacity of the vehicles. The availability of workers and machines, the deadlines, and the importance weights of each job. The skills, preferences, and availability of each worker. The skills required to complete each task. وزارة التعليم Ministry of Education 2024-1446 283

The objective function is formulated as a mathematical expression to be optimized (maximized or minimized) based on the relevant variables. This function represents the goal of the optimization problem, such as maximizing profit or minimizing costs. It is usually defined in terms of the decision variables and sometimes the state variables. Similarly, the constraints can be formulated using variables and mathematical inequalities. There are several types of mathematical programming, including Linear Programming (LP), Quadratic Programming (QP), and Mixed Integer Programming (MIP). This lesson focuses on MIP, which is used for problems where the decision variables are restricted to integers, such as scheduling or routing problems. The Knapsack Problem A simple example of using MIP to formulate the objective function and constraints is the 0/1 Knapsack Problem. The problem is defined as follows: you are given a knapsack with a maximum weight capacity of C and a set of items I. Each item i has two state variables: a weight w, and a value V.. The requirement is to fill the knapsack with the maximum possible value within the knapsack's capacity. A decision variable is also used to keep track of the combinations of items to be packed in a knapsack, where x = 1 if the item i is selected to be added to the knapsack and x, = 0 otherwise. The goal is to select a subset of the items from I such that: xi • Constraint: The sum of the weights of the selected items is not greater than the maximum capacity C. • Objective function: The sum of the values of the selected items is as high as possible (maximized). V=20 w=5 W₁ =8 V₁=15 V=10 w=10 --> V₂-23 W₂-19 --> W₁ =11 V₁ =7 Figure 5.6: The Knapsack problem W=2 V=7 A knapsack instance is illustrated in figure 5.6 with six items having specific weights and values and a maximum knapsack capacity of 40. The following code installs and uses the open-source Python library mip (mixed integer programming) to solve an instance of the 0/1 Knapsack problem and imports the necessary modules. !pip install mip #install the mip library # imports useful tools from the mip library from mip import Model, xsum, maximize, BINARY values = [20, 10, 23, 15, 7, 7] # values of available items Pilill weights = [5, 10, 19, 8, 11, 2] #weights of available items قرارة التعليم Mistry of Education 2024-1446

C = 40 # knapsack capacity I range(len(values)) # creates an index for each item: 0,1,2,3,... = solver Model("knapsack") # creates a knapsack solver = solver.verbose = 0 # setting this to 1 will print more information on the progress of the solver x = [] # represents the binary decision variables for each item # for each items creates and appends a binary decision variable for i in I: x.append(solver.add_var(var_type = BINARY)) # creates the objective function solver.objective = maximize (xsum(values[i] * x[i] for i in I)) # adds the capacity constraint to the solver solver + xsum (weights[i] += #solves the problem solver.optimize() * x[i] for i in I) <= C <Optimization Status. OPTIMAL: 0> The code creates a list x to store the binary decision variables for the items. The mip library provides: • the add_var(var_type=BINARY) tool for creating binary variables and adding them to the solver. ⚫ the maximize() tool for optimization problems that need to maximize an objective function, and the minimize() tool used for optimization problems that need to minimize an objective. the xsum() tool for creating mathematical expressions that include sums. In the above example, the tool is used to compute the total weight of the items in a solution and create the capacity constraint. ⚫ the optimize() tool for finding a solution that optimizes the objective function while respecting the constraints. The tool uses MIP to efficiently consider different combinations of values for the decision variables and find the one that optimizes the objective. ⚫ the + operator to add additional constraints to an existing solver. In the implementation below, the list x holds one binary variable for each item. After the solution has been computed, each variable will be equal to 1 if the item was included in the solution and equal to 0 otherwise. The mip library uses the x[i].x syntax to return the binary value for the item with index i. The solver computes the decision variable x, then finds the total value and weight of the selected items by iterating over the decision variable x, accumulating the weights and values for each selected item i, based on x[i], and displaying them as shown in the following code. • total_weight total weight = 0 # stores the total weight of the items in the solution total_value = 0 # stores the total value of the items in the solution وزارة التعليم Ministry of Education 2024-1446 285

for i in I: # for each item if x[i].x == 1: # if the item was selected print('item', i, 'was selected') # updates the total weight and value of the solution total_weight += weights[i] total_value += values[i] print('total weight', total_weight) print('total value', total_value) item 0 was selected item 2 was selected item 3 was selected item 5 was selected total weight 34 total value 65 Traveling Salesman Problem Another problem that can be solved via MIP is the Traveling Salesman Problem (TSP). It is a classic problem that seeks to determine the shortest possible route a salesman must take to visit a set of cities exactly once and then return to his starting point, without visiting any city twice. The figure 5.7 visualizes an instance of this problem. Each circle (node) represents a city or location that has to be visited. There is an edge connecting two locations if it is possible to travel between them. The number on the edge represents the cost (distance) between the two locations. In this example, the locations have been numbered according to their order in the optimal solution to the problem. The total cost of the route 1→2→4→3→1 is 10 + 25 +30 + 15 = 80, which is the shortest possible route that visits every city exactly once and returns to the starting point. TSP has practical applications in logistics, transportation, supply chain management, and telecommunications. It belongs to a broader family of routing problems that also includes other famous problems, which are presented below: TSP graph instances are fully connected; there is an edge connecting every pair of nodes. 1 20 20 10 15 4 25 30 2 35 3 Figure 5.7: Instance of Traveling Salesman problem • The Vehicle Routing Problem involves finding the optimal routes for a fleet of vehicles to deliver goods or services to a set of customers while minimizing the total distance traveled. Applications include logistics, delivery services, and garbage collection. • The Pickup and Delivery Problem involves finding the optimal routes for vehicles to pick up and deliver goods or people to different locations. Applications include taxi services, emergency medical services, and shuttle services. The Train Timetabling Problem involves finding the optimal train schedules in a railway network while minimizing delay percentage and ensuring efficient use of resources. Applications include .ailway transportation and schedulingرة التعليم Ministry of Education 286 2024-1446

The following code can be used to create an instance of the TSP. The function accepts the number of locations to be visited and the distance range (minimum and maximum distance) between two locations. It then returns: • a distance matrix that includes the distance between every possible pair of locations. • a set of numeric location ids (one for each location). ⚫ the location that serves as the start and end of the route. This is referred to as the 'startstop' location. import random import numpy from itertools import combinations def create_problem_instance(num_locations, distance_range): #initializes the distance matrix to be full of zeros dist_matrix = numpy.zeros((num_locations, num_locations)) # creates location ids: 0,1,2,3,4,... location_ids = set(range(num_locations)) # creates all possible location pairs location_pairs = combinations (location_ids, 2) for i,j in location_pairs: # for each pair distance = random.randint(*distance_range) # samples a distance within range # the distance from i to j is the same as the distance from j to i dist_matrix[j,i] = distance dist matrix[i, j] = distance # returns the distance matrix, location ids and the startstop vertix return dist_matrix, location_ids, random.randint(0, num_locations 1) - The following code used the above function to create a TSP instance with 8 locations and pairwise distances between 5 and 20: dist_matrix, location_ids, startstop = create_problem_instance(8, (5, 20)) print(dist matrix) print(startstop) وزارة التعليم Ministry of Education 2024-1446 [[O 19. 17. 15. 18. 17. 7. 15.] [19. 0. 15. 18. 11. 6. 20. 5.] [17. 15. 0. 17. 15. 7. 5. 11.] [15. 18. 17. Θ 19. 7. 7. 16.] [18. 11. 15. 19. Θ 17. 20. 17. ] [17. 6. 7. 7. 17. 0. 15. 14.] 7. 20. 15. 0. 14.] [ 7. 20. 5. [15. 5. 11. 16. 17. 14. 14. 0.11 Notice that the diagonal represents the distances from the nodes to themselves (dist_ matrix[i,i]), and thus is zero. 287

Creating a Brute-Force Solver for the Traveling Salesman Problem The following function uses brute force to exhaustively enumerate all possible routes (permutations) and return the shortest one. It accepts the distance matrix and the startstop location returned by the create_problem_instance() function. Note that a solution to a TSP instance is a permutation of cities, starting and ending at the startstop city. from itertools import permutations def brute force_solver (dist_matrix, location_ids, startstop): step # excludes the starstop location location_ids = location_ids - {startstop} # generate all possible routes (location permutations) all_routes = permutations (location_ids) best_distance = float('inf') # initializes to the highest possible number best_route = None # best route so far, initialized to None for route in all_routes: # for each route = distance curr_loc = 0 # total distance in this route startstop # current location for next_loc in route: distance += dist_matrix[curr_loc, next_loc] # adds the distance for this curr_loc = next_loc #goes to the next location distance += dist_matrix[curr_loc, startstop] # goes to the starstop location if distance < best_distance: # if this route has shorter distance than the best route best_distance = distance best_route = route # adds the startstop location at the beginning and end of the best route and returns return [startstop] + list(best_route) [startstop], best_distance + The brute force solver uses the permutations() tool to create all possible routes. Note that the startstop location is excluded from the permutations, as it must always appear at the start and end of each route. For example, if we have 4 locations 0,1,2,3 and 0 is the startstop location, then the list of possible permutations would be: for route in permutations({1,2,3}): print(route) (1, 2, 3) (1, 3, 2) (2, 1, 3) (2, 3, 1) • (3, 1, 2) (3, 2, 1) وزارة التعليم Magtry of Education 2024-1446

The brute force solver computes the total distance of each route and finally returns the one with the shortest distance. The following code applies the solver to the TSP instance generated above. brute force_solver (dist_matrix, location_ids, startstop) ([3, 5, 2, 7, 1, 4, 0, 6, 3], 73.0) Similar to the brute-force solvers that were described in the previous lessons, this solver is only applicable to small TSP instances. This is because the number of possible routes is a number that grows exponentially as N gets larger and is equal to (N-1)!. For example, for N=15, the number of possible routes is equal to 14! = 87,178,291,200. Using MIP to Solve the Traveling Salesman Problem To use MIP to solve the TSP, a mathematical formulation that covers both the objective function and the constraints of the TSP needs to be created. ij The formulation requires a binary decision variable x.. for every possible transition i→j from a location i to another location j. For a problem with N locations, the number of possible transitions is equal to Nx(N-1). If x is equal to 1, the solution includes a transition from location i to location j. Otherwise, if xij is equal to 0, then this transition is not included in the solution. ij Accessing elements in a 2-dimensional numpy array can be easily done via the [i,j] syntax. For example: arr numpy.full((4,4), 0) #creates a 4x4 array full of zeros = print(arr) arr[0, 0] arr[3, 3] = = 1 1 [[0000] Го 0 0 01 [ 0 0 0 0] [0 0 0 0]] [[1 0 0 0] print() print(arr) [0 0 0 0] [ 0 0 0 0 ] [0 0 0 1]] The code also uses the product() tool from 'itertools' to compute all possible location transitions. For example: ids = {0, 1, 2} for i, j in list (product(ids, ids)): print(i, j) 0 0 0 1 0 2 1 0 1 1 1 2 2 0 2 1 2 2 وزارة التعليم Ministry of Education 2024-1446 289

The following code uses the Python mip library to create an MIP solver. It then adds one binary decision variable for every possible transition in the TSP instance generated above. from itertools import product # used to generate all possible transitions from mip import BINARY from mip import Model, INTEGER solver = Model() # creates a solver solver.verbose = 0 # setting this to 1 will print info on the progress of the solver # 'product' creates every transition from every location to every other location transitions = list(product (location_ids, location_ids)) N = len(location_ids) # number of locations # creates a square numpy array full of 'None' values x numpy.full((N, N), None) = # adds binary variables indicating if transition (i->j) is included in the route for i, j in transitions: x[i, j] = solver.add_var(var_type = BINARY) The above code uses the numpy.full() tool to create NxN numpy array for storing the binary x variables. After adding the x decision variables, the following code can be used to formulate the objective function for the TSP. The function iterates over every possible transition i→j and multiplies its distance dist_matrix[i,j] with its decision variable x[i,j]. If the transition is included in the solution, then x[i,j]=1 and dist_matrix[i,j] will be considered. Otherwise, dist_matrix[i,j] will be multiplied by 0 and will be ignored. # the minimize tool is used then the objective function has to be minimized from mip import xsum, minimize # objective function: minimizes the distance solver.objective = minimize (xsum (dist_matrix[i, j]*x[i][j] for i,j in transitions)) minimize(xsum(dist_matrix[i,j]*x[i][j] The next step is to ensure that the solver only reports solutions that visit each location, except for the startstop, exactly once, as the TSP requires. Visiting each location once means that a valid route can only: • arrive at each location exactly once. • depart from each location exactly once. These arrive/depart constraints can be easily added as follows: # for each location id for i in location_ids: solver += xsum(x[i,j] for j in location_ids - {i}) ©solver += xsum(x[j,i] for j in location_ids - {i}) == 1 # exactly 1 arrival == 1 # exactly 1 departure وزارة التعليم Ministry of Education 290 2024-1446

The complete formulation of the TSP includes one more type of constraint to ensure the computation of connected routes. Consider the TSP instance in figure 5.8. Location O is assumed to be the startstop location. In this instance, the shortest possible route is 03→4→1→2, with a total travel distance of 24. However, without a connectivity constraint, a solution with two unconnected routes 0 3→4→0 and 1-2-1 is also valid. This 2-route solution satisfies the arrive/depart constraints defined in the code above, as every location is entered and exited exactly once. However, this is not an acceptable TSP solution. A solution with a single connected route can be enforced by adding the decision variable y for every location i. These variables will capture the order in which each location is visited in the solution. 0 5 2 50 1 6 4 7 5 3 3 Figure 5.8: TSP instance # adds a decision variable for each location y = [solver.add_var(var_type = INTEGER) for i in location_ids] For example, if the solution is 0⇒3→4→1→2→0, then the y values would be y3=0, y4=1, y1=2, y2=3. Location O is the startstop location, so its y value is not considered. These new decision variables can be used to ensure connectivity by adding a new constraint for each transition that does not include the startstop: # adds a connectivity constraint for every transition that does not include the startstop for (i, j) in product (location_ids - {startstop}, location_ids - {startstop}): # ignores transitions from a location to itself if i != j: solver += y[j] - y[i] >= (N+1) * x[i, j] - N If a transition i→j has x=1 and is included in the solution, then the above inequality becomes y[j] >= y[i] + 1. This states that locations that are visited later are required to have higher y values. Combined with the arrive/depart constraints, a route that does not include the startstop is valid only if: • If it starts and ends with the same location, to ensure that all locations have exactly one arrival and one departure. ⚫ It assigns higher y values to locations that are visited later since y[j] has to be greater than y[i] for all transitions included in the route. This also avoids adding the same edge from a different direction (e.g., i→j and j→i). However, if a location serves as both the start and end of a route, then it needs to have a y value that is both higher and lower than those of all the others in the route. Given that this is impossible, adding the connectivity constraint eliminates any solutions with routes that do not include the startstop. وزارة التعليم Ministry of Education 2024-1446 291

As an example, consider the 1→2→1 route in the 2-route solution of the TSP instance shown in the figure above. The connectivity constraint requires that y₂ ≥ y₁ + 1 and y₁ ≥ y2 + 1. This is not possible, so the solution would be eliminated. 2 1 1 4 y3 1 4 In constrast, the correct solution 0→3→4→1→2→0 requires that y₁ ≥ 13 + 1, y₁ ≥ y + 1, and Y₂ ≥ Y₁ + 1. These can be satisfied by setting y₁ = 0, y₁ =1, y₁ =2, y2=3. Connectivity constraints do not apply to transitions that include the startstop node. 2 1 The following function puts everything together to create a complete MIP solver for the TSP: from itertools import product from mip import BINARY, INTEGER from mip import Model from mip import xsum, minimize def MIP_solver(dist_matrix, location_ids, startstop): solver Model() # creates a solver = solver.verbose = 0 # setting this to 1 will print info on the progress of the solver # creates every transition from every location to every other location transitions = list(product(location_ids, location_ids)) N = len(location_ids) # number of locations # create an empty square matrix full of 'None' values x numpy.full((N, N), None) = # adds binary decision variables indicating if transition (i->j) is included in the route for i, j in transitions: x[i, j] solver.add_var(var_type = BINARY) # objective function: minimizes the distance solver.objective = minimize(xsum(dist_matrix[i,j]*x[i][j] for i,j in transitions)) # Arrive/Depart Constraints for i in location_ids: solver += xsum(x[i,j] for j in location_ids - {i}) solver += xsum(x[j,i] for j in location_ids - {i}) # adds a binary decision variable for each location == 1 # exactly 1 arrival == 1 # exactly 1 departure y = [solver.add_var(var_type=INTEGER) for i in location_ids] # adds connectivity constraints for transitions that do not include the startstop for (i, j) in product(location_ids - {startstop}, location_ids - {startstop}): if i != j: # ignores transitions from a location to itself solver y[j] - y[i] >=(N+1)*x[i,j] - N += solver.optimize() #solves the problem # prints the solution if solver.num_solutions: # if a solution was found best_route = [startstop] # stores the best route curr_loc startstop #23 the currently visited location while True: for next loc in location_ids: # for every possible next location if x[curr_loc, next_loc].x == 1: # if x value for the curr_loc->next_loc transition is 1 best_route.append(next_loc) #appends the next location to the route curr_loc=next_loc # visits the next location break if next_loc break == startstop: # exits if route returns to the startstop return best route, solver.objective_value #returns the route and its total distance وزارة التعليم Ministry of Education 292 2024-1446

The following code generates 100 TSP instances with 8 locations and a distance range between 5 and 20. It also uses the brute-force and the MIP solver to solve each instance and reports the percentage for which the two methods reported routes with the same distance. same_count = 0 for i in range(100): dist_matrix, location_ids, startstop-create_problem_instance(8, [5,20]) route1, dist1 = brute force_solver (dist_matrix, location_ids, startstop). route2, dist2 = MIP_solver(dist_matrix, location_ids, startstop) # counts how many times the two solvers produce the same total distance if dist1 == dist2: same_count += 1 print(same_count / 100) 1.0 The results verify that the MIP solver reports the optimal solution for 100% of the problem instances. The following code also demonstrates the speed of the MIP solver, by using it to solve 100 larger instances with 20 locations. import time start = time.time() # starts timer for i in range(100): dist_matrix, location_ids, startstop = create_problem_instance(20, [5,20]) route, dist = MIP_solver(dist_matrix, location_ids, startstop) stop time.time() # stops timer print(stop - start) # prints the elapsed time in seconds 188.90074133872986 Even though the exact execution time will depend on the processing power of the machine that you use to run this Jupyter notebook, it should generally take just a few minutes to compute the solution for all 100 datasets. This is quite impressive, considering that the number of possible routes for each of the 100 instances translates to 19! = 121,645,100,000,000,000 different routes. Such a large number of routes is far beyond the capabilities of the brute-force approach. However, by efficiently searching this massive space of all possible solutions, the MIP solver can find the optimal route quickly. Despite its advantages, mathematical programming also has its limitations. It requires a solid understanding of mathematical modeling and may not be suitable for complex problems where the objective function and constraints are hard to express via mathematical formulas. In addition, even though mathematical programming is much faster than the brute-force approach, it might still be too • slow for large datasets. In such cases, the heuristic approach demonstrated in the previous two lessons presents a much faster alternative. وزارة التعليم Ministry of Education 2024-1446 293

Exercises 1 Explain how mathematical programming can be used to solve complex optimization problems. 2 What are the advantages and disadvantages of the MIP approach for solving optimization problems? وزارة التعليم Ministry of Education 294 2024-1446

3 Analyze two optimization problems that can be solved with mathematical programming and outline their state and decision variables. 4 List three different optimization problems from the family of routing problems. وزارة التعليم Ministry of Education 2024-1446 295

5 You want to create a brute-force solver function for the Traveling Salesman Problem. Complete the following code so that the solver function returns the best route and best total distance. from itertools import permutations def brute force_solver (dist_matrix, location_ids, startstop): # excludes the startstop location location_ids = - · { # generates all possible routes (location permutations) all_routes = } ) best_distance = float('inf') #initializes to the highest possible number best_route = None # best route so far, initialized to None for route in all_routes: # for each route distance = 0 #total distance in this route curr_loc= for next loc in route: distance += curr_loc= distance += #current location [curr_loc, next_loc] # adds the distance of this step #goes the next location [curr_loc, ] # goes to back to the startstop location if distance < best_distance: # if this route has lower distance than the best route best_distance = distance best_route = route وزارة التعليم Ministry of Education 296 2024-1446 # adds the startstop location at the beginning and end of the best route and returns return [startstop] + list(best_route) + [startstop], best_distance

6 You want to create an MIP solver for the Traveling Salesman Problem. Complete the following code so that the decision variables and connectivity constraints are selected correctly. def MIP_solver(dist_matrix, location_ids, startstop): solver = () # creates a solver solver.verbose = 0 # setting this to 1 will print info on the progress of the solver # creates every transition from every location to every other location transitions = list( (location_ids, location_ids)) N = len(location_ids) # number of locations # creates an empty square matrix full of 'None' values x numpy.full((N, N), None) = #adds binary decision variables indicating if transition (i->j) is included in the route for i, j in transitions: x[i, j] = = solver. (var_type= (xsum(dist_matrix[i, j] * x[i][j] for # objective function: minimizes the distance solver.objective = i, j in transitions)) #Arrive/Depart Constraints for i in location_ids: solver += xsum( solver += xsum( # Adds a binary decision variable for each location y = [solver. location_ids] for j in location_ids - {i}) == 1 for j in location_ids - {i}) == 1 (var_type= ) for i in # Adds connectivity constraints for transitions that do not include the startstop for (i, j) in product(location_ids - {startstop}, location_ids {startstop}): if i != j # ignores transitions from a location to itself solver. solver + y[j] - y[i] >= (N + 1) * x[i, j] - N += () # solves the problem وزارة التعليم Ministry of Education 2024-1446 297