Abstract

This paper proposes to address the problem of the simultaneous optimization of the liner shipping route and ship schedule designs by incorporating port time windows. A mathematical programming model was developed to minimize the carrier’s total operating cost by simultaneously optimizing the port call sequence, ship arrival time per port of call, and sailing speed per shipping leg under port time window constraints. In view of its structure, the nonlinear nonconvex optimization model is further transformed into a mixed-integer linear programming model that can be efficiently solved by extant solvers to provide a global optimal solution. The results of the numerical experiments performed using a real-world case study indicated that the proposed model performs significantly better than the models that handle the design problems separately. The results also showed that different time windows will affect the optimal port call sequence. Moreover, port time windows, bunker price, and port efficiency all affect the total operating cost of the designed shipping route.

1. Introduction

Container lines transport containers following fixed sequences of port calls with a regular service frequency [1, 2]. A port, the core of the shipping logistics chain, needs to service a large number of ships on different routes according to predetermined schedules. In recent years, with the fast growth of container trade, ports are becoming more and more congested. In April 2017, Shanghai port, the world’s largest container port, was not immune from port congestion, which surprised the industry. Meanwhile, slow port expansion projects in some nondeveloped countries and regions have exacerbated this phenomenon. This shows that the daily service capacity of a port is limited, and, therefore, it cannot guarantee that a ship can be serviced as soon as it arrives at the port. As a result, the availability of ports, i.e., port time windows, is the first factor to be considered when designing a ship schedule [3].

However, in practice, it is not feasible to design a ship schedule under port time windows constraints for a given route without considering the optimization of the port rotation. We take the S2 service route operated by SITC Container Lines as an example to illustrate. The route, consisting of four ports, provides services following the port call sequence of Taicang-Shanghai-Osaka-Kobe-Taicang. Without the constraints of time windows, only one container ship is needed to maintain weekly service frequency on the route. Now, we assume that the time windows of each port call are as follows: Taicang (Tuesday), Shanghai (Monday), Osaka (Friday), and Kobe (Thursday). As shown in Figure 1(a), owing to the constraints of time windows, a ship needs to wait in the anchorage for nearly 6 days before arriving at Shanghai as well as at Kobe following the original port call sequence on the route. (It takes only a few hours for a ship from Taicang to Shanghai as well as from Osaka to Kobe.) Therefore, it would take 3 weeks for a ship to perform a round-trip journey, which means that three ships need to be deployed on the route. As a result, the weekly operating cost of the route will increase significantly. In addition, the longer transit time leads to delay in delivery, which also results in loss of sales, goodwill, and service satisfaction to the carrier. However, if we change the port call sequence to Shanghai-Taicang-Kobe-Osaka-Shanghai, as shown in Figure 1(b), it only takes 1 week for a ship to perform a round-trip journey. From the above comparison, it is easy to see that the adjusted port call sequence can greatly reduce the operating cost compared with the original one, subject to the constraints of time windows. The reason why we can make such an adjustment is the operational flexibility of the route. Taicang and Shanghai are both ports in China. For the carriers operating foreign trade, there is no cargo flow carried by the ships between the two adjacent O-D port pairs in the same country except for special reasons (empty container repositioning, transshipment operations due to skipping port calls, etc.). Moreover, the distance between the two ports is very short. In fact, such similar port pairs or groups are fairly common, such as the Kansai port group, which includes the Osaka and Kobe ports; the Kanto port group, which includes the Tokyo, Yokohama, and Kawasaki ports; the Manila north and south ports; and the Laem Chabang and Bangkok ports.

Hence, liner shipping route design and ship schedule design with time windows should be carried out at the same time rather than one after another. Otherwise, the designed shipping route and schedule will not be feasible in reality. This study provides a decision support system for container lines to improve their services and reduce their operating costs.

1.1. Literature Review

This study is related to two research topics: liner shipping network design and ship scheduling. The former has been a hot topic in the field of liner shipping research over recent decades. The work of Rana and Vickson [4] was the early contribution in this field. They formulated a mathematical model to determine the port rotation, the number of trips each ship makes in a planning horizon, and the amount of transported cargo for multiple ships by employing a Lagrangian relaxation method with the objective of maximizing the total profit. Shintani et al. [5] proposed a two-stage programming model to design a single shipping route by explicitly considering empty container repositioning. A genetic algorithm-based heuristic was employed to solve the problem. The work by Agarwal and Ergun [6] was the first academic study to consider the weekly frequency constraint in the design of liner transport service networks composed of two or more routes. However, the sailing speed in their model was predetermined, not variable. Álvarez [7] jointly optimized the optimal routing and fleet deployment by using a combined tabu search and a column generation-based heuristic method. Wang and Meng [8] proposed a novel aspect in the liner shipping network design problem by reversing the port rotation direction. They demonstrated that the operating cost of a container liner shipping network could be significantly reduced by reversing the port rotation direction of liner service routes. Song and Dong [9] proposed a shipping route design problem that considers ship deployment and empty container repositioning with the objective of minimizing the total cost. A three-stage solution method was proposed to deal with the problem. However, they failed to consider the speed variation in each leg. Wang and Meng [10] considered the delivery deadlines in the liner shipping network design problem. The problem was proved to be NP-hard and solved by a column generation-based heuristic method. Karsten et al. [11] examined the liner shipping network design problem that considers the coordination between vessels and cargo transit time restrictions. Wang et al. [12] designed a single intercontinental service by simultaneously optimizing the sailing speed, the port rotation of selected ports in a given set of ports, and the shipping demand for each type. A mixed-integer linear programming model was developed and solved by using exact algorithms. Wang et al. [13] proposed a port call adjustment problem to determine the new additional port calls to be inserted and the existing port calls to be removed for the routes in a shipping network. A two-phase solution method was proposed using two heuristic methods.

A considerable number of papers on vessel scheduling have been published in the past. Meng and Wang [14] jointly determined the sailing speed, service frequency, and fleet deployment for an operating strategy problem. An exact ε-optimal algorithm was proposed to solve the problem. Wang and Meng [15] investigated the optimal sailing speed on each shipping leg and the optimal fleet size in a liner shipping network. They reformulated the bunker consumption function from historical operating data, which was applied to a real-world case study. Wang et al. [16] simultaneously optimized the schedules of shipping routes and cargo allocation with the objective of reducing the demurrage cost. Cheaitou and Cariou [17] tried to optimize the sailing speed under a semielastic demand. Their results showed that the slow steaming strategy may not be optimal. Hvattum et al. [18] developed an exact algorithm to optimize the sailing speed for a fixed sequence of port calls with hard time windows. They proved that the optimal speeds can be obtained in quadratic time. Wang et al. [19] studied the liner ship schedule design problem with hard time windows. Wang et al. [3] studied the same problem with consideration of the availability of each berth with the objective of minimizing the bunker cost, operating cost, and inventory cost for a single shipping route. Alharbi et al. [20] examined the ship schedule design problem with time windows for a liner shipping network. In their model, the vessel could utilize a premium berth with high penalty when violating the berth time window. Aydin et al. [21] focused on the sailing speed optimization problem with time windows on a single voyage that considers the waiting cost of early arrivals. Dulebenets [22] investigated the vessel scheduling problem with a heterogeneous fleet for a single shipping route. Zhen et al. [23] proposed an integrated planning model to simultaneously determine the fleet size, ship schedule, sailing speed, and cargo allocation for a given shipping network. Kim et al. [24] designed a simple single route by considering the variable sailing speed on each leg and the fleet size to maximize a carrier’s profit. Koza et al. [25] proposed a mixed-integer programming model for the integrated liner shipping network design and scheduling problem that incorporates the transshipment times, together with a solution method based on column generation. Jiang et al. [26] developed a mixed-integer linear programming model for the near-sea liner shipping schedule design problem considering big customers’ time preferences.

From the abovementioned literature, we can see that the two key and interrelated problems, namely, liner shipping route and ship schedule, are usually treated separately by most existing studies. Although a few articles combined the two, they did not consider the port time window constraints and the variable sailing speed per shipping leg. This study tries to fill the gap. As mentioned in Section 1 of this study, it is highly desirable to simultaneously optimize port call sequence and ship schedule subject to the constraints of time windows.

1.2. Objectives and Contributions

This study focuses on jointly optimizing the port rotation, sailing speed on each leg, and ship arrival time at each port of call for a given set of ports with time windows while minimizing the sum of the vessel operating cost, fuel cost, and waiting cost. This problem has practical significance because it can help container lines improve service levels and reduce total operating costs.

The contribution of this study is threefold. First, it simultaneously determines the optimal port rotation, the optimal sailing speed per shipping leg, and the optimal arrival time per port of call under port time window constraints. Second, a tailored mixed-integer programming model was developed that considers the special structure of the problem. Finally, some useful management insights obtained from case studies can provide support for container lines.

2. Problem Description and Mathematical Formulation

We consider a group of ports, denoted by a set . The set has number of physical ports indexed by and . The voyage from port to port is called the shipping leg . We define the set of all shipping legs as . Since the shipping route is a closed circle, we can arbitrarily choose one port as the first port of call. In this study, we choose physical port 1 as the first port of call. The details on problem description and key constraints are provided in the following sections.

2.1. Port Rotation

The port rotation of a shipping route forms a loop in practice. We use the binary variables to indicate whether the shipping leg is used. In this study, we assume that each physical port can be visited only once during a round-trip journey; then, the constraint should hold for each port .

2.2. Bunker Consumption

A ship’s unit bunker consumption largely depends on its sailing speed. We denote by the variable the sailing speed of a ship on leg , and we define as the bunker consumption in tons per nautical mile at the sailing speed on leg . Based on the results of existing studies [15, 27], the fuel consumption is a power function of the sailing speed: , where and are two coefficients estimated from the practical data. Let be the unit bunker price (in US dollars per ton), and let be the length of shipping leg ; then, the weekly fuel cost of the designed shipping route can be calculated by .

2.3. Port Time Windows

Since the service provided by the ports has time windows, we set to represent the time windows at port . We define the time 00 : 00 of a certain Sunday as time 0. Let be the ship arrival time at port , and let be the time windows at port in the first cycle. For example, means that port is available on Monday, Thursday, Friday, and Saturday, and a ship can call the port for service in these periods. Otherwise, the ship must wait in the anchorage until the time windows open again. Owing to the characteristics of weekly service frequency of liner shipping, ships call port at the same time in different weeks. Hence, an arrival time is feasible if and only if . For example, a ship can call port for service at because .

It should be noted that, as a ship has a minimum sailing speed, it needs to wait in the anchorage when arriving at the port outside the time windows. On the one hand, it is well known that, when a ship is waiting in the anchorage, although its main engine is not running, its auxiliary engine still does. Thus, there is still fuel consumption on the ship. On the other hand, the inventory cost of containers also rises with the increase in waiting time. Hence, the weekly waiting cost of the designed route can be calculated by , where the variable is the waiting time for the ship in the anchorage at port after using shipping leg , and is the fixed waiting cost per hour in the anchorage at port .

2.4. Fixed Operating Cost

A fleet of homogeneous ships is deployed on the designed shipping route to achieve a weekly service frequency for each port of call. Let the integer variable represent the number of ships deployed on the route. Thus, the weekly operating cost of the designed route can be calculated by , where is the fixed weekly operating cost in US dollars per week per ship.

2.5. Mathematical Model

The problem of the simultaneous optimization of the liner shipping route and ship schedule designs with time windows aims to determine the optimal port call sequence, the optimal sailing of each shipping leg, and the optimal ship arrival time that satisfies the time window constraints at each port to minimize the total cost including the vessel operating cost, fuel cost, and waiting cost. Before presenting the model, we list the variables and parameters below.Decision variables:: binary variables indicating whether shipping leg is used on the designed shipping route: auxiliary variables for subtour elimination constraints: sailing speed of a container ship on leg : sailing time of a container ship on leg : waiting time for the ship to be in the anchorage at port after using leg : integer variable representing the number of ships deployed on the designed shipping route: ship arrival time at port : the time when the ship returns to the first port of call after visiting all ports of call of a round-trip journey: ship arrival time at port in the first cycle: numbers of containers of container shipping demand on leg Sets:: set of all ports: set of all shipping legs: set of all container shipping demands: set of nonnegative integersParameters:: number of ports on the designed shipping route: length of shipping leg : minimum sailing speed of a container ship on leg : maximum sailing speed of a container ship on leg : port time windows at port : fixed waiting cost per hour in the anchorage at port : fixed dwell time at port : fixed weekly operating cost of one container ship: maximum number of ships deployed on the designed shipping route: volume of container shipping demand : port of loading and port of discharging for container shipping demand : very large numbers: unit bunker price: a coefficient estimated from the practical data: a coefficient estimated from the practical data: capacity of one container ship

Thus, the problem of the simultaneous optimization of the liner shipping route and ship schedule designs with time windows can be formulated as follows:subject to

The objective function (1) minimizes the total weekly operating cost, which amounts to the sum of the vessel operating cost, fuel cost, and waiting cost. Constraint (2) ensures that each port can be visited only once. Constraint (3) is the Miller–Tucker–Zemlin (MTZ) subtour elimination constraints. Constraints (4) and (5) jointly guarantee that if leg is used, or otherwise. Constraint (6) defines the lower and upper bounds of the sailing speed on each shipping leg. Constraints (7) and (8) jointly impose the port time window restrictions at each port. Constraints (9) and (10) establish the relationship between the port call arrival times: if leg is used, the ship arrival time at port is equal to the sum of the ship arrival time at port , the fixed dwell time at port , the sailing time on shipping leg , and the waiting time in the anchorage at port . Constraints (11) and (12) define the time when the ship returns to the first port of call after a round-trip journey. Constraints (13) and (14) indicate that if leg is not used, the waiting time for the ship in the anchorage at port should be 0. Constraint (15) enforces a weekly service frequency on the designed shipping route. Constraint (16) indicates that the number of ships deployed on the shipping route is a positive integer. Constraints (17) and (18) indicate that all container shipping demands should be satisfied while ensuring the container flow balance. Constraint (19) indicates that the total number of containers loaded on each shipping leg should not exceed the capacity of one container ship.

3. Approximate MILP Model

To take advantage of state-of-the-art MILP solvers such as CPLEX to solve the developed model [M1], we need to deal with three difficulties: (i) the bunker consumption, which contains two variables and power functions in the objective function, is nonlinear; (ii) constraint (5) contains the reciprocal of the sailing speed; and (iii) both the “mod” operation and the set lead to a nonconvex domain. These three difficulties make the proposed model challenging to solve using the MILP solvers. Next, we deal with these difficulties and transform the model [M1] into an approximate mixed-integer linear programming model.

3.1. Linearization of the Bunker Consumptions

First, we define a new variable to be the reciprocal of the sailing speed

Then, constraints (5) and (6) can be transformed into linear constraints associated with the variable :

We introduce an auxiliary variable to replace the nonlinear combinational items and in the objective function. We define . Hence, we have if , i.e., leg is used; otherwise, . By introducing a very large number , we can replace the above if-then complementary conditions equivalently with the big-M reformulations of the following constraints:

By introducing the auxiliary variables , we define . Objective function (1) can be transformed by constructing its epigraph form as follows:subject toand to constraints (23) and (24).

We define . Because and , is convex in . We divide the interval into equal segments:where is a given integer number. As long as is large enough, the error of the approximate model can be controlled within an acceptable range. We use to represent any value of the equal segments. From the property of the first-order condition of the convex function, we can obtainwhere represents the tangent line of at the point . Thus, the function can be converted into the following piecewise function:

By combining constraints (26) and (29), we can obtain separate inequalities

3.2. Transformation of the Time Windows

Since the set may be a discrete domain, by transforming the hard time windows into soft time windows, we can obtain piecewise functions that can eliminate the discrete domain and can be linearized. For example, we assume that the time windows at port are Tuesday and Thursday. By introducing a very large number , we can express the penalty function at port as follows:

The penalty cost is high enough to limit a ship to call the port only within the time windows. For example, as shown in Figure 2, ships can only call port on Tuesdays and Thursdays.

Next, we linearize the penalty function . Since has turning points, let be the turning points and let be the set of all turning points at port . For example, as shown in Figure 2, there are 10 turning points: . By introducing the continuous auxiliary variables and the binary auxiliary variables , we can transform the piecewise function into the following linear constraints:where can realize the function linearization and can limit the value range of to limit the interval of .

Thus, the total penalty costs of time windows at all ports are as follows:

By introducing the auxiliary variables , we have

It is noticeable that constraints (42)–(44) are equivalent to constraint (7).

Hence, an approximated mixed-integer linear optimization model [M2] can be obtained, and it can be solved efficiently by using state-of-the-art solvers, subject to constraints (2)–(4), (9)–(19), (21)–(24), (30), and (32)–(44).

4. Case Study

This section presents the results of the numerical experiments we conducted on the PBT1 service operated by SITC Container Lines to assess the applicability of the proposed models and algorithms. We first tested the validity of the solution approach to determine how it can significantly reduce the total cost. Then, we conducted a sensitivity analysis to test some key input parameters in the solution approach.

4.1. Parameter Settings

The PBT1 service includes nine ports of call, and its port rotation is Qinhuangdao (QHD) ⟶ Xingang (TXG) ⟶ Dalian (DAL) ⟶ Nagoya (NAG) ⟶ Yokkaichi (YCI) ⟶ Shimizu (SZU) ⟶ Tokyo (TOK) ⟶ Yokohama (YOK) ⟶ Kawasaki (KAW) ⟶ Qinhuangdao (QHD). We assumed that 2000-TEU ships were deployed on the designed route, i.e., . The operating cost , the maximum speed , the minimum speed , the bunker price , the maximum number of ships deployed , and the fixed waiting cost in the anchorage at each port US$100/h. For ease of calculation, we take the tightest upper bounds; then, the very large numbers (, and ) can be set as , , , and , respectively. To avoid violating the time window constraints, the very large number was set to 300,000. The number of segments for the interval was . The bunker consumption function of each leg was set to . The port-to-port distance, the dwell time at each port, and the time windows are shown in Tables 13, respectively. CPLEX 12.8.0, programmed by the MATLAB toolbox YALMIP running on a 2.5 GHz Intel Core Duo laptop with 16 GB of RAM, is called to solve the MILP model. The MILP model inputting the above parameters has 3,473 continuous variables, 207 integers (198 binary) variables, and 21,466 constraints, and the optimal solution with 0 gap among objective bounds can be found in about 1 hour.

4.2. Comparison with the Original Port Rotation

In this section, we present the results of the validation of the proposed model and the solution approach. The results obtained using the parameters in Section 4.1 are shown in Tables 4 and 5. It is noticeable that the arrival times do not violate port time window constraints. The solution approach is feasible and optimal. We further compared our model [M2] with the model of [3], which optimized sailing speeds without optimizing port rotation, as shown in Table 6. We can see that it takes nearly 400 h for a container ship to wait in the anchorage in a round-trip journey using the model in [3]. Hence, five container ships need to be deployed on the route to meet the weekly service frequency, while only two container ships need to be deployed following the proposed model [M2]. The total cost is also reduced by 45% compared with that in [3]. Therefore, the proposed model and the solution approach can greatly reduce the total cost.

4.3. Effects of Port Time Windows

In this section, we discuss the impact of port time windows on the solution. First, we examine the effect of relaxing the time windows on the results. We set the time windows to three scenarios, where scenario 1 is the setting in section 4.1. As shown in Table 7, in scenario 2, we give an extra day for the time window of each port while keeping the time windows of scenario 1; we did the same in scenario 3 on the basis of scenario 2. The results are shown in Table 8. It is noticeable that the total cost decreases and that the port rotation is different. This is because, with the relaxation of the time windows, there will be a better port rotation to reduce the total cost.

Then, we examine the effect of different time windows on the results. Scenarios 4 and 5 are set to different time windows from that of scenario 1, and each of them has only 1 day available in a week on each port, as shown in Table 9. From the results in Table 10, we can see that, under different time window restrictions, both the optimal port rotation and the total cost are different. The results also confirm the applicability and flexibility of our proposed models and algorithms, which always seek the optimal port rotation subject to the constraints of the time windows.

4.4. Effects of Bunker Prices

The bunker price fluctuates, and the fuel cost constitutes the bulk of the total operating expenses. Hence, in this section, we examine the impact of bunker price on the solution. We set the bunker price from US$200 to US$600, and the other parameters are the same as those in Section 4.1. The results are shown in Table 11. We can see that, when the bunker price rises, the port rotation remains the same, while the total cost increases significantly. Moreover, with the increase in the bunker price, one more container ship is deployed to maintain the weekly service frequency. The reason is that when more ships are deployed, the round-trip time is longer, and hence ships can sail at a lower speed to reduce bunker consumption.

4.5. Effects of Port Efficiency

In this section, we discuss the effect of port efficiency on the solution. We varied the dwell time at each port by −50%, −25%, +25%, +50%, and +100% of the benchmark value, but the other parameters remained the same. The results are shown in Table 12. We can see that, with the decrease in port efficiency, the total cost increased significantly. The port rotation changed only when the dwell time was +100% of the benchmark value, and in other cases, it remained the same. This indicates that the port rotation under the benchmark value cannot meet the limit of the time windows. Hence, the optimal solution approach gives a new port call order.

5. Conclusions and Future Work

In this paper, we investigated a new and practical tactical-level decision support system that simultaneously optimizes the liner shipping route and ship schedule designs with time windows for the liner shipping industry. Owing to the constraint of the time windows, the liner shipping route design and the schedule design cannot be dealt with separately. Otherwise, the designed schedule may result in a delay in delivery and a huge increase in operating costs. To deal with the problem, we developed a mathematical programming model to simultaneously determine the optimal port rotation, the optimal arrival time at each port of call, and the optimal sailing speed on each leg, with the objective of minimizing the total cost of the shipping route. The nonconvex and nonlinear model is then transformed into a mixed-integer linear programming model by reformulating the piecewise functions from the reconstruction of the hard time windows and by approximating the convex nonlinear function through an outer linear approximation technique. The new model can be solved efficiently by using the current mainstream optimization solver to obtain the optimal solution, which proves the applicability of the model to the problem.

The model was applied to a real-world case study involving SITC Container Lines. A number of numerical experiments were conducted, and a sensitivity analysis was performed on the bunker price, port efficiency, and port time windows. First, a higher bunker price and a lower port efficiency result in a higher total cost, and the changes barely affect the port call sequence. Second, with the relaxation of the time windows, a better port rotation with a lower total cost is obtained.

Our model can also accommodate the possibility of violating time windows with only minor modifications. Actually, can be deemed as a high penalty cost using a premium berth. If is set to be an appropriate penalty cost, a ship can arrive at the port outside time windows by paying a penalty cost that is equal to . Then, the objective is to minimize the sum of the vessel operating cost, fuel cost, waiting cost, and penalty cost.

In the future, first, we will extend the proposed model by incorporating the inventory cost. Second, as this study focuses on a single shipping route, future research could extend the investigation of the problem to a liner shipping network composed of multiple shipping routes, on each of which a port can be visited twice a week. Third, future research could incorporate the proposed problem into a fleet deployment problem while considering decisions making for coordinating the contractual and spot voyages assignment [28, 29]. Fourth, the convergent speed of the solvers is generally slow with large-scale cases, so we will attempt to use appropriate heuristic algorithms to accelerate solution efficiency.

Data Availability

All data included in this study are available from the corresponding author upon request.

Conflicts of Interest

The authors declare that they have no conflicts of interest regarding the publication of this article.

Acknowledgments

This research was supported by the National Key Research Development Program of China (2018YFB1600900), National Natural Science Foundation of China (72071041 and 71874067), Natural Science Foundation of Zhejiang Province (LGG18E080005), the “Six Industry Talent Peak Project” Fund of Jiangsu Province (RJFW-049 and JNHB-115), and the “Green and Blue Project” Fund of Jiangsu Province (2017SJB1641).