A Simulated Annealing Algorithm for Intermodal Transportation on Incomplete Networks

Growing competition in the world enforces the need for an efficient design of transportation networks. Furthermore, a competitive transportation network should also be eco-friendly. As road transportation is responsible for the largest quantities of CO2 emissions, Intermodal Transportation (IT) might be a potential alternative. From this perspective, intermodal terminals location is a cornerstone for building a sustainable transportation network. The purpose of this paper is to study and efficiently solve the Intermodal Terminal Location Problem on incomplete networks. We model this problem as a mixed integer linear program and develop a simulated annealing algorithm to tackle medium and large instances. The computational results show that the obtained solutions using simulated annealing are competitive and close to the exact solutions found by CPLEX solver for small and medium instances. The same developed algorithm outperforms the best found solutions from the literature using heuristics for larger instances.


Introduction
The growth in container-related activities over the last decades is having a remarkable impact on transportation demand. In fact, producers need more and more transportation services to move raw materials and manufactured products to final customers. Furthermore, worldwide economic growth has led to an increasing growth of domestic and international streams of goods and services. Accordingly, handling these flows of goods from production to transportation and storage induce major environmental concerns. Henceforth, in many countries, government, companies and Non-Government Organizations (NGO) are more and more inclined to purchase green products. They are also increasingly urged to incorporate green policy throughout their decision-making process. Therefore, service providers are motivated to provide solutions involving the integration of environmental aspects. It is well known that Intermodal Transportation (IT) is a cornerstone for achieving sustainable development through modal shift from road to eco-friendly modes. However, the IT competitiveness depends on several factors. Among the key elements for building a successful IT system, we cite: (i) optimal use of the most attribute of each transportation mode included in the system; (ii) overall network design; and (iii) optimal location of intermodal terminals. The European Conference of Ministers of Transport (ECMT) defined IT "as the carriage of goods by at least two different modes of transport in the same loading unit (an Intermodal Transport Unit or ITU) without stuffing or stripping operations when changing modes". This type of transport is the door-to-door chain from the initial shipper to the final client. A typical component of IT is the container-based system. In fact, container flows are strongly intermodal. Overall, an intermodal transportation chain is composed by three segments: (i) the first mile segment, which is, in most cases, a short segment (less than 50 km) that uses using road mode (or river); (ii) the main transportation segment, which represents the main segment and uses a consolidation mode such as rail, river, or maritime mode; and (iii) the last mile segment, which is similar to the first mile.
IT is a competitive alternative for door-to-door road transport when the unit cost of the consolidation mode (the second segment) is smaller than in road mode and compensates the additional fees of transshipment operations. IT is very successful and even crucial when we have the two following situations: (i) when the transport is not time critical; and (ii) when there is a large volume of goods. In fact, despite the great speed of the carriage (e.g., rail mode), the overall journey of an intermodal transport chain takes more than door-to-door road. Accordingly, intermodal transportation is not suitable for urgent deliveries. In addition, having a large amount of goods is necessary to benefit from the scale economy generated by the consolidation mode. Intermodal terminals are equipped centers for the transshipment of goods between different modes of transportation. The location of such terminals is a strategic decision that may have a strong and direct influence on the overall operational costs. This problem is addressed in the scientific literature by analogy to well-known Facility Location Problems (FLP) and to Hub Location Problems (HLP). The pioneering paper of [1] is the starting point of Intermodal Terminal Location Problem. The majority of later studies were based on this paper. The authors either proposed a modified version of the models or proposed solving approaches. The proposed works commonly assume that all located terminals are fully connected in a generated solution.
The main contributions of the current paper are the following: • We propose a generalization of a well-known problem in the literature, namely the Intermodal Terminal Location Problem (ITLP) by introducing the incomplete version of the problem relaxing the completeness of the inter-terminal network. • We build a mathematical MIP (Mixed Integer Program) model to formulate the problem and we discuss the meaning of each components. • We solve small and medium instances to optimality using CPLEX solver. • We develop an efficient simulated annealing is adopted to solve real life size instances. • We present extensive numerical experiments and analyze the obtained results.
The remainder of the paper is structured as follows. We review related works in Section 2. The problem description and formulation are detailed in Section 3. We develop the simulated annealing method in Section 4. Numerical results are reported in Section 5. We provide results analysis and comments in Section 6. We give concluding remarks and perspectives in the last section.

Related Works
Although the research in the ITLP is a very young field, the number of papers related to this subject is steadily increasing. This reveals the importance of the Intermodal Terminal Location for both academics and industry. A review synthesis of works done in the field of IT can be found in the paper of [2]. Specifically, the Network Design Problem (NDP) is a well-studied problem in the framework of intermodal transportation. A stochastic approach for robust intermodal transportation plans in a real-world network was proposed by [3]. The relationship between IT and supply chain was discussed by [4]. They developed a scenario-based estimation to quantify the demand for the Italian IT Network. Owing to its safety, IT is increasingly considered for hazardous materials (Hazmat) distribution. For instance, Ref. [5] suggested using scheduled direct and faster trains to reduce the transport risk of hazmats. The proposed solution is as a bi-objective model solved by a tabu-search heuristic. Ref. [6] studied optimization challenges regarding rail processes at maritime terminals. They provided a survey of the related works in this important node of IT. Ref. [7] proposed a fuzzy version of the Intermodal terminal Location Problem. He studied the case in which the decision-maker may suggest locating a number of rail terminals around a fixed number "p". This scenario is modeled by adding a fuzzy constraint to the developed model. In [8], a compromise programming technique combined with fuzzy modeling is proposed to tackle the location problem of real-road terminals. The developed Multi-Criteria Decision Making Method (MCDM) is considered to model the conflicting objectives of different stakeholders involved in the intermodal transportation chain. Ref. [9] modeled the intermodal container terminal location problem in urban areas as a non-linear Mixed Integer Program (NLMIP). They considered an entropy maximization approach and a Lagrangian relaxation technique to solve the problem in real-life case. A combined Delphi and Analytical Network Process (ANP) is proposed in [10] to solve the Intermodal Terminal location problem in accordance with sustainable requirements. Ref. [11] proposed a modeling approach for Intermodal Terminal Hub Location Problem with transfer modes and considering two-types of hubs. They proved that the problem is NP-hard and reported numerical results on real-world container distribution networks. In addition, multi-agent simulation tools are used in [12] to simulate the seaport-dry port system as intermodal terminal. The developed model is validated on the Ningbo-Zhouhan port. A combined SWOT analysis and Analytic hierarchy process (AHP) and Preference Ranking Organization Method for Enrichment of Evaluations (PROMETHEE) is proposed in [13] to study the location of intermodal terminals. Ref. [14] studied user-specified reliability levels using a model that provide the optimal route through road links, rail links and via intermodal terminals. A Markov Chain Monte Carlo method is proposed in [15] to solve the intermodal terminal location problem. The proposed approach is a two-layer approach that efficiently solve the problem with routing decisions. Computational results are reported. Ref. [16] proposed a Best-Worst Method (BWM) with a consensus model from the perspective of shipping line companies. The objective is to select the inland terminal locations from potential sites. A sensitivity analysis was performed by authors to assess the impact of the increasing volume of containers on the desirable locations. In [17], a quantitative method is proposed to locate inland rail terminals based on cost savings. The model was validated in real-world railway network from Croatia. Table 1 summarizes the main research works related to the intermodal transportation. In addition to review papers that survey research works in this field, terminal location problems, route selection, network design, transport policy, and integrated location and routing problems are the main topics of related works. Terminal location problems are the strategic planning problem aiming to open a fixed number of rail terminals among a set of potential sites given by a decision-maker. The route selection problem is a tactical decision problem where the objective is to assign the best combination of transportation modes for each commodity. The network design is another example of strategic planning that aims to determine the initial network to use for the transshipment of demands. Some papers deal with transport policy issues where the general management rules of transport are addressed. Moreover, few research papers considered the integrated location and routing problem in the frame of the intermodal transportation. The particularity of the last problem lies in the fact that the two decisions of terminal location and the routing of demands are tackled simultaneously.  [46,47] We can clearly remark that terminal location problems are steadily attracting researchers attention. This proves the importance of locating rail terminal in the design of the intermodal transportation networks. In fact, while we can find different types of research papers in the field of intermodal transportation, various aspects of terminal location problem are studied. We can cite the strategic choice of potential sites, the study of uncertainty, the use of multi criteria decision making, etc. Furthermore, we provide in Table 2 an overall classification of the models used for ITLP and a summary of works using either exact or heuristics methods. The first papers dealing with the studied problem use deterministic models to formulate the problem and exact methods to solve it. The use of such formulations enables applying the rich analytical deterministic tools to solve the problem. More recent papers use either stochastic or fuzzy modeling to capture the uncertainty, vagueness, and complexity of the problem. In such case, heuristics are often used as no exact method is able to solve non-linear or large size problems.

The Complete and Incomplete Version of the ITLP
The ITLP aims to determine which of a set of potential sites locations to open as terminals and how to route the supply and demand of a set of customers (representing zones of supply and demand) through the network. The transportation of these demands may be done using a unimodal link (road mode) or a combination of two modes (rail-road). The objective is to minimize the total cost of the transportation as well as the installation cost of terminals. Up to now, the studied version of the ITLP introduced by [1,25], studied by [27,29], etc. assumes that the inter-terminals network is a complete graph. Figure 1 shows a small complete intermodal terminal network with three terminals T 1 -T 3 and 5 customers (1-5). As in hub problems, the assumption that the induced graph by the located terminals must be a fully-connected graph (complete graph) is a restrictive constraint that may increase the total investment cost. Moreover, this restriction is non-realistic since, in real life cases, a rail-road network is usually incomplete due to the heavy cost incurred for the construction of its infrastructure. The current work relax this assumption and proposes a general model in which the located rail-road terminals are partially connected. Figure 2 shows a small incomplete intermodal terminal network. The induced graph by terminals is not incomplete. For instance, there is no rail link between rail terminals T 1 and T 2 .

Mathematical Model
Let us introduce the following parameters: I : Set of senders/receivers. K : Set of candidates sites for intermodal terminals.
l : Fixed number of links to be established. q ij : The amount of goods to be shipped from the customer i to customer j. c km ij : The intermodal unit cost for transportation from sender i to receiver j using the two rail terminals k and m.
c ij : The unimodal unit cost for transporting goods by truck from sender i to receiver j.
f k : The investment cost for location of a terminal on site k. C k : Fixed capacity of a located rail terminal k.
Let us consider the following decision variables: l km : A binary variable equal 1 if the rail link between k and m is established, 0 otherwise. Accordingly, l kk equal 1 if the node k is a terminal, 0 otherwise. r ij : The amount of goods routed from the sender i to receiver j. s km ij : The amount of goods routed from the sender i to receiver j through the two rail terminals k and m.
The incomplete location problem of intermodal terminals may be formulated as follows: The linear objective function (1) minimizes the total cost for routing goods and the total cost for locating rail terminals. The first term expresses the cost for routing goods through the intermodal chain: from the initial sender i to the first rail terminal k using trucks, then from the terminal k to the second terminal m by rail mode, and later from the terminal m to the final customer j using road mode. The second term is the cost for routing goods directly using trucks from i to j while the third term is the overall investment cost for locating rail terminals. ∑ k,m∈K s km ij + r ij = q ij , ∀i, j ∈ I The constraint (2) states that the sum of goods transported via the intermodal chain through the terminals k and m and that routed directly from i to j equal exactly the total goods to be transported from i to j.
The constraint (3) imposes to respect the terminals capacities.
The two inequalities (4) and (5) guarantee that a rail link is established if the two adjacent terminals are located. In fact, we have the following implications: if l kk = 0 or l mm = 0 =⇒ l km = 0 Equation (6) guarantees that a rail link is opened in the two directions.
∑ k,m∈K l km = l Equation (7) guarantees the establishment of exactly the given number of rail links. The decision variable l km controls the rail links to be opened. Thus, if due to economic, geographic or environmental constraints a railway link is impossible between two terminals, such scenario is forbidden by l kl = 0. s km ij ≤ l km , ∀k, m ∈ K, ∀i, j ∈ I Equation (8) forbids the rail flows between terminals if they are not located. The developed model is a Mixed binary Integer Program (MIP). If we denote by p the number of potential sites (the cardinal of K) and by n the numbers of senders/receivers (the cardinal of I), then the program has n 2 p 2 + 3p 2 + n 2 + p + 1 constraints and n 2 p 2 + n 2 + p 2 variables. Proposition 1. If the triangular inequality holds for the transportation unit cost (assumed in the literature and holds for dataset used in the numerical experiments), then s kk ij = 0, ∀i, j ∈ I, ∀k ∈ K.
This proves that the constraint s kk ij = 0, ∀i, j ∈ I, ∀k ∈ K, used in several mathematical models in recent papers in the literature, is an unnecessary constraint [27].

Proof.
To express the fact that rail mode generate a scale economy (consolidation mode), the intermodal unit cost may be written as follows: c km ij = c ik + αc km + c mj where c ik is the cost between i and k, c km the cost between the two rail terminals, and c mj is the cost between the rail terminal m and the final receiver j with α < 1. We denote by α the discount coefficient expressing the scale economy generated by the rail mode. This parameter is supposed to be equal to 0.5 in [27]. Since c kk ij = c ik + αc kk + c kj = c ik + 0 + c kj = c ik + c kj ≥ c ij (the last inequality is the triangular inequality), and as we aim to minimize the total cost expressed by Equation (1), given the conservation flow constraint of Equation (2), it follows that all goods from i and destined to j will be routed directly using trucks, and then r ij = q ij , which implies that ∀i, j ∈ I, s kk ij = 0 Proposition 2. Let l be the number of rail links to be located and |K| = p the cardinal of the candidates sites set. The Incomplete ITLP modeled by the aforementioned MIP is unfeasible when l > p(p−1)

.
Proof. It is an immediate result of the handshaking lemma (a classical lemma in graph theory: an incomplete graph with p vertices has at most p(p−1) 2 edges).

Model Variants
The basic model introduced in the previous subsection may be enriched by adding other constraints to match other realistic situations.
For instance, one may model the situation where both the number of rail links l and the number of rail terminals t are given by the decision-maker by adding the constraint: In such case, as the number of rail terminals is known, the objective function becomes: Furthermore, the basic formulation may be extended to consider potential handling operations in the located rail terminals. In fact, let us denote by h km the handling cost of a unit of goods transhipped from the rail terminal m to the rail terminal k (symmetry may be considered by assuming that h km = h mk ). Therefore, the objective function considering asymmetric handling costs becomes:

General Scheme
The Simulated Annealing (SA) algorithm is a local search metaheuristic used for finding an approximate global optimum. The SA is based on the Metropolis-Hastings criterion enabling to escape from being trapped in local minima. We used this single solution metaheuristic for the two reasons: (1) to solve the intractable instances for the CPLEX solver; and (2) because this method was successfully applied to similar location problems as in [51,52]. We adopted the classical general scheme of the SA, as shown in Algorithm 1.

Initial Solution and Neighbors Generation
To generate the initial solutions, we use Algorithm 2.

Algorithm 2: Generation of initial solutions
Randomly generate rail links ; if l kk = 1 and l mm = 1 then To generate neighbors of a solution, we define three types of moves: • Permutation of rail links: We modify a given initial solution by interchanging two non-connected rail terminals with two connected ones. This move has the particularity to generate new feasible solution since the number of rail links remains the same. • Permutation of terminal locations: This move interchanges a terminal node with a non-terminal one. As the generated solution may be infeasible, a correction phase is conducted to guarantee the respect of the model constraints. In fact, if the capacity of the new terminal is less than the old one, a new rail link is added. Otherwise, the flows are routed trough the new terminal. • Changing the number of terminals: From a given solution with t number of terminals, we generate neighbor solutions by adding or deleting one or several rail terminals. The correction phase is made to respect links and capacities of terminals.
All this moves are used in a sequential way in the implemented algorithm. The best solution from the initial iteration to the last one is saved in a specific variable called BestSol and reported in the end of the SA execution.

Parameters Tuning of the SA Method
Since the simulated annealing algorithm is a parametric approach, the setting of values of these parameters has a significant effect on the objective function. For this, we ran a set of experiments to determine the most suitable values to reach good solutions. The controller parameters of the developed simulated annealing are the initial temperature T 0 , the cooling rate δ, the neighbor moves rate γ and the minimal temperature T min . The correspondent intervals of these parameters are presented in Table 3. The Tagushi design of MINITAB software [53] was used for parameters tuning of the simulated annealing. Equation (12) defines the signal-to-noise ratio (S/N) of Tagushi method. This ratio expresses the effect of parameters on the objective function. A maximum value of this ratio minimizes the effect of uncontrollable factors. This ratio is given by: The constant z represents the total number of experiments and y i is the result of the experiment i. The results of parameters tuning using Tagushi design for a typical set of training instances are illustrated by charts of Figures 3 and 4.

Numerical Results
All experiments were done on a notebook Intel core i5, 2.53 GHz with 8 GB RAM memory under Windows 10. The SA algorithm is coded using Python 3.

SA Solutions vs. Optimal and Best Found Solutions
To assess the efficiency of the developed Simulated Annealing (SA) algorithm, we used the dataset introduced in [27]. In this dataset, customers and potentials sites coordinates are randomly generated between (0, 0) and (10 4 , 10 4 ). The goods demands are generated from the interval [0, 500]. The investment cost f k in the interval [0, 5.10 5 ] and potentials sites capacities are drawn from [0, 10 4 ]. After that, the direct cost c ij is equal the Euclidean distance between the customers i and j, while the intermodal cost is computed using the following formula c km ij = c ik + 1 2 c km + c mj . The objective function values of the studied model are detailed in [31]. For this, we lunched the SA method for 10 runs and we choose the best solutions. As shown in Figure 5, saving best solutions in an archive set enables to get a decreasing values of the objective function. Accordingly, the developed SA method is a descent strategy. Table 4 reports the gaps of the SA solutions to CPLEX solutions. The second column gives the gap of the best obtained solution from the CPLEX solution, while the third column gives the average gap of the conducted 10 runs.
Note that we used the following formula to compute the gap: The obtained results show that our method is efficient and may even reach the optimal solutions for the small size instances and reaches best solutions for medium instances in short times.
To further evaluate the performance of the developed SA, we compare our results with those of [27]. Note that, for a meaningful comparison, since our problem and that studied in [27] are not the same, we made a transformation of our problem to get the same structure of [27]. In fact, since our problem is a generalization of [27], we added the two following constraints to our formulation: These two equations transform our incomplete (general case) to complete problem. The gaps of the last column in Table 5 show that our solutions by SA approach outperform best found solutions by [27] in 11 out of 15 instances.

Completeness Effect
To study the completeness effect on the solution cost, let t = ∑ k∈K l kk which denotes the number of terminals in a solution. Then, we define the "Completeness rate" of inter-terminals induced graph as follows: We modify the linear program introduced in the Section 3 by changing Equation (7) to become: ∑ k,m∈K l km ≤ l (14) As shown in Figure 6, the completeness rate of the SA solutions is overall between 50% and 100%. This means two important things: (1) a sufficient number of rail links is necessary to reduce the total transportation cost as the rail mode enables to achieve scale economy by transporting huge volumes of flows; (2) a fully connected terminal network is unnecessary and increases the total transportation cost.
To further study the completeness rate effect on the generated solutions, we plot in Figure 7 the completeness rate over iterations of the archived best solutions for an instance with 40 nodes. This figure highlights that, while a sufficient number of rail links is necessary to reduce the total transportation cost (more than 50%), a fully connected rail network may not achieve an efficient transportation plan.

Concluding Remarks and Perspectives
Research in the freight terminal location problems is still in its beginning. However, this domain is attracting the attention of scholars and policy makers. The current paper contributes to this knowledge area by optimizing the location of intermodal terminals. Such terminals are the cornerstone in each intermodal freight system. Considering the incomplete inter-terminal network is relatively ignored in the literature, a mixed integer model for the intermodal terminal location problem is proposed in the current work. Due to its NP-hardiness, we develop an efficient Simulated Annealing Algorithm to tackle medium and large instances. Our solving approach allows us to retrieve the optimal solutions for small instances considering CPLEX 12.10. solver. In addition, our best found solutions are competitive in the medium and larger sized instances as the worst CPU time required appears to be less than 10 s. Specifically, on larger instances, our algorithm outperforms the best integer solutions returned by CPLEX, in both the objective cost and CPU time.
The numerical results show that considering the general case where all located terminals are not necessarily fully connected may substantially reduce the overall costs of location and transportation. Accordingly, this leads to the reduction of committed investments costs.
While we provide an efficient modeling and solving of the incomplete intermodal terminal location problem, an important extension of the current work is to integrate the routing decision in the mathematical formulation. In fact, the proposed formulation assumes that routing is efficiently managed once locations of terminals are fixed. However, routing and location are interdependent. Therefore, the current paper may be a first step to integrate routing decisions for the incomplete version of the ITL.
From methodological perspective, using other solving methods may be explored in a future study. Especially, exact methods such branch and bound or cutting planes methods may be adapted to the studied problem.
Another potential opportunity for future research may be the integration of the hub location formulation in the framework of the intermodal location considering the hubs as intermediate warehouse for the shipped commodities. Furthermore, it may be interesting to deal with a case study in the framework of intermodal transportation to consider more realistic situations.