List Of Pareto Optimal Solutions of a Biobjective Shortest Path Problem

Many applications in practice involve the search for a shortest path in a network by optimizing two conflicting objective functions. Such problems often are referred to as biobjective optimization problems. Their goal is to find special optimal paths that are nondominated and are also known in the specialized literature as to as Pareto optimal. While most of the existing methods aim to find the minimum complete set of Pareto optimal paths, we propose an approach that is able to generate a list of all Pareto optimal solutions in a given network. The described method solves the biobjective optimization problem in the case in which the first objective function is a linear (MINSUM), while the second objective function is from the “bottleneck” type (MAXMIN). The presented approach is based on two modifications of the Dijkstra's shortest path algorithm that solve the MINSUM and the MAXMIN problems respectively. We prove the correctness and the computational complexity of the presented algorithms. Also, we provide detailed numerical examples that illustrate their execution.


I. INTRODUCTION
P ARETO optimal solutions of biobjective (bicriterion) op- timization problems are a subject of extensive research in combinatorial optimization and operation research disciplines, and in particular the biobjective shortest path problems [1].These special type of shortest path problems can arise in numerous applications in practice including transportation problems, computer networking, robot motion, and many others.
The search of biobjective Pareto optimal solutions in a shortest path problem is an optimization problem that is a subject of two objective functions, lets say f and g.The Pareto optimal paths (also called nondominated) are a set of paths, such that for any path α in it, it is impossible to improve either f or g criterion, without getting worse the other [2].
By determining the objective functions f and g, different types of bicriterion path problems can be defined.The first notable work that examines these types of problems is [3] in which Hansen defines ten types of bicriterion path problems, and also introduces their abbreviations.In particular, the MINSUM-MAXMIN problem is solved with an algorithm with a polynomial complexity O(m 2 log n), where n is the number of vertices of the network, and m is the number of directed edges.In the MINSUM-MAXMIN problem the first objective function is a liner one, while the second objective function is from the bottleneck type.
In the literature actually there are quite few works that focus on the solution of the MINSUM-MAXMIN problem.One of them is [4] which proposes an extension of the Martin's algorithm [5] for a multiobjective shortest path problem with a MAXMIN objective function.Most of the methods that can be found in the literature focus on the combination of two linear functions, for example [2], [6], [7], [8] solve the MINSUM-MINSUM bicriteria path problem.Other works focus on MINMAX-MINSUM problem (see [9] and [10]), however in their case the authors do not consider the Pareto optimality, rather they aim to define a singe objective function by combining the MINMAX and MINSUM criteria.
Another subject that is rarely considered in the literature is the calculation of all Pareto optimal paths from a source vertex v 0 and a destination (target, terminal) vertex v t .Most of the existing algorithms aim to find the minimal complete set of Pareto optimal paths, which means that from each equivalent set of Pareto optimal paths a single path is discovered (see [3]).In [11] the authors look for a set of alternative Pareto optimal paths in a method that solves a concrete practical problem for routing of Hazardous materials and propose shortest path algorithm on a network with two criteria: one that corresponds to road length, and the other that corresponds to a risk measure.
The work [2] is one of the few in the literature that pays special attention to the calculation of all Pareto optimal paths.The authors propose two algorithms depending whether the paths may contain or may not contain loops, both of them based on k shortest paths algorithms in graphs.However, as mentioned above, in this work the two objective functions are linear, and the algorithms do not cover the case in which one of the functions is a bottleneck function.
The exact methods that are present in the literature, are generally classified into labeling and ranking paths algorithms [12].Labeling algorithms can be split into two categories: label setting [3], [5], [8]; and label correcting [6], [7], [11].In the category ranking paths we can classify methods that are based on the k shortest paths algorithms, for example [2].
The other major branch of methods are based on heuristics approaches.For example, in [13] the authors propose a so-lution to the biobjective shortest path problem that is based on a genetic algorithm for which the authors report to find the Pareto optimal set in 77% of the instances.Also, the probabilistic technique ant colony optimization (ACO) and its modifications are adopted in the solution of various complex combinatorial problems (see for example [14]).Heuristic methods often are adopted in various practical problems, like electric vehicle shortest path problem [15].
In this paper we propose an exact method for calculation of all Pareto optimal paths for the MINSUM-MAXMIN problem in a network.We define two helper problems, MINSUM list and MAXMIN list, and we provide two algorithms that solve them, which are based on generalization of the Dijkstra's algorithm [16].We use the solutions of the two helper problems to formulate the method for general problem solution.The correctness of all algorithms is proved, and their computational complexity is shown.Also, we illustrate the algorithms with detailed examples that show their execution.
The paper is organized as follows.In Sec.II we introduce the notations and problems formulation.In Sec.III we describe the two algorithms that solve the two helper problems.In Sec.IV we show how the Pareto optimal solutions list is constructed based on the solution of the two helper problems.Finally, Sec.V contains conclusions and discussions.

A. Notations
number of vertices and m = |E| number of directed edges.Without loss of generality we will assume that V = {1, 2, . . ., n} and E ⊆ V 2 .
We define the following two functions on the set of edges of the digraph, f : E → R + and g : E → R + .The function f assigns to each edge (i, j) ∈ E the positive number f (i, j), which we call the length of the edge e = (i, j).The function g assigns to each edge (i, j) ∈ E: For convenience, we will call the value g(i, j) the capacity of the edge e = (i, j).
The digraph G together with the functions f and g defines the network G = (V, E, f, g) (see [17]).The network is represented by the adjacency list of the outgoing neighbors [18] that is augmented with the length and capacity of the edges in the following way: where Adj(i) = {(j, f (i, j), g(i, j)) : (i, j) ∈ E}, ∀i ∈ V .In this way, if q = Adj(i, k) for some i ∈ V , and a positive integer k, q(1) will denote the k-th outgoing neighbor of the vertex i, q(2) will denote the length of the edge (i, q(1)), and q(3) will denote the capacity of the edge.To denote the adjacency list of outgoing neighbors of a given network G, we will use the notation G.Adj.
Path in the network G is the finite sequence of the type where v j ∈ V, ∀j ∈ {0, 1, . . ., t} are distinct vertices, and e i is an edge with starting vertex v (i−1) and ending vertex v i , that belongs to E for all i ∈ {1, 2, . . ., t}.The path consists of (t + 1) vertices and t edges, the vertex v 0 is the source of the path, and the vertex v t is the destination of the path.The path (2) with a source v 0 and a destination v t connects v 0 with v t and is called a (v 1 , v t )-path, which we will denote with an ordered sequence α of vertices: For each path α = (v 0 , v 1 , . . ., v t ) we define two functions: We will call the number x(α) the length of the path, and the number y(α) the capacity of the path α.Both functions x(α) and y(α) define the objective functions of the problems that we will discuss.We denote all (1, j)-paths with the shorter W j .Then, we will call the number a distance between the vertex 1 and the vertex j.Also, we will call the number the capacity of the vertex j.
For each W j the term Pareto optimal path is defined as follows.
Definition 1.We call the path α ∈ W j Pareto optimal when there does not exist another path β ∈ W j , for which any of the following two conditions is fulfilled: We say that α and β are equivalent (α ∼ β), when x(α) = x(β) and y(α) = y(β).
Besides that, we will denote the distance from vertex 1 to any vertex v with r(v).

B. Problems formulation
Based on the above definitions, we formulate the main problem considered in this paper: Problem 1 (List of Pareto optimal solutions).Compute a list of all Pareto optimal solutions for W n .
To solve the List of Pareto optimal solutions problem, we will use the solutions of the following two helper problems.
The solution of the first helper problem requires the definition of a function minsum(G.Adj) that computes a list of all shortest paths in the network.
Problem 2 (MINSUM list).Compute a list of all (1, n)-paths with minimal length, given by: The solution of the second helper problem requires the definition of a function maxmin(G.Adj), and it is a version of the first helper problem that computes a list of all maximum capacity paths in the network.
Problem 3 (MAXMIN list).Compute a list of all (1, n)-paths with maximal capacity, given by:

III. SOLUTION OF THE TWO HELPER PROBLEMS
To solve the two helper problems we propose two modifications of the Dijkstra's algorithm [16], in which the results hold as well in the case in which the source vertex is selected i 0 ̸ = 1.In both modifications we assume that the network G = (V, E, f, g) is defined using the adjacency list of the outgoing neighbors.
In the computer program implementation of the modified versions of Dijkstra's algorithm we apply the Fibonacci heap data structure [19] for all priority queue operations.Even though the relative complexity of its implementation, this advanced data structure introduces a significant speedup of the algorithm to O(n log n + m), which is proved based on amortized analysis [18].

A. List of all shortest paths
We will solve the problem of computing of a list of all (1, n)shortest paths by finding the subnetwork G = (V, E, f, g) of the shortest paths.Definition 2. We will say that G = (V, E, f, g) is a subnetwork of the shortest paths in the network G = (V, E, f, g), if the following two properties hold: The solution of the MINSUM list helper problem is given by the definition of a function minsum(G.Adj) (Alg.2), which for a given network G calculates the adjacency list of the outgoing neighbors of the shortest paths subnetwork G.The implementation of the minsum(G.Adj) function is based on a modification of the Dijkstra's algorithm [16], as follows.
The algorithm splits the set of network vertices in into subsets.The first subset V 0 denotes the vertices that are not yet traversed by the algorithm.The second subset U = V \ V 0 Algorithm 1 Function relaxS(u, v, d, pr) Input: vertices u, v, and vectors d, pr Output: stores the traversed vertices.The procedure that traverses the network guarantees that for each vertex v ∈ V 0 and each vertex u ∈ U .Initially, V 0 = V and U = ∅.In each of n consecutive iterations of execution the function minsum(G.Adj) a selected vertex is transferred from V 0 into U .The algorithm uses two vectors d and pr, both of them with n components.Initially, d(1) = 0, and all other components of d are ∞.The initial values of pr are equal to the empty set ∅.After the completion of the algorithm, d will store the distances from the source vertex to each of the other vertices in the network, in other words d(j) = r(j), ∀j ∈ V , and pr will store the adjacency list of the ingoing neighbors of the digraph (V, E).The last step of the minsum(G.Adj) function composes the adjacency list of the outgoing neighbors of the shortest paths subnetwork G.
In the proposed variant of the Dijkstra's algorithm that solves Prob.2, the function that implements the relaxation procedure is modified.We define the function relaxS(u, v, d, pr) that performs relaxation of the edge (u, v) by changing the current state of d and pr, as it is given in Alg. 1.
We use two more helper functions in our modified Dijkstra's implementation of the minsum(G.Adj) procedure: extract(V 0 ) and outadj(pr, G.Adj).
The input of extract(V 0 ) is the subset V 0 ⊂ V , and the output is The purpose of the outadj(pr, G.Adj) function is to build the adjacency list of the outgoing neighbors of the network G out of the digraph presented by the adjacency list of the ingoing neighbors pr.The function composes each edge e = (i, j) from the digraph given by pr, and takes the corresponding edge length f (i, j) and capacity g(i, j) from the outgoing adjacency list of the original input network G.Adj. {u, V 0 } ← extract(V 0 ) for each q ∈ G.Adj(u) do 6: {d, pr} ← relaxS(u, q(1), d, pr) end for 8: end while G.Adj ← outadj(pr, G.Adj) 10: return {d(n), G.Adj} algorithm (see [18]).We will only note that after the each iteration of the while loop the following properties hold for each vertex v ∈ V 0 : Besides that, when following the proof of the Dijkstra's algorithm, we find out that d(v 1 ) = r(v 1 ).Then in the set V 0 does not exist a vertex that is the one before the last one in a (i 0 , v 1 )-path with length d(v 1 ).
The while loop of the algorithm stops when V 0 = ∅, and then d(j) = r(j), ∀j ∈ V .
We define the network G = (V, E, f, g), where E = {(i, j) ∈ E : i ∈ pr(j)}.In this way, G is the shortest paths subnetwork of the network G, and pr is the adjacency list of the ingoing neighbors of the digraph (V, E), which is verified directly, by using the fact that pr(v) is the set of those vertices u ∈ V for which there exists a path α = (i 0 , i 1 , . . ., u, v) such that x(α) = r(v).
It is clear that using the adjacency list G.Adj, it is easy to compose the list of all (1, n)shortest paths.The following example illustrates this observation.
Example 1.Let G 1 is the network given on the Figure 1.We will find all (1, 5)shortest paths in G 1 .

Solution:
The outgoing adjacency list of G 1 is given by: The function minsum(G 1 .Adj) produces the following result: where The proof follows directly from the computational complexity of the Dijkstra's algorithm, in the case in which the Fibonacci heap data structure is used for the implementation of the priority queue operations (refer to [18]).

B. List of all maximum capacity paths
We denote the capacity of each vertex of the network v ∈ V with c(v).Definition 3. We say that G = (V, E) is a maximal capacity digraph of the network G = (V, E, f, g), if the following two properties hold: In this section we will define the function maxmin(G.Adj) which calculates the adjacency list of the outgoing neighbors G.Adj of the maximal capacity graph G = (V, E).For its implementation, first we will formulate the function capacity(G.Adj) (Alg.3) which calculates the capacity of a vertex n in the network G, again following the Dijkstra's algorithm.
We denote with d a vector of n components that initially has all its elements equal to −∞, except the first element d(1) = ∞.During the calculations of the function capacity(G.Adj) the components of d are growing, and when the function relax(u, q(1), G.Adj, d) end for 8: end while return d(n) completes, the j-th component of d stores the capacity of the vertex j, or in other words d(j) = c(j).Note that by definition c(1) = ∞.We will use the following helper functions in the implementation of capacity(G.Adj).
The function extract(V 0 , d) for an arbitrary subset V 0 = {j 1 , j 2 , . . ., j k } ⊆ V and using the vector d, calculates the pair {j s , V 1 }, where The function relax(i, j, G.Adj, d) for each two vertices i and j, and based on G.Adj and d, performs the following calculations: 1) Defines r = min{d(i), g(i, j)}.
Proposition 3. The function capacity(G.Adj) is correctly defined.
Proof: Let 1 ≤ k < n.We assume that after k iterations of the while loop (see Alg. 3), from the initial set V 0 = V = {1, 2, . . ., n} are excluded k number of vertices u s , and for the resulting set V 0 = V \ {u 1 , . . ., u k } the following properties hold: 1) For each u The above three properties are true for k = 1.After the (k + 1)-st iteration of the while loop, the vertex v k is excluded from v 0 and we get V ′ 0 = V 0 \ {v k }.We will prove that the three properties hold for the set V ′ 0 .From the definition of v k and the fact that the three properties hold for V 0 , it follows directly that properties 1) and 2) are fulfilled for Also, let α be an arbitrary (1, v k+1 )-path and α = (1, u 1 , . . ., u s , v j1 , . . ., v jr , v k+1 ), Algorithm 4 Function maxmin(G.Adj) Input: G.Adj Output: the capacity c of the vertex n and G.Adj c ← capacity(G.Adj) 2: for i ← 1 to n do for each q ∈ Adj(i) do if q(3) ≥ c then pushback( G.Adj, q(1)) 6: end if end for 8: end for return {c, G.Adj} where u i ∈ U ′ , ∀i ∈ {1, . . ., s} and v j1 ∈ V ′ 0 .From the definition of d(v j1 ) it follows that for the capacity y(β) of the path β = (1, u 1 , . . ., u s , v j1 ) it is fulfilled that y(β) ≤ d(j 1 ).Then, for that capacity y(α) of the path α we have that Having the function capacity(G.Adj), we can define the function maxmin(G.Adj), given in Alg. 4.

Proposition 4. The function maxmin(G.Adj) is correctly defined.
Proof: Since the function capcity(G.Adj) that is triggered on the first line of the Alg. 4 is correct, it follows that c is the capacity of the vertex n.The outer for loop on the line 2 defines that adjacency list G.Adj of the graph G = (V, E), where E = {(i, j) ∈ E : g(i, j) ≥ c}.Now we will prove that the graph G = (V, E) is the maximum capacity graph.
If α is a (1, n)-path in the network G with capacity y(α) = c, then for each edge (i, j) of α it holds that g(i, j) ≥ c, and hence, (i, j) ∈ E. This shows that α is a (1, n)-path of G. Now, let β is a (1, n)-path in G.This means that β is a (1, n)-path in the network G, and for each its edge (i, j) it is fulfilled that g(i, j) ≥ c.From here it follows that y(β) ≥ c.However, since c is the capacity of the vertex n, then y(β) = c.

Proposition 5. The computational complexity of the function
The proof follows directly from the complexity of the Dijkstra's algorithm implemented with Fibonacci heap, since the function capacity(G.Adj) repeats exactly it steps.
We will introduce the function list(m, G.Adj) that will help us to clarify the following examples.The argument of the function m is a natural number or ∞.The function maps a list S of (1, n)-paths in G that satisfies the following properties: 2) If m is a natural number, S contains all (1, n)-paths of G if their number is not greater than m.Otherwise, S contains m number of all (1, n)-paths of G.
Example 2. For the network G 1 (see Fig. 1) we will compose the list S y (9) of all (1, 5)-paths with maximal capacity.
On the first step, using the function maxmin(G 1 .Adj), we calculate that:

IV. PARETO OPTIMAL SOLUTIONS SET
We denote with P the set of Pareto optimal paths in the network G = (V, V, f, g).It is clear that where P i are the classes of Pareto equivalent paths.We will find the set P by composing a list Q of all classes of Pareto equivalent paths P i .We will compose the list Q using the following procedure, which we will call Pareto Optimal Paths (POP): Lemma 1.The POP procedure correctly composes the list Q.
Proof: We will prove the correctness of the POP procedure by induction.
Base case.The first iteration of POP defines: The following properties hold: 1) From the definitions of the sets P 1 , Y 1 and Z 1 it follows: If Y 1 = ∅, then W n = P 1 ∪ Z 1 and P 1 is the set of Pareto optimal paths.This means that in (14) the constant k 0 = 1, and the calculations of the POP procedure will stop.
If Y 1 ̸ = ∅, then the constant k 0 > 1.In this case, for each β ∈ Y 1 the inequalities are fulfilled: As a result, also in this case P 1 is a set of equivalent Pareto optimal paths α, and it is correctly included in the list Q.
Here, we set W = Y 1 , and we go back to step 2 of the second iteration of the procedure.
The second iteration of the POP procedure defines: In analogy with the first iteration of POP, the following properties are proved.1) From the definitions of the sets P 2 , Y 2 and Z 2 it follows: 2) If α ∈ P 2 , then by the definition of P 2 the equalities x(α) = d 2 and y(α) = c 2 hold, and hence the paths that belong to P 2 are equivalent.Besides that, from (15) follows that Then, for β one of the following statements hold: Therefore, each path β ∈ Z 2 is dominated by any path from P 2 .Let α ∈ P 2 , and β is such a (1, n)-path, so that β / ∈ Y 1 .Then y(β) ≤ c 1 < y(α).The following two cases are possible: In this case α and β cannot be compared.That proves that if Y 2 = ∅, then the elements of P 2 are Pareto optimal and the equation ( 14) has the form P = P 1 ∪ P 2 , in other words, k 0 = 2. From here it follows that P 2 is correctly included in the list Q, and the termination of the computation of POP procedure is correct.If Y 2 ̸ = ∅, then k 0 > 2. In this case, for each β ∈ Y 2 the following inequalities are fulfilled.
The first inequality follows from the definition of Y 2 , and the second one -from the definition of d 2 and c 2 .Therefore, also in the case in which Y 2 ̸ = ∅, the set P 2 is a set of Pareto optimal paths, and it is correctly included in the list Q.In this case, we set W = Y 2 and we go back to step 2 in the procedure to start its third iteration.
It is clear that in the third iteration of the procedure Inductive step.We assume that after k ≥ 2 iterations of the POP procedure, the following components are defined.
• The sets P i , Z i , Y k , i ∈ {1, 2, . . ., k}, that have no common elements.• The numbers d i and c i i ∈ {1, 2, . . ., k}, for which the following five properties are fulfilled.

It is clear that if
If Y k ̸ = ∅, we implement the (k + 1)-st iteration of the POP procedure.Using steps from 2 to 5, we define: Repeating the proof of the second iteration, we find out that for Since W n has finite number of elements and P i ̸ = ∅, ∀i, then after a finite number of k 0 iterations the procedure POP stops, and we prove that: 2) P 1 , P 2 , . . ., P k0 are sets of equivalent Pareto optimal paths for which a) d i = x(α) and c i = y(α), ∀α ∈ P i and ∀i ∈ {1, 2, . . ., k 0 }; b) d i < d i+1 and c i < c i+1 , ∀i ∈ {1, 2, . . ., (k 0 −1)}.
Corollary 1.For the classes P i of Pareto optimal paths the following holds for each i ∈ {1, 2, . . ., k 0 }.Corollary 2. For P k0 the following equality hods: Remark 1.The list Q can be composed by: first apply the Corollary 2 to separate the set P k0 ; after that, consecutively separate the sets P k0−1 , P k0−2 , and so on, until P 1 is separated.
Every subset P i turns out to be a set of all (1, n)-paths in a special digraph G i .The algorithm that finds a List of Pareto Optimal Paths which we will call LPOP (given in Alg.5), composes a list S of the adjacency lists G i .Adj for each i ∈ {1, 2, . . ., k 0 } by implementing the POP procedure.
We will note that once we have the list S, we can easily obtain a list P ′ of Pareto optimal solutions by taking predefined number of elements for each class P i , as well we can obtain a list P of all optimal solutions.Besides the previously defined functions minsum(G.Adj), maxmin(G.Adj) and capacity(G.Adj), in the formulation of the LPOP algorithm (Alg.5), we will use the function restrict(G.Adj, c), that is defined as follows.
The function restrict(G.Adj, c) takes as an input the outgoing adjacency list G.Adj and the number c.It calculates the adjacency list of those edges (i, j) from G.Adj, for which g(i, j) > c.
Proof: The correctness of the Alg. 5 follows from the correctness of the functions capacity(G.Adj), minsum(G.Adj), maxmin(G.Adj), and from Lemma 1.It is easily verified by proving that the while loop of the lpop(G.Adj) function implements the POP procedure using these functions.
We will examine the first iteration of the while loop  we know that α is a (1, n)-path in G, exactly when α is a (1, n)-path in the network G and x(α) = d 0 .Hence, the set X from the POP procedure is the set of all (1, n)-paths of the subnetwork, defined by the adjacency list G.Adj.By using the definition of the function maxmin(G.Adj), we will prove that its function call implements step 3 of the POP procedure.
Indeed, when applied on the adjacency list of the subnetwork G, the function maxmin(Adj) calculates the maximal capacity c 1 of a (1, n)-path in G, and defines the adjacency list G.Adj of the maximal capacity digraph G of the subnetwork G.This means that α is a (1, n)-path in the digraph G, if and only if it is a (1, n)-path in the network G and y(β) = c 1 .In the notations of the POP procedure, this means that c 1 = max β∈X {y(β)} and α is a (1, n)-path in the digraph G if and only if α ∈ P 0 .For that reason the adjacency list G.Adj is included in the list S (line 10 of Alg. 5).
The step 6 from the POP procedure is implemented in the alternative branch of the if statement which verifies whether the set Y , defined by the step 4 of the procedure, is the empty set.If c 0 = c 1 , then Y = ∅, and the algorithm is terminated.Otherwise, we define the adjacency list RAdj of the subnetwork for which β is a (1, n)-path if and only if it is a (1, n)-path in G, and y(β) > c 1 .
We will illustrate the above proof with the following example.

Fig. 1 .
Fig. 1.Example network G 1 composed by five vertices with length and capacity of each edge given next to it

Algorithm 3
Function capacity(G.Adj) Input: G.Adj Output: the capacity c n of the vertex n
LASKO LASKOV, MARIN MARINOV: LIST OF PARETO OPTIMAL SOLUTIONS OF A BIOBJECTIVE SHORTEST PATH PROBLEM 609 Function lpop(G.Adj) Input: G.Adj Output: the list S with k 0 number of elements RAdj ← G.Adj 2: c 0 ← capacity(G.Adj) more ← true 4: while more = true do {d 0 , G.Adj} ← minsum(G.Adj) Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.Algorithm 5