1. Introduction

Given a bipartite network whose node set IJ is partitioned into m origin nodes (I={i 1,…, i m }) and n destination nodes (J={j 1,…, j n }), the transportation problem (TP) aims at satisfying the demand for a commodity that is supplied from the origin nodes, at minimum cost. Its mathematical model is as follows:

whose dual problem is as follows:

Nothing essential is lost if we assume c ij ⩾0 for all i, j, a i >0 for all i, b j >0 for all j, and

The last condition ensures that supply and demand are equal in order that the TP be feasible.

The TP is well known in the operational research community from theoretical (Intrator and Szwarc, 1988; Arsham, 1992; Dubuc et al, 1999; Sultan and Goyal, 1988) as well as the computational point of view (Glover et al, 1974; Goyal, 1984; Arsham and Kahn, 1989; Bertsekas and Castanon, 1989; Gass, 1990; Wilsdon, 1990; Ji and Chu, 2002). It is well solved (polynomially in the size of the input data) and has many interesting features, the most important of them is that the constraint matrix is totally unimodular.

There are many algorithms to deal with the TP. Most of them are simplex-like such as the Stepping-Stone, the Simplex-Type algorithm in (Arsham and Kahn, 1989), and the Dual-Matrix method in (Ji and Chu, 2002). All of them are pseudo-polynomial methods since their complexity depends on the size of certain data coefficients, often on m, n and C=max i,j c ij . There are also other paradigms, such as the Auction algorithm (Bertsekas and Castanon, 1989) and a well-known and widely used network flow algorithm, the Primal–Dual algorithm (PDA, see eg, Gondran and Minoux, 1979). These two are also pseudo-polynomial in m, n and C.

Let S be a subset of I × J. Given dual variables u i , i∈I, and v j , j∈J, if the following holds (ie, the dual variables are feasible through the set S)

the set S is called a set of admissible cells.

In this paper, we shall present an improvement of the primal–dual method. First, we give a short description of this method. The PDA begins with an initialization phase by setting: u i =0, i∈I, and v j =min i∈I c ij , j∈J (as in Ji and Chu, 2002). This ensures an admissible cell in every row and in every column. After that the algorithm alternates primal and dual phases. The primal phase solves a maximum flow problem through the arcs of S. If the maximum flow obtained is a feasible solution for TP (which is equivalent to say that its maximum value is equal to ∑ i∈I a i ), then it is optimal because of the dual feasibility, and the algorithm stops. If the maximum flow is infeasible for TP (which means that its maximum value is <∑ i∈I a i ), using information from the primal phase, the dual phase comes in play to modify the dual variables, and therefore the set S of admissible cells (see Gondran and Minoux, 1979).

2. A necessary condition for obtaining a feasible solution for TP through a set of admissible cells

The primal phase in the PDA is useless and unproductive, since there is no hope for finding a feasible solution through the admissible cells, unless the following condition holds.

Claim 1

  • If there exists a feasible solution through a set S of admissible cells then

Proof

  • Suppose there exists a feasible solution (x ij ) through S. We have

    since for fixed j, x ij a i for all i∈I and for fixed i, x ij b j for all j∈J. □

Obvious as it may seem, the conditions of the claim are very useful as we shall see.

3. An algorithm for finding a set of admissible cells satisfying the necessary condition

The algorithm we are about to describe does not consider the primal solution. It acts only on the basis of the information obtained from the dual variables and the conditions (1) and (2) when they do not hold. It consists of changing iteratively the dual variables (without regard for the primal), while maintaining their feasibility, until they satisfy the required conditions.

Step 0::

Initialization phase

u i ←min j∈J c ij for all i∈I;

v j ←min i∈I (c ij u i ) for all j∈J;

Repeat Steps 1–3 until conditions (1) and (2) hold

(ie, γ i , δ j ⩾0 for all i,j).

Step 1::

Compute

γ i =∑ j∣(i,j)∈S b j a i for all i∈I and

δ j =∑ i∣(i,j)∈S a i b j for all j∈J.

If there is an i 0 ∈I for which γ i0<0 go to Step 2.

If there is a j 0J for which δ j0<0 go to Step 3.

Step 2::

Compute

Update dual variables

For every j∈J such that set v j v j α;

Go to Step 1.

Step 3::

Compute

Update dual variables

For every iI such that set u i u i α;

Go to Step 1.

This is a real improvement. Every iteration (one execution of the repeat loop) in our algorithm involves only O(mn) operations, whereas the primal phase of the PDA solves a maximum flow in a bipartite network which needs much more computing time (eg, the ‘push-relabel’ algorithm with FIFO vertex selection rule whose complexity is O(m+n)3).

We point out that the condition of the claim is not sufficient. This means that, in general, after our algorithm finishes, the maximum flow need not be feasible. If it is not feasible, we have the recourse to the PDA. On the other hand, as it is known for the PDA, we are never concerned by the degeneracy phenomenon as we are not seeking basic feasible solutions, but admissible cells.

4. A numerical example

The data are given in Table 1. We compute in Step 0 the dual variables u i and v j . The admissible cells are put in rectangular boxes. We compute the values γ i and δ j , the non-negative values of which are marked with the symbol √. Adopting, for example the rule of the smallest negative value γ i or δ j , we choose the second row with γ 2=−25. We compute the value α=6−1−3=2, which corresponds to the cell with (2, 3). We put u 2=1+2=3. We search in this row for admissible cells. Whenever we find one, we change the corresponding v j by the value v j α. The admissible cell in Row 2 is (2, 2). Therefore, we put v 2=0−2=−2 (see Table 2). We begin Iteration 2 by computing the values γ i and δ j (observe that the set of admissible cells has changed). We choose the first row from which we find α=9−5−3=1. So u 1=5+1=6. Since the admissible cell in Row 1 corresponds to j=4 we put v 4=0−1=−1 (see Table 3). In Iteration 3 we compute the values γ i and δ j in connection with the new set of admissible cells. We choose the first column where we find α=8−3−0=5. So v 1=0+5=5. Since the admissible cell in Column 1 corresponds to i=3 we put u 3=1−5=−4 (see Table 4). We compute, in the last Iteration 4, the values γ i and δ j and find all of them non-negative. The algorithm stops.

Table 1 Iteration 1
Table 2 Iteration 2
Table 3 Iteration 3
Table 4 Iteration 4

Now it is time to apply the primal phase of the PDA, that is to search for a maximum flow through the admissible cells (not marked with ◊. See Table 5). The flow is feasible for the TP, and therefore, optimal.

Table 5 The maximum flow through the admissible cells

5. Computational experience

The PDA as well as the primal–dual algorithm with preprocessing (PDAWP) is coded in C (with the GNU gcc compiler) on a Pentium IV 466 MHz. The codes are available from the first author. The two algorithms are run and compared on a set of 30 randomly generated instances. The data are generated by fixing the values of m and n. The costs are uniformly distributed between 1 and 99 and the demands between 50 and 999. The supplies are uniformly distributed between 50 and 999 and then multiplied by the factor n/m. The results are shown in Table 6. The column labeled ‘Dual changes’ gives the number of times the dual solution is changed in the PDA (which is the number of times the maximum flow problem is solved). The label ‘NC iterations’ indicates the number of iterations of the algorithm described in Section 3 for satisfying the condition of the claim. Computing times (including input/output operations) are expressed in seconds.

Table 6 Computational results

The computation and the comparison of the two algorithms show two important facts. The first is predictable. It confirms the advantage of our approach from computing time point of view. The second is rather surprising. Although the size of the instances increases, the algorithm PDA necessitates less iterations and the algorithm PDAWP requires less iterations for satisfying the necessary condition of the claim.

6. Conclusion

Our contribution consists in speeding up the PDA. Indeed, there is no hope to find a feasible solution for the TP through the admissible cells unless the conditions of the claim hold. Therefore, instead of solving in vain, at each iteration, a maximum flow problem, we continue searching for dual variables satisfying our claim. Once such dual variables are found, we search for a maximum flow through the admissible cells they define. If the flow is feasible for the TP, it is optimal and we are done. If it is not, we start the PDA.