CG Project 2019

Abstract

This project covers Nimrod Megiddo's paper [1] that presents a generalisation of how parallelism can benefit the design of serial algorithms in computational problems. In many cases, a good parallel algorithm for one problem may turn out to be useful for designing an efficient serial algorithm for another problem.

We start by (visually) exploring Meggido's preliminary example that demonstrates his general method, known as parametric search. The technique is frequently used to solve optimization problems in computational geometry. The method uses a decision algorithm, i.e. an algorithm that checks if a given condition holds or not, to solve an optimization problem, i.e. an algorithm that determines an optimal solution for a problem. By cleverly running the decision problem as a concurrent process, Megiddo shows that parallelism can achieve a lower theoretical bound in some cases.

Next, a more advanced setting is explored in the context of the Minimum Ratio Cycle problem (MRC), where the goal is to find a cycle in a graph with edge-costs and edge-times, such that the ratio of the total cost to the total time of edges on the cycle is minimized. This problem was originally coined by Dantzig et al. [2] in an attempt to solve a ship routing problem.

Finally, we present an attempt at an implementation of the parameterization method for the MRC, inspired on Megiddo's earlier work [3] and conclude with some remarks and outlooks.

My presentation that covers this project in more depth can be found here.

Megiddo's Preliminary Example

Let $ \{f_{1}(\lambda), \dots, f_{n}(\lambda)\} $ be a set of pairwise distinct functions of $ \lambda $ such that $ f_{i}(\lambda) = a_{i} + b_{i}\lambda \ | \ b_{i} > 0 $. Next, let $ F(\lambda) $ denote the median of the set $ \{f_{1}(\lambda), \dots, f_{n}(\lambda) \} $ for every $ \lambda \in \mathbb{R} $. Finally, let $ \lambda_{ij} $ denote the intersection of two distinct functions $ f_{i} $ and $ f_{j} $ in the set $ \{f_{1}(\lambda), \dots, f_{n}(\lambda) \} $, such that $ a_{i} + b_{i}\lambda_{ij} = a_{j} + b_{j}\lambda_{ij} $ with $ i \neq j$.

Function Visualisation

To showcase how functions from the set $ \{f_{1}(\lambda), \dots, f_{n}(\lambda) \} $ behave, consider an interactive visualisation of three functions $ f_{1}(\lambda) = a_{1} + b_{1}\lambda, \ f_{2}(\lambda) = a_{2} + b_{2}\lambda $ and $ \ f_{3}(\lambda) = a_{3} + b_{3}\lambda \ $ together with their corresponding $ F(\lambda) $ evaluated on $ \{f_{1}(\lambda), f_{2}(\lambda), f_{3}(\lambda) \} $. Notice that $ b_{i} $ is constrained to be positive and that the intersections marked in black $ \lambda_{12}, \ \lambda_{13} $ and $ \lambda_{23} $ are critical points for $ F(\lambda) $, because the evaluation of $ F(\lambda) $ can only change at these points.

$ a_{1} = $ 2

-10 10

$ a_{2} = $ 10

-10 10

$ a_{3} = $ -10

-10 10

$ b_{1} = $ 0

0 20

$ b_{2} = $ 4

0 20

$ b_{3} = $ 9

0 20

It's easy to see now that $ F(\lambda) $ is a piecewise linear function that lies on the line of the function that evaluates to the median for a particular $ \lambda $. As the $ b_{i} $'s are positive, this means that each segment of the function is monotone-increasing.

In piesewise linear functions, the points where the slope of the function changes are known as breakpoints. Given that we know that a function of $ n $ lines has at most $ \frac{n^{2} - n }{2} $ intersections, the upper bound for breakpoints is $ O(n^{2}) $. It is straightforward to see now that if we are given $ \lambda $, we can evaluate $ F(\lambda) $ using a linear-time selection algorithm [4] such that $ F(\lambda) $ can be evaluated with $ O(n) $ comparisons when all $ f_{i}(\lambda) $'s have been computed.

A first serial attempt

If we were to try to find a $ \lambda $ such that $ F(\lambda) = 0 $, the visualisation above might suggest one possible way to do that:

Identify the intersection points $ \lambda_{ij} $ where $ a_{i} + b_{i}\lambda_{ij} = a_{j} + b_{j}\lambda_{ij}$ and $ i \neq j $. We know that every breakpoint of $ F $ is equal to one of the $ \lambda_{ij}$'s.
Search the set of identified $ \lambda_{ij} $'s for two values $ \lambda^{1} $ and $ \lambda^{2} $, such that there are no other $ \lambda_{ij} $'s in the open interval $ (\lambda^{1}, 0) \cup (0, \lambda^{2}) $.

The set of intersection points can be sorted and binary-searched for an interval $[\lambda^{i}, \lambda^{i+1}]$ that requires $ O(log(n)) \ F$ evaluations in which the cardinalities of the searched sets are $ n, \frac{n}{2}, \frac{n}{4}, \dots \ $. As $ F(\lambda) $ is linear over the interval $[\lambda^{1}, \lambda^{2}]$, the intersection with the origin can be found in $ O(1) $ once $ \lambda^{1} $ and $ \lambda^{2} $ are known. The total running time of this algorithm is $ O(n^{2} + n \ log(n)) + log(n)) $ = $ O(n^{2})$, where the dominating factor is the evaluation of all the intersection points $ \lambda_{ij} $,

A second serial attempt

An alternative approach suggested by Megiddo [3] yields the same serial running time of $ O(n^{2}) $ but is shown to benefit from parallel computation. Remember that we are looking for for two values $ \lambda^{1} $ and $ \lambda^{2} $, such that there are no other $ \lambda_{ij} $'s in the open interval $ (\lambda^{1}, 0) \cup (0, \lambda^{2}) $. With this interval, we can solve for the origin in a single step.

Let $ \lambda^{*} $ be the unknown solution to $ F(\lambda^{*}) = 0 $. In order to evaluate $ F $ by a linear-time median algorithm, a value $ \lambda $ needs to be known to compare each $ f_{i}(\lambda) $ in the set $ \{f_{1}(\lambda), \dots, f_{n}(\lambda) \} $. Nevertheless, we know that only intersection points $ \lambda_{ij} \ | \ i \neq j $ really matter for the comparison of the median algorithm as they are critical values. Calculating the intersection point $ \lambda_{ij} $ itself does not require an input value $ \lambda $ and we know that for any input value $ \lambda $: if $ \lambda \geq \lambda_{ij} $, then $ f_{i} (\lambda) \geq f_{j}(\lambda) $ whereas for $ \lambda \leq \lambda_{ij} $ we know that $ f_{i} (\lambda) \leq f_{j}(\lambda) $, or vice versa.

So when we consider two functions $ f_{i} (\lambda) $ and $ f_{j} (\lambda) $, all we really need is the intersection point $ \lambda_{ij} $ to compare $ F(\lambda_{ij}) $ and $0$ to solve for $ F(\lambda^{*}) = 0 $. As each function is monotone increasing, we know that if the sign of $ F(\lambda_{ij}) $ is negative, then $ \lambda^{*} \geq \lambda_{ij} $ and if the sign is positive, then $ \lambda^{*} \leq \lambda_{ij} $.

The remark above suggests a parallel algorithm, as each intersection point $ \lambda_{ij} $ of two functions in the set $ \{f_{1}, \dots, f_{n} \} $ can be computed and compared independently with a shared variable that keeps track of the smallest interval bounds so far. If an intersection point $ \lambda_{ij} $ lies between the current bounds $ [\lambda_{min}, \lambda_{max}]$ then we can update that interval to be either $ [\lambda_{ij}, \lambda_{max}]$ or $ [\lambda_{min}, \lambda_{ij}]$ depending on the sign of $ F(\lambda_{ij}) $. If an intersection point $ \lambda_{ij} $ lies outside the current bounds $ [\lambda_{min}, \lambda_{max}]$ then the result of the comparison is uniform over the interval and the interval can remain as is.

In serial time, the total running time remains $ O(n^{2}) $ as we have to gather all intersections of the set $ \{f_{1}, \dots, f_{n} \} $ in $ O(n^{2}) $ and out of all intersections, only $ O(n) $ critical points may require a $ O(n) $ median finding to evaluate $ F $, leaving us with $ O(n^{2} + n^{2}) = O(n^{2}) $.

Speeding Up with Parallelism

Before making any claims about the parallel running time, it's important to note that Megiddo does not account for overhead resulting from memory conflicts, memory syncronization or deadlocks that are typically associated with parallelism. Rather, Megiddo simulates parallelism serially by running each independent process after one another. This reduces the complexity of his analysis and results.

Let $ P $ denote the number of processors and let $ T(n, P) $ denote the number of comparison steps for sorting $ n $ values on a machine with $ P $ processors.

As already mentioned, the intersection $ \lambda_{ij} $ of two functions $ f_{i} $ and $ f_{j} $ from the set $ \{f_{1}, \dots, f_{n} \} $ can be computed in $ O(1) $ such that one step in the simulated multiprocessor yields $ P $ critical values. As suggested in the explanation of the previous algorithm, these comparisons can happen independently of each other.

Without loss of generality, we can say that there exists a labeling such that the critical values resulting from a multiprocessor step can be written as $ \lambda_{1} \leq \lambda_{2} \leq \dots \leq \lambda_{P}.$ Megiddo does not assume that these values are sorted in his running time analysis, but with a parallel sorting algorithm like Valiant's [5] or Preparata's [6], a sorting time of $ T(n, P) = O(log(n)) $ or $ T(n, P) = O(log(n)log(log(n)) $ can be achieved when $ P = n $ and $ P = log(n) $ respectively.

With the set of $ P $ critical values sorted, an interval $ [\lambda_{i}, \lambda_{i+1}] \ | \ i \in [0, P] $ can be identified with a $ O(log(P)) $ binary search, using a linear-time median algorithm to pivot, such that $ \lambda_{i} \leq \lambda^{*} \leq \lambda_{i+1} $. In each step, we take an intersection between the found interval and the interval of the previous step to get the smallest possible interval that satisfies $ \lambda_{i} \leq \lambda^{*} \leq \lambda_{i+1} $. This repeats $ T(n, P) $ times until an interval is found that does not intersect with the previous interval and the median function $ F $ crosses the zero level over the interval. Here, the equation can be solved in $ O(1) $.

The running time of the algorithm presented above yields $ O(nlog(n)^{2}) $ for Preparata's scheme and $ O(nlog(n)^{2}log(log(n))) $ for Valiant's. Both bounds being superior to the previous $ O(n^{2}) $.

Minimum Ratio Cycle problem (MRC)

The MRC problem was coined by Dantzig et al. [2] in the context of a ship routing problem posed by the American Office of Naval Research. Practically, the problem concerns a ship owner that wants to maximize his mean daily profit over time while making a round trip through multiple ports. The solution to the problem gives exactly the path that the skipper should take to maximize his profit. As these kinds of problems are usually posed in a minimization setting, the problem can be transformed to minimize the negative profit over time ratio, which is exactly the MRC problem.

Representation

More abstractly, each vertex $ i $ within a directed graph $ \mathcal{G} = (V, E)$, with $V$ a set of vertices and $E$ a set of edges, has an associated cost $ c_{ij} $ and travel time $ t_{ij} $ to reach a vertex $ j $. Self-loops are not allowed in $ \mathcal{G} $.

The problem can be represented as an adjacency matrix $ A \in \{0, 1\}^{|V| \times |V|} $ determining the directed graph, a cost matrix $ C \in \mathbb{R}^{|V| \times |V|} $ and a time matrix $ T \in \mathbb{R}^{|V| \times |V|} $. If there is a directed edge from vertex $ i $ to $ j $, then $ A_{ij} = 1 $, if not then $ A_{ij} = 0 $. As there are no self-loops in $ \mathcal{G} $, the diagonal elements of the matrix are always $ 0 $. An interactive example of the setup is given below.

$$ A = \begin{bmatrix} 0 & 1 & 0 & 0 & 0 \\ 1 & 0 & 1 & 0 & 1 \\ 0 & 0 & 0 & 1 & 0 \\ 1 & 0 & 1 & 0 & 0 \\ 1 & 1 & 0 & 0 & 0 \\ \end{bmatrix}\ C = \begin{bmatrix} 0 & 8 & 0 & 0 & 0 \\ 3 & 0 & 9 & 0 & 4 \\ 0 & 0 & 0 & 7 & 0 \\ 2 & 0 & 1 & 0 & 0 \\ 6 & 11 & 0 & 0 & 0 \\ \end{bmatrix}\ T = \begin{bmatrix} 0 & 3 & 0 & 0 & 0 \\ 1 & 0 & 4 & 0 & 1 \\ 0 & 0 & 0 & 2 & 0 \\ 3 & 0 & 1 & 0 & 0 \\ 2 & 9 & 0 & 0 & 0 \\ \end{bmatrix}\ $$

Megiddo's Theorem

To continue, we need to get this problem in a formulation that expresses the minimization. Numerous combinatoral optimization problems can be formulated as linear minimization problems subject to certain constraints. Megiddo highlights two classes of problems [3] and proves an important relationship between them.

$$ \textbf{Problem A. } \ \text{Minimize } \ c_{1}x_{1} + \dots + c_{n}x_{n} \ \text{ subject to } \ x = (x_{1}, \dots, x_{n}) \in D $$ $$ \textbf{Problem B. } \ \text{Minimize } \ \frac{a_{0} + a_{1}x_{1} + \dots + a_{n}x_{n}}{b_{0} + b_{1}x_{1} + \dots + b_{n}x_{n}} \ \text{ subject to } \ x = (x_{1}, \dots, x_{n}) \in D $$ $$\text{With } \ D \ \text{ a set of conditions or constraints to which } \ x = (x_{1}, \dots, x_{n}) \\ \text{ must adhere in order to be a valid solution.} $$

Examples satisfying the definition of problem $ A $ are the shortest (simple) path, the traveling salesman, the Chinese postman, the minimum spanning tree and various other scheduling problems. Notice that problem $ B $ is a generalization of problem $ A $ that also includes MRC as a valid problem instance.

Before getting to any concrete algorithm for MRC, it is beneficial to look at a more general setting. The goal is to relate problem $ A $ with problem $ B $, such that we can use algorithms for problem $ A $ to solve problem instances of problem $ B $. Megiddo's main result goes as follows.

$\textbf{Theorem. } $ If problem $A$ is solvable within $O(p(n))$ comparisons and $ O(q(n)) $ additions, then problem $ B $ is solvable in time $ O(p(n)(q(n) + p(n)))$.

Megiddo starts by showcasing a standard trick to solve ratio minimization problems. Given a problem $B $ with values for the coefficients, pick a fixed number $t \in \mathbb{R} $ and solve problem $A $ for $ \{ c_{i}(t) = a_{i} - tb_{i} \ | \ i = 1, \dots, n \} $ under the same constraints $D $.

Now suppose that $ v $ is an optimal value of problem $ A $ such that $ v = \min \big( (a_{1} - tb_{1})x_{1} + \dots + (a_{n} - tb_{n})x_{n} \big) $ with $ x \in D $. If $ v $ can be written as $ v= tb_{0} - a_{0}$, then $ t $ must be the optimal value for problem $ B $, such that $ t = \min \big( \frac{a_{0}}{b_{0}} + \frac{a_{1}}{b_{1}}x_{1} + \dots + \frac{a_{n}}{b_{n}}x_{n} \big) $. On the other hand, the solution $ x = (x_{1}, \dots, x_{n}) \in D $ of problem $ A $ is optimal with respect to $ t $ because of the assumptions.

In the case that $ v < tb_{0} - a_{0}$, a smaller $t$ should be tested and if $ v > tb_{0} - a_{0}$, a bigger $t$ should be tested. Essentially, Megiddo proofs his theorem by parameterizing problem $ A $ with a linear function $ c_{i}(t) $ for which the optimal $ t^{*} $ needs to be found. Finding $ t^{*} $ is done by comparing each minimum value $v(t) = \min \big( (a_{1} - tb_{1})x_{1} + \dots + (a_{n} - tb_{n})x_{n} \big) $ of problem $A$ and gradually adjusting the interval $[e, f]$ that contains $ t^{*} $. This parameterization and convergence to an interval in allow Megiddo to relate problem $ A $ with problem $ B $ and derive his theorem.

Preparing the MRC Implementation

We are now ready to understand a practical example that uses this result to find a minimum ratio cycle. We can see that a matrix $ X \in \{0, 1\}^{|V| \times |V|} $ could encode the valid paths that are available to take, to determine what cost-time ratios are accounted for. Megiddo suggests that $ D $ in this case is the set of all possible zero-one matrices $ X $ that encode a simple directed cycle. An example of an $ X \in D $ based on the problem defined by adjacency matrix $ A $ is given below. The numerator and denominator of problem $ B $ correspond with the previously discussed $ C $ and $ T $ matrix respectively.

$$ A = \begin{bmatrix} 0 & 1 & 0 & 0 & 0 \\ 1 & 0 & 1 & 0 & 1 \\ 0 & 0 & 0 & 1 & 0 \\ 1 & 0 & 1 & 0 & 0 \\ 1 & 1 & 0 & 0 & 0 \\ \end{bmatrix} \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ X = \begin{bmatrix} 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 \\ 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 \\ \end{bmatrix} $$

The first step towards a MRC solution is an algorithm that finds the shortest simple cycle, i.e. no repetition of vertices and edges, in a directed graph where each edge has an associated distance $ c_{ij}$ that may take a negative value. This setting can be formulated as a minimization problem for the sum $ \sum{c_{ij}x_{ij}} $ such that it becomes a valid instance of the problem $ A $ class that we discussed earlier. Additionally, allowance of negative distance values requires negative cycle detection to stop the algorithm from running forever.

The Floyd-Warshall algorithm (FW), originally solving the all-pairs-shortest-path problem (APSP), can be turned into a shortest simple cycle finder and can also detect negative cycles by running it again, and checking if the values remain the same for each element in the APSP matrix. If they do not, then there is a negative cycle in the graph, otherwise, there are no negative cycles. More background on using FW for negative cycle detection can be found in Hougard's paper [7].

Before we can apply Floyd-Warshall as a (negative) cycle detector, we need to encode the cost-time ratios in our adjacency matrix. Aditionally, we encode the absence of a path between two nodes with $ \infty $. Note that the difference between Floyd-Warshall as a solution to the APSP problem and the cycle detection problem is quite subtle. To apply it as a APSP, the diagonal elements are encoded with $ 0 $. To use FW as a simple cycle detector, the elements on the diagonal are encoded with $ \infty $.

$$ X^{'} = \begin{bmatrix} \infty & \frac{C_{0, 1}}{T_{0, 1}} & \infty & \infty & \infty \\ \frac{C_{1, 0}}{T_{1, 0}} & \infty & \frac{C_{1, 2}}{T_{1, 2}} & \infty & \frac{C_{1, 4}}{T_{1, 4}} \\ \infty & \infty & \infty & \frac{C_{2, 3}}{T_{2, 3}} & \infty \\ \frac{C_{3, 0}}{T_{3, 0}} & \infty & \frac{C_{3, 2}}{T_{3, 2}} & \infty & \infty \\ \frac{C_{4, 0}}{T_{4, 0}} & \frac{C_{4, 1}}{T_{4, 1}} & \infty & \infty & \infty \\ \end{bmatrix}\ = \begin{bmatrix} \infty & \frac{\Large 8}{\Large 3} & \infty & \infty & \infty \\ \frac{\Large 3}{\Large 1} & \infty & \frac{\Large 9}{\Large 4} & \infty & \frac{\Large 4}{\Large 1} \\ \infty & \infty & \infty & \frac{\Large 7}{\Large 2} & \infty \\ \frac{\Large 2}{\Large 3} & \infty & \frac{\Large 1}{\Large 1} & \infty & \infty \\ \frac{\Large 6}{\Large 2} & \frac{\Large 11}{\Large 9} & \infty & \infty & \infty \\ \end{bmatrix} $$

Floyd-Warshall Algorithm as APSP

The basic idea of the Floyd-Warshall algorithm comes from the observation that an edge $ e_{ij} $ is the shortest path from vertex $i$ to $j$ if and only if there is no other route $ \{e_{ik},e_{kj} \} $ such that $ \delta(e_{ij}) > \delta(e_{ik}) + \delta(e_{jk}) $ with $\delta(x) $ a distance function of choice. The input $ m $ of the FW algorithm corresponds with the matrix $ X^{'} $ from above. Note that the resulting matrix only gives the shortest distance between two nodes and that reconstructing a shortest path requires an additional 2-dimensional matrix that keeps track of the shortest subsequent node for a given target node.


      function basic_floyd_warshall(m) {
        for(var k = 0; k < n; k++) {
          for(var i = 0; i < n; i++) {
            for(var j = 0; j < n; j++) {
              if(m[i][k] + m[k][j] < m[i][j]) {
                m[i][j] = m[i][k] + m[k][j]; // Detour is shorter
              }
            }
          }
        }
      }

Now, let $ u_{ij}^{(m)}(t)$ be the length of a shortest simple path from $i $ to $j $, only using nodes from the set $ \{1 \dots m-1 \} \cup \{i, j\} $ and using a distance function $ c_{ij}(t) = a_{ij} - tb_{ij} $, the parameterization that is a core element for algorithm $ B $. In our implementation, this is represented as a function $ f(t) $ that is returned from a call to the function $parameterize(A, C, T)$, where $A$ is an adjacency matrix, $C$ is a cost matrix and $T$ a time matrix for the given problem.

With the $t$-parameterized function, the first step is to find a $ t' \in [e,f] $ such that $ u_{ij}^{(m)}(t) = u_{im}^{(m)}(t) + u_{mj}^{(m)}(t)$. An interactive example for the parameterization of the $ A $ algorithm for a fully connected $ K_{4} $ graph, where the minimum ratio cycle is colored in blue, and its corresponding $ u_{ij}^{(m)}(t),\ u_{im}^{(m)}(t),\ u_{mj}^{(m)}(t) $ and $ u_{im}^{(m)}(t) + u_{mj}^{(m)}(t) $ functions is given below.

$i = $ 2

0 3

$ j = $ 2

0 3

$ m = $ 1

0 3

In the visualisation above, the red line represents $ u_{ij}^{(m)}(t) $, the blue lines represent $ u_{im}^{(m)}(t), u_{mj}^{(m)}(t) $ and the green line represents $ u_{im}^{(m)}(t) + u_{mj}^{(m)}(t) $. When $ t = 0$, a clear nod can be observed in all functions due to the sign change that transforms the distance functions from $ c_{ij}(t) = a_{ij} - tb_{ij} \ $ to $\ c_{ij}(t) = a_{ij} + tb_{ij} $. This causes the shortest path with respect to $c_{ij}(t) $, i.e. the minimum path distance value, to jump to the smallest, most negative path. As Megiddo does not cover how this should be handled in his paper, we assume this to be correct. Another obvious observation is all functions are linear over $ (-\infty, 0) $ and $ (0, \infty) $ and that the value $ t' $ for which $ u_{ij}^{(m)}(t') = u_{im}^{(m)}(t') + u_{mj}^{(m)}(t')$ is the intersection of the red and yellow line, marked with a black point.

In a similar way as the first example, we use the function $ c_{ij}(t) $ to guide the search and refine the solution interval.

Intialise the search interval $[e, f]$ to $[-\infty, \infty]$ and set $ i = j = m = 0 $ where $ i $ and $ j $ are constrained such that $ 0 \leq i \leq j \leq n $, with $n$ the total number of vertices in the graph.
Solve for $ t $ in $ u_{ij}^{(m)}(t) = u_{im}^{(m)}(t) + u_{mj}^{(m)}(t)$:

If $ u_{ij}^{(m)}(t) - u_{im}^{(m)}(t) - u_{mj}^{(m)}(t) = 0$ has a unique solution $ t' \in [e, f]$, then go to the next step.
If $ u_{ij}^{(m)}(t) - u_{im}^{(m)}(t) - u_{mj}^{(m)}(t) = 0$ has multiple solutions , go to step 4.

Check if the graph with distances $ c_{kl}(t') = a_{kl} - t'b_{kl}$ has negative cycles and / or zero cycles.

If there is a zero cycle and no negative cycle, then terminate with the zero cycle as minimum ratio cycle.
If there is a negative cycle, then update the upperbound $ f $ of the search interval to $ t' $.
If all cycles are positive, then update the lowerbound $ e $ of the search interval to $ t' $.

Set $ u_{ij}^{(m+1)}(t) = \min( u_{ij}^{(m)}(t), u_{im}^{(m)}(t) + u_{mj}^{(m)}(t)) $
Update parameters such that:

If $ j < n $, update $ j \leftarrow j + 1 $ and go back to step 1.
If $ j = n $ and $ i < n $, update $ i \leftarrow i + 1 $ and go back to step 1.
If $ i = j = n $ and $ m < n $, update $ i \leftarrow 1, j \leftarrow 1, m \leftarrow m + 1 $ and go back to step 1.
If $ i = j = n $ and $ m = n $, then go to step 5.

Find a $ k$ such that $ u_{kk}^{n+1}(f) < 0$. The minimum ratio cycle is the shortest simple cycle with respect to the distance metric $ c_{ij}(f)$ that contains $ k $.

Our attempt at the implementation of a parametric MRC algorithm can be found here and is made available as a global function that can be run in the browser console, as below. $ A_4 $ represents the adjacency matrix, $ C_4 $ the cost matrix and $ T_4 $ the time matrix.


    // Open browser console and run.
    mrc(A_4, C_4, T_4);

We note however that this implementation is not fully waterproof due to unclarity in the algorithm specification. For example, we do not know what should be done in the case that the upper bound of the search interval $ f $, remains $ \infty $. This can be the case when no single unique solution $ t' $ gives rise to a negative cycle during the execution of the MRC algorithm in step 2.

As the final step (5) consists of evaluating $u_{kk}^{(n+1)}(f)$, it is crucial for $ f $ to evaluate to an evaluable value, different from $ \infty $. Nevertheless, we can still proceed to discuss the parametric properties of the algorithm.

MRC Parameterization

We are now ready to see that Megiddo's preliminary example is in fact very similar to the MRC algorithm presented above. The idea of the $ B $ algorithm, i.e. the minimum ratio cycle, is solved by parametrically employing algorithm $ A $, i.e. the shortest simple cycle, over an interval $ [e, f] $.

Whenever we try to solve for $ t' $ in step 1, a comparison between two linear functions $ c_{1}(t) = u_{ij}^{(m)} $ and $ c_{2}(t) = u_{im}^{(m)} + u_{mj}^{(m)} $ has to be made to consider an update of the interval $ [e, f] $. When the algorithm gets to step 3, the result of this comparison could be a piecewise linear function $ c_{3}(t) = \min{(c_{1}(t), c_{2}(t))} $, with one breakpoint at most. This breakpoint is then tested by running the $ A $ algorithm on it, i.e. the shortest simple cycle.

In the same fashion as the preliminary example, postponing the evaluation of each critical point by algorithm $ A $, results in a set of critical points $ e = t_{0} < t_{1} < \dots < t_{k} = f $, where each subinterval corresponds with a piecewise linear function over $ [t_{i-1}, t_{i}] $. Provided that the critical values $ t_{0}, t_{1}, \dots t_{k}$ are sorted, finding a subinterval $ [t_{i-1}, t_{i}] $ can be done by employing a binary search as in the first example.

Conclusion and Outlook

In this work we have studied two in-depth examples of parametric search. The first being a more artificial example, and the second a formulation of a minimization problem. While both problems seem different at first glance, the parametric search formulation of both problems shows to be surprisingly similar. The general nature and applicability of parametric search contribute to its success as an optimization technique for a wide range of problems [1] in combinatorial optimization and computational geometry.

In a nutshell [8], parametric search uses a decision problem $ P(\lambda) $ that is monotone in parameter $ \lambda \in \mathbb{R} $, meaning that if $ P(\lambda) = v$ then $ P(\lambda') = v $ for all $ \lambda' < \lambda $. The goal is then to solve the optimization problem of finding an optimal value $ \lambda^{*} $ that is the largest value satisfying $ P(\lambda^{*}) = v $.

Despite the success, algorithms based on parametric search are considered difficult to implement, making it primarily of theoretical interest. We found that the formulations of problems could benefit from more visual representations of the setting, yielding faster intuitions of how parametric search can be used in different settings. We would also incorporate more practical implementations, as the gap between theoretical algorithms and their implementation is often bigger than one might expect. The benefit of using JavaScript in this case is that imlementations can easily be run in the browser, without installations of any kind. With technologies like web workers, even parallelism could be incorporated in the algorithms, making it easier for other researchers to quickly explore and test these algorithms.

References

[1]

Megiddo, N. (1981, October). Applying parallel computation algorithms in the design of serial algorithms. In 22nd Annual Symposium on Foundations of Computer Science (sfcs 1981) (pp. 399-408). IEEE.

[2]

Dantzig, G. B., Blattner, W. and Rao, M. R. (1967). Finding a Cycle in a Graph with Minimum Cost to Time Ratio with Application to a Ship Routing Problem. In Theory of Graphs, P. Rosenstiehl, ed. Dunod, Paris, and Gordon and Breach, New York. 77-84.

[3]

Megiddo, N. (1978, May). Combinatorial optimization with rational objective functions. In Proceedings of the tenth annual ACM symposium on Theory of computing (pp. 1-12). ACM.

[4]

Blum, M., Floyd, R. W., Pratt, V. R., Rivest, R. L., & Tarjan, R. E. (1973). Time bounds for selection. J. Comput. Syst. Sci., 7(4), 448-461.

[5]

Valiant, L.G. Parallelism in comparison problems SIAM J. Comput 4 (1975), 348-355.

[6]

Preparata, F.P. New parallel-sorting schemes. IEEE Trans Comput. C-27 (1978), 669-673.

[7]

Hougardy, Stefan. (2010). The Floyd-Warshall algorithm on graphs with negative cycles. Inf. Process. Lett.. 110. 279-281. 10.1016/j.ipl.2010.02.001.

[8]

Ilinkin, I. (2013, June). Parametric search visualization. In Proceedings of the twenty-ninth annual symposium on Computational geometry (pp. 343-344). ACM.