Brain.ai Workshop tasks.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Brain.ai Workshop with Visma"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Prerequisites\n",
    "* Have Python 3+ installed on your machine\n",
    "* Install jupyter and check that you can open jupyter with 'jupyter notebook'\n",
    "* If you're able to, work on two monitors"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Acknowlegdements\n",
    "\n",
    "This notebook has been inspired by: https://gist.github.com/hollemawoolpert/612b51f5a489d276a4240f298a3cc9ed#file-postbuild"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Setup"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Install libraries\n",
    "!pip install ortools\n",
    "!pip install gmaps\n",
    "!pip install googlemaps\n",
    "!pip install geopy\n",
    "!pip install numpy\n",
    "!pip install matplotlib"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Enable gmaps extension\n",
    "!jupyter nbextension enable --py --sys-prefix gmaps"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Imports\n",
    "from __future__ import print_function\n",
    "from ortools.constraint_solver import routing_enums_pb2\n",
    "from ortools.constraint_solver import pywrapcp\n",
    "\n",
    "import gmaps\n",
    "import geopy\n",
    "import googlemaps\n",
    "import numpy as np\n",
    "import random\n",
    "import math\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Sections\n",
    "1. Backstory\n",
    "2. Preparation\n",
    "3. Generate Random Baseline Solution\n",
    "5. Find solution manually\n",
    "4. Design and use a Genetic Algorithm\n",
    "5. Solve using OR Tools"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Backstory\n",
    "\n",
    "The rollout of the Covid-19 vaccine in Norway has finally started to pick up pace. A large batch of vaccines has been delivered to St. Olavs Hospital, and they need to be distributed to the vaccination sites. The sites include your favorite examination sites!\n",
    "\n",
    "The vaccines have to be kept at 203 Kelvin, so special pupose vehicles are needed. There is a limited amount of these vehicles, and the vaccines need to be distributed swiftly.\n",
    "\n",
    "Your have been assigned with the task of determining which route each vehicle should drive in order to minimize the time used to distribute the vaccines!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Preparation\n",
    "\n",
    "#### Create account and get API key to use Google Distance Matrix API and Maps API\n",
    "\n",
    "We wish to use Google's Distance Matrix API to efficiently get distances between places. We also need some other APIs. To use these APIs, you have to register your credit card information. You will get enough free credits for the purpose of this workshop, and Google promises you won't be billed without explicitly accepting that. However, if you are not comfortable with this, you can use the provided distance matrix and skip the map parts.\n",
    "\n",
    "If you want to use the APIs, you need an account and an API key:\n",
    "\n",
    "1. If you have never used the Google Cloud Console to create a billing account or a project, go to https://developers.google.com/maps/gmp-get-started and click \"Get Started\". Follow the steps.\n",
    "2. Check the \"Routes\" and \"Places\" products and click \"Enable\".\n",
    "3. Copy your API key. \n",
    "4. Try opening this in your browser:\n",
    "[api-request](https://maps.googleapis.com/maps/api/distancematrix/json?origins=studentersamfundet+trondheim|gl%C3%B8shaugen+trondheim|downtown+trondheim&destinations=studentersamfundet+trondheim|gl%C3%B8shaugen+trondheim|downtown+trondheim&mode=walking&key=API_KEY), replacing `API_KEY` at the end with your API key. You should get a json distance matrix as response.\n",
    "5. Assign your API key to `API_KEY` in the cell below.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "API_KEY = None"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Using the Maps API\n",
    "\n",
    "Let's set some locations to visit and plot them in a map. Feel free to change locations, and run the cell below.\n",
    "\n",
    "NB: The distance matrix API used later has a limit of 10 locations per request."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "location_names = [\n",
    "    \"St. Olavs Hospital\",\n",
    "    \"Nidaros legekontor\",\n",
    "    \"Møllenberg legesenter\",\n",
    "    \"Sit Gløshaugen legesenter\",\n",
    "    \"Trondheim spektrum\",\n",
    "    \"Byåsen legesenter\",\n",
    "    \"Valentinlyst legesenter\",\n",
    "    \"Utleirahallen\",\n",
    "    \"Dragvoll\",\n",
    "    \"Sluppen\"\n",
    "]\n",
    "depot = 0  # start and end location index\n",
    "\n",
    "gmaps.configure(api_key=API_KEY)\n",
    "geocoder = geopy.geocoders.GoogleV3(api_key=API_KEY)\n",
    "\n",
    "locations = [geocoder.geocode(name + \" Trondheim\") for name in location_names]\n",
    "coordinates = [(loc.latitude, loc.longitude) for loc in locations]\n",
    "\n",
    "depot_layer = gmaps.symbol_layer(\n",
    "    [coordinates[depot]],\n",
    "    hover_text=\"Start: \" + location_names[depot],\n",
    "    fill_color='white',\n",
    "    stroke_color='red',\n",
    "    scale=8\n",
    ")\n",
    "\n",
    "stops_layer = gmaps.symbol_layer(\n",
    "    coordinates[0:depot] + coordinates[depot+1:],\n",
    "    hover_text=location_names[0:depot] + location_names[depot+1:],\n",
    "    fill_color='white',\n",
    "    stroke_color='black',\n",
    "    scale=4\n",
    ")\n",
    "\n",
    "fig = gmaps.figure()\n",
    "fig.add_layer(depot_layer)\n",
    "fig.add_layer(stops_layer)\n",
    "\n",
    "fig"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Using the Distance Matrix API\n",
    "\n",
    "Next, we use the Distance Matrix api to get our distance matrix."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def build_distance_matrix(coordinates, travel_mode, measure='duration'):\n",
    "    \"\"\"Get distance matrix for coordinates.\n",
    "\n",
    "    Element (i,j) is the distance between the i-th and j-th coordinate.\n",
    "\n",
    "    Args:\n",
    "        coordinates: List of two-tuples of latidude/longitude-coordinates.\n",
    "        measure (string): 'duration' for seconds or 'distance' for meters.\n",
    "\n",
    "    Returns:\n",
    "        Distance matrix (list of lists) of size n*n, where n==len(coordinates).\n",
    "    \"\"\"\n",
    "    client = googlemaps.Client(key=API_KEY)\n",
    "    dm_response = client.distance_matrix(origins=coordinates, destinations=coordinates, mode=travel_mode)\n",
    "    dm_rows = [row['elements'] for row in dm_response['rows']]\n",
    "    distance_matrix = [[item[measure]['value'] for item in dm_row] for dm_row in dm_rows]\n",
    "    return distance_matrix\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "travel_mode = 'driving'\n",
    "distance_matrix = build_distance_matrix(coordinates, travel_mode)\n",
    "distance_matrix"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Baseline Solution\n",
    "We start by creating a random solution. This will serve as a baseline; we should be able to do much better than random. A random solution will also solve as a starting point when we solve the problem using a genetic algorithm later."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create a random solution\n",
    "def create_random_solution(locations, depot, num_vehicles):\n",
    "    \"\"\"Create a random solution.\n",
    "\n",
    "    Create a random solution to our vehicle routing problem.\n",
    "\n",
    "    Args:\n",
    "        locations: List of location objects.\n",
    "        depot (int): Index of the depot.\n",
    "        num_vehicles (int): Number of vehicles.\n",
    "\n",
    "    Returns:\n",
    "        list of ints represeting a solution (see \"Chromosome Representation\" below).\n",
    "    \"\"\"\n",
    "    num_locations = len(locations)\n",
    "    solution = list(range(num_locations))  # 1 depot\n",
    "    for _ in range(num_vehicles - 2):\n",
    "        solution.append(depot)  # num_vehicles-1 depots\n",
    "    random.shuffle(solution)  # note that this is in-place\n",
    "    return solution"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "random_solution = create_random_solution(locations, depot, num_vehicles=3)\n",
    "random_solution"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Evaluate the solution\n",
    "Let's create some methods for helping us evaluate solutions to the VRP."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_routes_dict(solution, depot):\n",
    "    \"\"\"Get a dictionary representing the solution.\n",
    "\n",
    "    Get a representation of the routes that is simpler to read and work with.\n",
    "\n",
    "    Args:\n",
    "        solution (list of ints): A solution chromosome.\n",
    "        depot (int): Index of the depot.\n",
    "\n",
    "    Returns:\n",
    "        A dictionary with one element of ordered indices for each vehicle route.\n",
    "    \"\"\"\n",
    "    routes = dict()\n",
    "    route_num = 0\n",
    "    route = [depot]\n",
    "    for i in solution:\n",
    "        route.append(i)\n",
    "        if i == depot:\n",
    "            routes[route_num] = route\n",
    "            route_num += 1\n",
    "            route = [depot]\n",
    "    route.append(depot)\n",
    "    routes[route_num] = route\n",
    "    return routes\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "random_routes = get_routes_dict(random_solution, depot)\n",
    "random_routes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_route_distance(route, distance_matrix):\n",
    "    \"\"\"Get the distance of a route.\n",
    "\n",
    "    Args:\n",
    "        route (list of ints): A vehicle route.\n",
    "        distance_matrix (list of lists of ints): A distance matrix.\n",
    "\n",
    "    Returns:\n",
    "        int: Distance of the route.\n",
    "    \"\"\"\n",
    "    distance = 0\n",
    "    for k in range(len(route) - 1):\n",
    "        i_from = route[k]\n",
    "        i_to = route[k + 1]\n",
    "        distance += distance_matrix[i_from][i_to]\n",
    "    return distance\n",
    "\n",
    "\n",
    "def get_max_route_distance(routes, distance_matrix):\n",
    "    \"\"\"Get the longest route distance.\n",
    "\n",
    "    Get the distance of the longest route in a set of routes.\n",
    "\n",
    "    Args:\n",
    "        routes (dict): Vehicle routes.\n",
    "        distance_matrix (list of lists of ints): A distance matrix.\n",
    "\n",
    "    Returns:\n",
    "        int: Distance of the longest route.\n",
    "    \"\"\"\n",
    "    return max([get_route_distance(route, distance_matrix) for route in routes.values()])\n",
    "\n",
    "\n",
    "def pretty_print(routes, location_names, distance_matrix):\n",
    "    \"\"\"Create a nice printout of a set of routes.\n",
    "\n",
    "    Args:\n",
    "        routes (dict): Vehicle routes.\n",
    "        location_names (list of str): Names of the locations.\n",
    "        distance_matrix (list of lists of ints): A distance matrix.\n",
    "    \"\"\"\n",
    "    for vehicle, route in routes.items():\n",
    "        route_distance = get_route_distance(route, distance_matrix)\n",
    "        if route_distance == 0:\n",
    "            print(f\"Vehicle {vehicle} not used\\n\")\n",
    "            continue\n",
    "        route_str = str(location_names[route[0]])\n",
    "        for i in route[1:]:\n",
    "            route_str += \" -> \" + location_names[i]\n",
    "        print(\n",
    "            f\"Route for vehicle {vehicle}:\\n\" +\n",
    "            route_str + \"\\n\" +\n",
    "            f\"Distance of route: {route_distance}\\n\"\n",
    "        )\n",
    "    print(f\"Maximum of the route distances: {get_max_route_distance(routes, distance_matrix)}\")\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pretty_print(random_routes, location_names, distance_matrix)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def map_solution(routes, coordinates, depot):\n",
    "    \"\"\"Create a map figure showing vehicle routes.\n",
    "\n",
    "    Args:\n",
    "        routes (dict): Vehicle routes.\n",
    "        location_names (list of str): Names of the locations.append\n",
    "        distance_matrix (list of lists of ints): A distance matrix.\n",
    "\n",
    "    Returns:\n",
    "        gmaps figure object.\n",
    "    \"\"\"\n",
    "    fig = gmaps.figure()\n",
    "    colors = ['blue','red','green','#800080','#000080','#008080']\n",
    "    assert len(colors) >= len(routes), \"More colors needed\"\n",
    "    for vehicle_id in routes:\n",
    "        waypoints = []\n",
    "        depot_coordinates = coordinates[depot]\n",
    "        \n",
    "        # skip depot (occupies first and last index)\n",
    "        for coordinate_index in routes[vehicle_id][1:-1]:\n",
    "            waypoints.append(coordinates[coordinate_index])\n",
    "        \n",
    "        if len(waypoints) == 0:\n",
    "            print('Empty route:', vehicle_id)\n",
    "        else:\n",
    "            route_layer = gmaps.directions_layer(\n",
    "                depot_coordinates, waypoints[-1], waypoints=waypoints[0:-1], show_markers=True,\n",
    "                stroke_color=colors[vehicle_id], stroke_weight=5, stroke_opacity=0.5,\n",
    "                travel_mode=travel_mode.upper())\n",
    "            fig.add_layer(route_layer)\n",
    "            \n",
    "            # complete the route from last location to depot\n",
    "            return_layer = gmaps.directions_layer(\n",
    "                waypoints[-1], depot_coordinates, show_markers=False,\n",
    "                stroke_color=colors[vehicle_id], stroke_weight=5, stroke_opacity=0.5,\n",
    "                travel_mode=travel_mode.upper())\n",
    "            fig.add_layer(return_layer)\n",
    "    \n",
    "    return fig\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "map_solution(random_routes, coordinates, depot)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Does this look like a good solution? In that case, try again...\n",
    "\n",
    "Combinatorial (or discrete) optimization problems, such as this one, have large solution spaces even for small inputs. So a random solution will most likely not be a good one.\n",
    "\n",
    "BONUS TASK: Quantify the complexity of the VRP."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. Manually solve the problem\n",
    "Let's try to solve the problem manually by inspecting the map and decide which routes the vehicles should make.\n",
    "\n",
    "TASK: \n",
    "\n",
    "1. Click the following link to open up a Google Drawing of the VRP instance. https://docs.google.com/drawings/d/1L0bQ9KUoX3CBnvLtdrDxZ2nC4Wq50ZfdqhvpKKT10mc/edit?usp=sharing\n",
    "2. Click 'File' in the top left corner and then 'Make a copy'. \n",
    "3. Inspect the problem and try to draw the best route for each of the three vehicles (Select the dropdown menu at 'Line' and select 'Scribble' to draw freely. \n",
    "4. Insert your solution into manual_solution below. This represents the routes taken by each vechicle. \n",
    "5. Evaluate your solution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Insert your solution here. \n",
    "manual_routes = {0: [0, 0], 1: [0, 0], 2: [0, 0]}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pretty_print(manual_routes, location_names, distance_matrix)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "map_solution(manual_routes, coordinates, depot)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. Solving with a genetic algorithm\n",
    "\n",
    "### Chromosome Representation\n",
    "\n",
    "A solution to this problem can be represented in a one-dimensional matrix of size (1, n) where n = (k-1) + c, where k is the number of vehicles and c is the number of destinations. \n",
    "\n",
    "A solution is represented as a permutation of integers, where each integer represents a location (destination or depot). The depot values is used to delimit different routes. \n",
    "\n",
    "Example: \n",
    "\n",
    "[2,1,0,4,3,0,5] \n",
    "\n",
    "Vehicle 1: 0-2-1-0\n",
    "Vehicle 2: 0-4-3-0\n",
    "Vehicle 3: 0-5-0\n",
    "\n",
    "Each solution will have to end with a zero (go back to depot), so we can omit this from the solution representation. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "class Chromosome:\n",
    "    \"\"\"The class that contains a chromosome, i.e. a solution to the problem\n",
    "    \n",
    "    Attributes:\n",
    "        route_representation (list): List representing the route of the vehicles\n",
    "        route_dict (dict): A dictionary with one element of ordered indices for each vehicle route.\n",
    "        longest_route_distance (float): Longest route a vehicle must make in a solution\n",
    "        fitness (float): Fitness value of the chromosome\n",
    "    \"\"\"\n",
    "    \n",
    "    def __init__(self, route_representation, route_dict):\n",
    "        \"\"\"Initializes a chromosome.\n",
    "        \"\"\"\n",
    "        self.route_representation = route_representation\n",
    "        self.route_dict = route_dict\n",
    "        self.longest_route_distance = None\n",
    "        self.fitness = None\n",
    "        \n",
    "    def __str__(self):\n",
    "        \"\"\"For outputting a string of the Chromosome attributes when printing.\n",
    "        \"\"\"\n",
    "        output_string = str(self.route_representation) + '\\n' + str(self.route_dict) + '\\n'\n",
    "        output_string += 'fitness: ' + str(self.fitness) + '\\nlongest_route_distance: ' + stl(self.longest_route_distance)\n",
    "        return output_string\n",
    "        "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Algorithm\n",
    "Implement the GA algorithm. Study the flow chart below and the comments in the code and make sure you understand how the algorithm works."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"ga.png\">"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def ga_solver(\n",
    "    locations,\n",
    "    depot,\n",
    "    num_vehicles,\n",
    "    distance_matrix,\n",
    "    population_size=50,\n",
    "    num_generations=50,\n",
    "    elite_size=2,\n",
    "    tournament_k=3,\n",
    "    tournament_p=0.95,\n",
    "    mutation_rate=0.8,\n",
    "    mutation_swaps=1,\n",
    "):\n",
    "    \"\"\"Main GA algorithm, performing all the steps of the genetic algorithm\n",
    "    \n",
    "    Args: \n",
    "        locations: List of location objects.\n",
    "        depot (int): Index of the depot.\n",
    "        num_vehicles (int): Number of vehicles.\n",
    "        distance_matrix (list of lists of ints): A distance matrix.\n",
    "        population_size (int): Size of the population.\n",
    "        num_generations (int): Number of generations to evolve.\n",
    "        elite_size (int): Number of best solutions to bring directly into next generation.\n",
    "        tournament_k (int): Number of competitors in the tournament selection.\n",
    "        tournament_p (float): Probability of the best competitor winning the tournament.\n",
    "        mutation_rate (float): Probability of performing mutation on a Chromosome.\n",
    "        mutation_swaps (int): Number of swaps to do in swap mutation.\n",
    "        \n",
    "    Returns:\n",
    "\tdict: The chromosome with the best solution.\n",
    "    \"\"\"\n",
    "    # Create initial random population\n",
    "    population = initialize_population(locations, depot, num_vehicles, population_size)\n",
    "    \n",
    "    # Calculate fitness of population\n",
    "    for chromosome in population:\n",
    "        calculate_fitness(chromosome, distance_matrix)\n",
    "\n",
    "    generation_bests = {'fitness': [], 'max_distance': []}\n",
    "\n",
    "    # Loop for evolving the populations through generations\n",
    "    for generation in range(num_generations):\n",
    "        \n",
    "        # Generate new generation\n",
    "        population = evolve_population(\n",
    "            population, elite_size, tournament_k, tournament_p, mutation_rate, mutation_swaps, depot\n",
    "        )\n",
    "        \n",
    "        # calculate fitness\n",
    "        for chromosome in population:\n",
    "            calculate_fitness(chromosome, distance_matrix)\n",
    "        \n",
    "        # store best values\n",
    "        max_distances = np.array([chromosome.longest_route_distance for chromosome in population])\n",
    "        best_i = max_distances.argmin()\n",
    "        chromosome = population[best_i]\n",
    "        generation_bests['fitness'].append(chromosome.fitness)\n",
    "        generation_bests['max_distance'].append(chromosome.longest_route_distance)\n",
    "\n",
    "    return population[best_i], generation_bests\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Initialize population\n",
    "Random initial population - Use the same approach as in generating a random solution.\n",
    "\n",
    "TASK: Implement the content of the for loop in the method below. Hint: Use `create_random_solution()` and `get_routes_dict()`. You might want to test that the method works as expected."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def initialize_population(locations, depot, num_vehicles, population_size):\n",
    "    \"\"\"Initialize population of Chromosome instances\n",
    "    \n",
    "    Args:\n",
    "        locations (list): Location objects.\n",
    "        depot (int): Index of the depot.\n",
    "        num_vehicles (int): Number of vehicles.\n",
    "        population_size (int): Size of the population.\n",
    "        \n",
    "    Returns\n",
    "        List of Chromosome instances.\n",
    "    \"\"\"\n",
    "    population = []\n",
    "    for i in range(population_size):\n",
    "\n",
    "        ## START CODE HERE ##\n",
    "\n",
    "        # Create a random solution\n",
    "        route_representation = None\n",
    "\n",
    "        # Convert to dictionary representation\n",
    "        route_dict = None\n",
    "\n",
    "        # Create a chromosome of the random solution and append it to the solution\n",
    "\n",
    "        ## END CODE HERE ##\n",
    "\n",
    "    assert len(population) == population_size\n",
    "    return population"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Evolve population\n",
    "\n",
    "Evolving the population involves some steps we have not implemented yet: getting elites, (tournament) selection and mutating the new generation. We will get to these methods below. You can input the missing code immediately by just looking at what parameters are needed for the mentioned methods, or you can skip the implementation task here for now and come back to it later.\n",
    "\n",
    "For now, the most important thing is that you read through `evolve_population()` to get an understanding of the main algorithm.\n",
    "\n",
    "TASK: Implement the indicated steps in `evolve_population()` below. Hint: Use `get_elites()`, `tournament_selection()` and `mutate_population()`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def evolve_population(\n",
    "    population, elite_size, tournament_k, tournament_p, mutation_rate, mutation_swaps, depot\n",
    "): \n",
    "    \"\"\"Take in the population and perform all the steps to return a new generation. \n",
    "    \n",
    "    1. Add elites\n",
    "    2. Select parents\n",
    "    3. Crossover\n",
    "    4. Mutation\n",
    "    5. Return new generation\n",
    "    \n",
    "    Args: \n",
    "        population (list): List of Chromosome instances\n",
    "        depot (int): Index of the depot.\n",
    "        elite_size (int): Number of best solutions to bring directly into next generation.\n",
    "        tournament_k (int): Number of competitors in the tournament selection.\n",
    "        tournament_p (float): Probability of the best competitor winning the tournament.\n",
    "        mutation_rate (float): Probability of performing mutation on a Chromosome.\n",
    "        mutation_swaps (int): Number of swaps to do in swap mutation.\n",
    "        \n",
    "    Returns: \n",
    "        List of Chromosome instances of the next iteration.\n",
    "    \"\"\"\n",
    "    next_generation = []\n",
    "    \n",
    "    # Add elites\n",
    "    if elite_size > 0:\n",
    "\n",
    "        ## START CODE HERE ##\n",
    "        elites = None\n",
    "        next_generation = None\n",
    "        ## END CODE HERE ##\n",
    "\n",
    "    # Parent (tournament) selection\n",
    "    mating_pool_size = len(population) - elite_size\n",
    "\n",
    "    ## START CODE HERE ##\n",
    "    mating_pool = None\n",
    "    ## END CODE HERE ##\n",
    "    \n",
    "    # Generate children\n",
    "    # NOTE: Children are just a list of the route_representations, not Chromosome instances\n",
    "    children = breed_parents(mating_pool, depot)\n",
    "    \n",
    "    # perform mutation\n",
    "    \n",
    "    ## START CODE HERE ##\n",
    "    mutated_route_representations = None\n",
    "    ## END CODE HERE ##\n",
    "    \n",
    "    # transform the mutated route representations into Chromosome instances\n",
    "    for route_representation in mutated_route_representations:\n",
    "        route_dict = get_routes_dict(route_representation, depot)\n",
    "        next_generation.append(Chromosome(route_representation, route_dict))\n",
    "\n",
    "    assert len(next_generation) == len(population)\n",
    "    return next_generation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Fitness evaluation\n",
    "Find the longest route in the solution. \n",
    "\n",
    "Here, if we get imbalanced routes, we can add a overcapacity penalty for the number of nurses that has to travel to more than roof(c/k) destinations. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def calculate_fitness(chromosome, distance_matrix):\n",
    "    \"\"\"Calculate the fitness and assign it to the fitness attribute of the chromosome.\n",
    "    \n",
    "    Args:\n",
    "        chromosome (Chromosome): The Chromosome instance.\n",
    "        distance_matrix (list of lists of ints): A distance matrix.\n",
    "    \"\"\"\n",
    "    chromosome.longest_route_distance = get_max_route_distance(chromosome.route_dict, distance_matrix)\n",
    "    chromosome.fitness = 1/chromosome.longest_route_distance"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Elites"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_elites(population, elite_size): \n",
    "    \"\"\"Get the best chromosomes from the population.\n",
    "\n",
    "    Args:\n",
    "        population (list): List of the Chromosome instances.\n",
    "        elite_size (int): Number of best solutions to bring directly into next generation.\n",
    "        \n",
    "    Returns: \n",
    "        List of the elite Chromosomes.\n",
    "    \"\"\"\n",
    "    # sort population descendingly in-place according to fitness value\n",
    "    population.sort(key=lambda x: x.fitness, reverse=True)\n",
    "    elites = population[:elite_size]\n",
    "    return elites"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Selection\n",
    "Select mating pool using tournament selection. See the illustration below and the comments in the code for an explanation.\n",
    "\n",
    "TASK: Finish the tournament selection method below. Hint: `random.sample()` and `np.random.choice()` are methods that might come in handy."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"tournament_selection.png\">"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def tournament_selection(population, tournament_k, tournament_p, mating_pool_size): \n",
    "    \"\"\"Takes in the entire population and returns a mating pool of the same size. \n",
    "    \n",
    "    The mating pool is the parents selected for breeding.\n",
    "    \n",
    "    Args:\n",
    "        population (list): List of the Chromosome instances.\n",
    "        tournament_k (int): Number of competitors in the tournament selection.\n",
    "        tournament_p (float): Probability of the best competitor winning the tournament.\n",
    "        mating_pool_size: (int): The size of the mating pool, i.e. the number of parents selected for breeding.\n",
    "    \n",
    "    Returns: \n",
    "        List of Chromosome instances that form the mating pool.\n",
    "    \"\"\"\n",
    "    mating_pool = []\n",
    "    population.sort(key=lambda x: x.fitness, reverse=True)\n",
    "    \n",
    "    # While the mating pool is not full\n",
    "    while len(mating_pool) < mating_pool_size:\n",
    "        \n",
    "        ## START CODE HERE ##\n",
    "\n",
    "        # Select competitors at random\n",
    "        competitors = None\n",
    "        \n",
    "        # Sort competitors according to fitness\n",
    "        competitors = None\n",
    "        \n",
    "        # Get the probabilities for the competitors. probability: p * (1-p)**i\n",
    "        probabilities = np.empty(tournament_k)\n",
    "        for i in range(tournament_k):\n",
    "            probabilities[i] = None\n",
    "        \n",
    "        # Normalize the probabilities (if you feel like it is necessary...)\n",
    "        probabilities = None\n",
    "        \n",
    "        # Choose the winner according to the propabilities\n",
    "        winner = None\n",
    "\n",
    "        ## END CODE HERE ##\n",
    "        \n",
    "        mating_pool.append(winner)\n",
    "        \n",
    "    assert len(mating_pool) == mating_pool_size\n",
    "    return mating_pool"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Crossover\n",
    "Use ordered crossover. Our chromosome representation makes this relatively straight-forward, but we still need to watch our step when implementing it. We have done this part for you.\n",
    "\n",
    "BONUS TASK: Try using egde combination crossover."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"ordered_crossover.png\">\n",
    "\n",
    "From: https://towardsdatascience.com/evolution-of-a-salesman-a-complete-genetic-algorithm-tutorial-for-python-6fe5d2b3ca35"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def breed_parents(mating_pool, depot):\n",
    "    \"\"\"Performs the crossover and returns the list of children chromosomes.\n",
    "    \n",
    "    Args:\n",
    "        mating_pool (list): List of Chromosome instances.\n",
    "        depot (int): Index of the depot.\n",
    "    \n",
    "    Returns:\n",
    "        List of Chromosome instances of the children produced by crossover.\n",
    "    \"\"\"\n",
    "    children = []\n",
    "    \n",
    "    # shuffle the order of the parents in-place\n",
    "    random.shuffle(mating_pool)\n",
    "    \n",
    "    # each crossover produces two children, so the crossover should be executed (size/2) times\n",
    "    num_iterations = math.ceil(len(mating_pool)/2)\n",
    "    for i in range(num_iterations):\n",
    "        parent1 = mating_pool[i]\n",
    "        parent2 = mating_pool[len(mating_pool)-i-1]\n",
    "        child1, child2 = ordered_crossover(parent1, parent2, depot)\n",
    "        children.append(child1)\n",
    "        children.append(child2)\n",
    "\n",
    "    # If mating pool is of odd length, we have one excess child. So we leave it behind...\n",
    "    if (len(children) > len(mating_pool)):\n",
    "        children.pop()\n",
    "    \n",
    "    return children"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def ordered_crossover(parent1_chromosome, parent2_chromosome, depot):\n",
    "    \"\"\"Performs ordered crossover and returns the route_representations of the two children.\n",
    "    \n",
    "    Args:\n",
    "        parent1_chromosome (Chromosome): Parent 1 Chromosome.\n",
    "        parent2_chromosome (Chromosome): Parent 2 Chromosome.\n",
    "        depot (int): Index of the depot.\n",
    "        \n",
    "    Returns: \n",
    "        Two lists with route_representations of the children (not Chromosome instances).\n",
    "    \"\"\"\n",
    "    # copy the route representations from the parents\n",
    "    parent1 = parent1_chromosome.route_representation.copy()\n",
    "    parent2 = parent2_chromosome.route_representation.copy()\n",
    "    \n",
    "    # Replace the depot values with strings to make them unique (ordered crossover requires unique cells)\n",
    "    def replace_depot_values(route):\n",
    "        counter = 1\n",
    "        for idx, value in enumerate(route):\n",
    "            if value == depot:\n",
    "                route[idx] = 'd' + str(counter)\n",
    "                counter = counter + 1\n",
    "    replace_depot_values(parent1)\n",
    "    replace_depot_values(parent2)\n",
    "\n",
    "    # generate random index slicing points\n",
    "    length = len(parent1)\n",
    "    start_slice = random.randint(0, length - 1)\n",
    "    end_slice = random.randint(start_slice + 1, length)\n",
    "    \n",
    "    # get the values from start (included) to stop (excluded)\n",
    "    child1_middle = parent1[start_slice:end_slice]\n",
    "    child2_middle = parent2[start_slice:end_slice]\n",
    "    \n",
    "    # get the values to fill in from the other parent in each child\n",
    "    child1_fill = [i for i in parent2 if i not in child1_middle]\n",
    "    child2_fill = [i for i in parent1 if i not in child2_middle]\n",
    "    \n",
    "    # put together the fill values and the unchanged part between the slice points\n",
    "    child1 = child1_fill[:start_slice] + child1_middle + child1_fill[start_slice:]\n",
    "    child2 = child2_fill[:start_slice] + child2_middle + child2_fill[start_slice:]\n",
    "    \n",
    "    # transform the depot values back to zeros\n",
    "    child1 = [0 if str(value).startswith('d') else value for value in child1]\n",
    "    child2 = [0 if str(value).startswith('d') else value for value in child2]\n",
    "    \n",
    "    return child1, child2\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Mutation\n",
    "Swap two places (ore more) with each other in one individual.\n",
    "\n",
    "TASK: Finish `swap_mutate()` below."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"swap_mutation.jpg\">\n",
    "\n",
    "From: https://www.theprojectspot.com/tutorial-post/applying-a-genetic-algorithm-to-the-travelling-salesman-problem/5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def mutate_population(children, mutation_rate, mutation_swaps):\n",
    "    \"\"\"Mutate population.\n",
    "\n",
    "    Takes in children which is a list of route_representations of the children\n",
    "    (because we don't bother with computing fitness before we have done mutation).\n",
    "    \n",
    "    Args: \n",
    "        children (list): List of lists containing the route representations of the children.\n",
    "        mutation_rate (float): Probability of performing mutation on a Chromosome.\n",
    "        mutation_swaps (int): Number of swaps to do in swap mutation.\n",
    "    \n",
    "    Returns:\n",
    "        List of lists with the route representations for the next generation.\n",
    "    \"\"\"\n",
    "    mutated_population = []\n",
    "    \n",
    "    for child in children:\n",
    "        if (random.random() < mutation_rate):\n",
    "            mutated_child = swap_mutate(child, mutation_swaps)\n",
    "            mutated_population.append(mutated_child)\n",
    "        else:\n",
    "            mutated_population.append(child)\n",
    "    \n",
    "    return mutated_population"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def swap_mutate(route_representation, number_of_swaps):\n",
    "    \"\"\"Run swap mutation on a route_representation list.\n",
    "    \n",
    "    Args:\n",
    "        route_representation (list): List of the route representation for a solution.\n",
    "        number_of_swaps (int): Number of swaps to do in swap mutation.\n",
    "    \n",
    "    Returns: \n",
    "        List containing the route representation of the mutated child.\n",
    "    \"\"\"\n",
    "    for _ in range(number_of_swaps): \n",
    "        \n",
    "        ## START CODE HERE ##\n",
    "\n",
    "        # choose two random indexes to swap\n",
    "        # INSERT CODE HERE\n",
    "        \n",
    "        # temporarily store the values in these indices\n",
    "        # INSERT CODE HERE\n",
    "        \n",
    "        # switch the values in the two indices\n",
    "        # INSERT CODE HERE\n",
    "\n",
    "        ## END CODE HERE ##\n",
    "        \n",
    "    return route_representation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Solve using the genetic algorithm"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "ga_solution, history = ga_solver(locations, depot, num_vehicles=3, distance_matrix=distance_matrix)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot how the solution improves through the generations\n",
    "plt.figure()\n",
    "plt.plot(history['fitness'])\n",
    "plt.xlabel('generation')\n",
    "plt.ylabel('fitness')\n",
    "plt.figure()\n",
    "plt.plot(history['max_distance'], )\n",
    "plt.xlabel('generation')\n",
    "plt.ylabel('max distance')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ga_routes = ga_solution.route_dict\n",
    "pretty_print(ga_routes, location_names, distance_matrix)\n",
    "map_solution(ga_routes, coordinates, depot)"
   ]
  },
  {
   "source": [
    "How does this solution compare to your manual solution? Is there anything you could do to increase the likelihood of finding a good solution - and at what cost?"
   ],
   "cell_type": "markdown",
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6. Solving with OR Tools\n",
    "\n",
    "It this section we will solve the same problem as before, but with the \"plug-and-play\" Google OR-tools python package. \n",
    "\n",
    "Comparing the amount of code and work we need when using OR-tools, compared to when designing an algorithm ourselves, illustrates the importance of not reinventing the wheel. More often than not the problems you are trying to solve are solved by someone else before you. This has to be taken into account before putting in a huge effort to build something from scratch.\n",
    "\n",
    "### Model\n",
    "\n",
    "First, we build the OR-tools model. See the [OR-tools documentation](https://developers.google.com/optimization/routing/vrp) for a more detailed explanations of the steps."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def distance_callback(from_index, to_index):\n",
    "    \"\"\"Returns the distance between the two nodes.\"\"\"\n",
    "    # Convert from routing variable Index to distance matrix NodeIndex\n",
    "    from_node = manager.IndexToNode(from_index)\n",
    "    to_node = manager.IndexToNode(to_index)\n",
    "    return distance_matrix[from_node][to_node]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "num_vehicles = 3\n",
    "\n",
    "# Create the routing index manager\n",
    "manager = pywrapcp.RoutingIndexManager(len(distance_matrix), num_vehicles, depot)\n",
    "\n",
    "# Create Routing Model\n",
    "routing = pywrapcp.RoutingModel(manager)\n",
    "\n",
    "# Register the distance callback. The solver will use this to get distances\n",
    "transit_callback_index = routing.RegisterTransitCallback(distance_callback)\n",
    "\n",
    "# Define cost of each arc - which is the distance\n",
    "routing.SetArcCostEvaluatorOfAllVehicles(transit_callback_index)\n",
    "\n",
    "# Add duration contraint as a dimension\n",
    "dimension_name = 'Duration'\n",
    "max_distance = max([max(distance_row) for distance_row in distance_matrix])\n",
    "max_vehicle_travel_distance = max_distance * 2 * len(distance_matrix) // num_vehicles\n",
    "routing.AddDimension(\n",
    "    transit_callback_index,\n",
    "    0,  # no slack\n",
    "    max_vehicle_travel_distance,\n",
    "    True,  # start cumul to zero\n",
    "    dimension_name,\n",
    ")\n",
    "distance_dimension = routing.GetDimensionOrDie(dimension_name)\n",
    "distance_dimension.SetGlobalSpanCostCoefficient(100)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Solve using OR-tools\n",
    "\n",
    "We find a solution using the model defined above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Set a baseline solution with heuristic\n",
    "search_parameters = pywrapcp.DefaultRoutingSearchParameters()\n",
    "search_parameters.first_solution_strategy = (\n",
    "    routing_enums_pb2.FirstSolutionStrategy.PATH_CHEAPEST_ARC)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Solve the problem, i.e. minimize duration of longest route\n",
    "or_solution_object = routing.SolveWithParameters(search_parameters)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Evaluate the solution\n",
    "\n",
    "TASK: Fill in the missing code for printing and mapping the solution below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def extract_routes(num_vehicles, manager, routing, solution):\n",
    "    \"\"\"Extract routes dictionary from solution object.\n",
    "    \"\"\"\n",
    "    routes = {}\n",
    "    for vehicle_id in range(num_vehicles):\n",
    "        routes[vehicle_id] = []\n",
    "        index = routing.Start(vehicle_id)\n",
    "        while not routing.IsEnd(index):\n",
    "            routes[vehicle_id].append(manager.IndexToNode(index))\n",
    "            previous_index = index\n",
    "            index = solution.Value(routing.NextVar(index))\n",
    "        routes[vehicle_id].append(manager.IndexToNode(index))\n",
    "    return routes\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "or_routes = extract_routes(num_vehicles, manager, routing, or_solution_object)\n",
    "\n",
    "## START CODE HERE ##\n",
    "\n",
    "# Print and map the solution\n",
    "\n",
    "## END CODE HERE ##"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Some obligatory reflection:\n",
    "- How does this solution compare to you GA solution?\n",
    "- Are any of them optimal?\n",
    "- Why does/doesn't OR-tools provide on optimal solution?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Compare solutions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "random_max_dist = get_max_route_distance(random_routes, distance_matrix)\n",
    "manual_max_dist = get_max_route_distance(manual_routes, distance_matrix)\n",
    "ga_max_dist = get_max_route_distance(ga_routes, distance_matrix)\n",
    "or_max_dist = get_max_route_distance(or_routes, distance_matrix)\n",
    "results = [random_max_dist, manual_max_dist, ga_max_dist, or_max_dist]\n",
    "\n",
    "ind = range(4)\n",
    "_ = plt.bar(ind, results)\n",
    "_ = plt.ylabel(\"Longest route distance\")\n",
    "_ = plt.xticks(ind, ['Random', 'Manual', 'GA', 'OR-tools'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "How do the solutions compare? Which approach is the most scalable? Below we have outlined some extra constraints/dimensions complicating things further. Which approach would you go for if you were to take these into account?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 7. Extra challenges\n",
    "\n",
    "Here are some ideas for extra challenges to have a go at if you finished early, or if you want to play some more around with different parts of what we have gone through. They should also illustrate well that there is almost no limit for how complex route optimization can become.\n",
    "\n",
    "- Host this notebook online with [binder](https://mybinder.org/)\n",
    "- Add stop/offload times. Hint: Start with the distance matrix.\n",
    "- Add more locations. Hint: The free distance matrix API has a limit for how many distances you can get in a single request, but we are free to post multiple requests :)\n",
    "- Add capacity constraint: Each location needs a number of vaccines, and each vehicle has a limited capacity. Use OR-tools and/or add support to the genetic algorithm. \n",
    "- Maybe the vehicles can refill when they are empty?\n",
    "- Time windows: Some locations can only receive the vaccines at certain times.\n",
    "- Handling the vaccines is complex, so at some locations two vehicles must be there together for the entire stop time before they can move on."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}