#!/usr/bin/env python
# coding: utf-8

# # Routing, Arbitrage
# 
# Chapter 8 of [Real World Algorithms](https://mitpress.mit.edu/books/real-world-algorithms).
# 
# ---
# 
# > Panos Louridas<br />
# > Athens University of Economics and Business

# ## The Bellman-Ford(-Moore) Algorithm
# 
# The Bellman-Ford(-Moore) algorithm is an alternative to Dijkstra's algorithm for finding shortest paths in graphs.
# 
# The basic idea is that we find the shortest paths from a start node to other nodes by using 1, 2, $n-1$ links, where $n$ is the number of nodes in the graph.
# 
# So we start with shortest paths that have only one link, then with two links, and so on.

# We will use the following function to read weighted graphs from files containing one line per edge.
# 
# This is familiar, it is the same that we used with Dijkstra's algorithm.

# In[1]:


def read_graph(filename, directed=False):
    graph = {}
    with open(filename) as input_file:
        for line in input_file:
            parts = line.split()
            if len(parts) != 3:
                continue # not a valid line, ignore
            [n1, n2, w] = [ int (x) for x in parts ]
            if n1 not in graph:
                graph[n1] = []
            if n2 not in graph:
                graph[n2] = []
            graph[n1].append((n2, w))
            if not directed:
                graph[n2].append((n1, w))
    return graph


# As a first example we will use the traffic grid graph again.
# 
# <img width="400" src="traffic_grid_graph.png"/>

# We go ahead and read it.

# In[2]:


import pprint

g = read_graph('traffic_grid_graph.txt')
pprint.pprint(g)


# We need a substitute for $\infty$ in the algorithm.
# 
# We will use `sys.maxsize` as a substitute for $\infty$ .
# 
# In Python there is no predefined maximum integer; `sys.maxsize` gives the maximum 32 bit or 64 bit integer, depending on the machine.
# 
# If that is not enough, we could simply assign a larger value.

# In[3]:


import sys 

MAX_INT = sys.maxsize


# The Bellman-Ford(-Moore) algorithm in Python is as follows; note that we just call it `bellman_ford(g, s)`.

# In[4]:


def bellman_ford(g, s):
    nodes = g.keys()
    num_nodes = len(nodes)
    # Initialize array holding path lengths.
    dist = [ MAX_INT ] * num_nodes
    dist[s] = 0
    # Initialize array holding node predecessors.
    pred = [ -1 ] * num_nodes

    # Try using paths of length up to num_nodes
    for i in range(num_nodes - 1):
        # Try all edges: get all nodes...
        for u in nodes: 
            # ...then for each node get all edges.
            for v, w in g[u]:
                # If the path to node v is bigger than the path to
                # node u and then the edge to v, update it.
                # If the distance to u is MAX_INT, it means we
                # have not reached it yet, so we should ignore it.
                if dist[u] != MAX_INT and dist[v] > dist[u] + w:
                    dist[v] = dist[u] + w
                    pred[v] = u
    
    return (pred, dist)


# We can apply it directly on our traffic grid graph:

# In[5]:


pred, dist = bellman_ford(g, 0)
print('pred', pred)
print('dist', dist)


# We can improve on our implementation by noting that at each iteration we do not need to check all edges.
# 
# In particular, at each iteration of the algorithm, instead of checking all edges, we need to check only the edges of the nodes whose estimates we updated at the previous iteration.
# 
# To do that, we will use a First-In First-Out (FIFO) queue, which we'll get from Python's [`deque`](https://docs.python.org/3/library/collections.html#collections.deque).
# 

# Each time we update an edge we will add the target node in the queue.
# 
# We will get the edges to check in each iteration by getting the edges of the first node in the queue each time.
# 
# We will kick-off the process by adding the source node into the queue.
# 
# We'll call the implementation `bellman_ford_qb(g, s)`, where `qb` stands for "queue-based".
# 

# In[6]:


import collections

def bellman_ford_qb(g, s):
    nodes = g.keys()
    num_nodes = len(nodes)
    # Initialize array holding path lengths.
    dist = [ MAX_INT ] * num_nodes
    dist[s] = 0
    # Initialize array holding node predecessors.
    pred = [ -1 ] * num_nodes
    # Initialize queue.
    q = collections.deque()
    # We'll use a list to check whether something
    # is already in the queue, so that we won't
    # add it twice.
    in_queue = [ False ] * num_nodes
    # We start by putting the starting node in the queue.
    in_queue[s] = True
    q.append(s)

    # While the queue is not empty:
    while len(q) != 0:
        u = q.popleft()
        in_queue[u] = False
        # For every edge of the current node, check
        # and update if needed.
        for (v, w) in g[u]:
            if dist[u] != MAX_INT and dist[v] > dist[u] + w:
                dist[v] = dist[u] + w
                pred[v] = u
                # Add the node of the updated path in the queue
                # if necessary.
                if in_queue[v] == False:
                    q.append(v)
                    in_queue[v] = True

    return (pred, dist)


# Here are the results of the algorithm on the traffic grid example:

# In[7]:


pred, dist = bellman_ford_qb(g, 0)
print('pred', pred)
print('dist', dist)


# The Bellman-Ford(-Mooore) algorithm is generally slower than Dijkstra's, but it handle graphs with negative weights.
# 
# Consider again the [negative_weights_graph.txt](negative_weights_graph.txt):
# 
# <img width="150" src="negative_weights_graph.png"/>

# We go ahead and read it:

# In[8]:


g = read_graph('negative_weights_graph.txt', directed=True)
pprint.pprint(g)


# Then if we try `bellman_ford(g, s)` on it, we'll see that we get the right results:

# In[9]:


pred, dist = bellman_ford(g, 0)
print('pred', pred)
print('dist', dist)


# The same with the queue-based version:
#  

# In[10]:


pred, dist = bellman_ford_qb(g, 0)
print('pred', pred)
print('dist', dist)


# What about negative cycles?
# 
# Consider the [negative_cycle_graph.txt](negative_cycle_graph.txt):
# 
# <img width="220" src="negative_cycle_graph.png"/>
# 
# 

# We go ahead and read it:

# In[11]:


g = read_graph('negative_cycle_graph.txt', directed=True)
pprint.pprint(g)


# Then we try bellman_ford(g, s) on it, and see what we get:

# In[12]:


pred, dist = bellman_ford(g, 0)
print('pred', pred)
print('dist', dist)


# The algorithm terminates, but the results do not make much sense.
# 
# For example, the distance to node 1 is -3, but in fact it is $-\infty$, as we can go round and round the negative cycle.
# 
# The real problem is that we did not get any indication that the results are bogus because of negative cycles.

# To fix this problem, we will be adding to the queue a special, *sentinel* node.
# 
# That is a fictitious node that does not exist in the graph.
# 
# It will demarkate each set of neighbours that we add in the queue.
# 
# Each time we meet it in the queue we will know that we have handled the complete set of neighbours of one node.
# 
# This cannot happen more than $n$ times, where $n$ is the number of nodes; so if it does, we have reached a negative cycle.

# In[13]:


def bellman_ford_nc(g, s):
    nodes = g.keys()
    num_nodes = len(nodes)
    # The sentinel is equal to the number of nodes,
    # as the nodes start from 0.
    sentinel = num_nodes
    # Initialize array holding path lengths.
    dist = [ MAX_INT ] * num_nodes
    dist[s] = 0
    # Initialize array holding node predecessors.
    pred = [ -1 ] * num_nodes
    q = collections.deque()
    # We'll use a list to check whether something
    # is already in the queue, so that we won't
    # add it twice.
    in_queue = [ False ] * num_nodes
    in_queue[s] = True
    # We start by putting the starting node and the
    # sentinel in the queue.
    q.append(s)
    q.append(sentinel)
    
    i = 1 # number of iterations
    # Repeat as long as the queue contains more than one
    # element (the sentinel) and we have not handled
    # the neighbours of all nodes.
    while len(q) != 1 and i < num_nodes:
        u = q.popleft()
        # If we have reached the sentinel, update the
        # nodes count and put the sentinel back in the
        # queue.
        if u == sentinel:
            i += 1
            q.append(sentinel)
        else:
            in_queue[u] = False
            # For every edge of the current node, check
            # and update if needed.
            for (v, w) in g[u]:
                if dist[u] != MAX_INT and dist[v] > dist[u] + w:
                    dist[v] = dist[u] + w
                    pred[v] = u
                    # Add the node of the updated path in the queue
                    # if necessary.
                    if in_queue[v] == False:
                        q.append(v)
                        in_queue[v] = True

    return (pred, dist, i < num_nodes)


# So now we are able to handle graphs with negative weights and detect cycles:

# In[14]:


pred, dist, no_negative_cycle = bellman_ford_nc(g, 0)
print('pred', pred)
print('dist', dist)
print('no_negative_cycle', no_negative_cycle)